Computer ‐ Aided Design of Peptidomimetic Inhibitors of Falcipain ‐ 3: QSAR and Pharmacophore Models

: In this work, antiparasitic peptidomimetics inhibitors (PEP) of falcipain ‐ 3 (FP3) of Plas ‐ modium falciparum ( Pf ) are proposed using structure ‐ based and computer ‐ aided molecular design. Beginning with the structure of Pf FP3 K11017 complex (PDB ID: 3BWK), three ‐ dimensional (3D) models of FP3 ‐ PEPx complexes with known activities ( IC (cid:2873)(cid:2868)(cid:2915)(cid:2934)(cid:2926) ) were prepared by in situ modi ‐ fication, based on molecular mechanics and implicit solvation to compute Gibbs free energies (GFE) of inhibitor ‐ FP3 complex formation. This resulted in a quantitative structure–activity rela ‐ tionships (QSAR) model based on a linear correlation between computed GFE ( ∆∆ G (cid:2913)(cid:2925)(cid:2923) ) and the experimentally measured IC (cid:2873)(cid:2868)(cid:2915)(cid:2934)(cid:2926) . Apart from the structure ‐ based relationship, a ligand ‐ based quantitative pharmacophore model (PH4) of novel PEP analogues where substitutions were di ‐ rected by comparative analysis of the active site interactions was derived using the proposed bound conformations of the PEPx. This provided structural information useful for the design of virtual combinatorial libraries (VL), which was virtually screened based on the 3D ‐ QSAR PH4. The end results were predictive inhibitory activities falling within the low nanomolar concentration range.


Introduction
Malaria is a widespread disease, with causative agent Plasmodium falciparum (Pf), transmitted mainly by female Anopheles mosquito bites [1]. The disease has been declared a public health concern by the World Health Organization (WHO) in many developing countries [2,3]. Additionally, since the implementation of artemisinin-combined therapy (ACT) in 2006, resistance cases have been recorded [4][5][6][7]. Meanwhile, the treatment of malaria mainly depends on ACT, despite resistance to this combination. This suggests the need for industry-academia partnerships for the search of new antimalarials which act via alternative modes of action. Two strategic approaches have been suggested in the search for new remedies against malaria; one focused on eliminating the parasite or preventing its contact with potential human hosts, and a second aimed at developing efficacious drugs to treat infected patients [1]. The latter is often aimed at the inhibition of a therapeutic target, often a vital enzyme involved in the parasite's life cycle. This often requires the search for or the design of new molecules capable of binding in a specific manner to known parasite vital enzymes.
During the last two decades, the identification of drug targets against Pf has increased tremendously [4][5][6][7], thus favouring the second approach. This is known as "rational drug design and discovery". As an example, the parasite breaks down a large amount of haemoglobin (Hb) from human red blood cells in order to obtain the required nutrients for its growth during the blood stage [8]. This involves several proteases, known as validated drug targets in [9]. These drug targets could be divided into two major groups: (i) those which are directly involved in the invasion and rupture of the red blood cells, and (ii) those dedicated to the breakdown of Hb [10]. Two protease families are involved in Hb breakdown by hydrolysis. These include aspartic proteases (plasmepsins) and cysteine proteases (falcipains, FPs) [10]. One metalloprotease called falcilysin [11], and one dipeptidyl aminopeptidase [12] are also involved.
Previous studies have focused on the search for inhibitors of falcipains 2 and 3 (FP2 and FP3), respectively [13], even though FP3, shown to be expressed later in the parasite life cycle, appeared to be a more efficient haemoglobinase than FP2 [14]. This indicates that FP3 inhibition is lethal to the parasite and, therefore, constitutes an attractive target in Pf drug discovery. Several FP3 inhibitors have been identified and described in the literature, which are capable of blocking the enzyme's activity by forming reversible or irreversible covalent bonds within the enzyme active site [12]. These inhibitors could be sub-classified into three categories: peptide-based, non-peptidic, and peptidomimetic inhibitors [15,16], although preference has been given to those known to be reversible and, hence, considered to be potentially more effective than irreversible ones [17,18]. The most promising inhibitors so far are those discovered by chemical synthesis [19][20][21][22][23][24], by molecular docking [25] and virtual screening studies [18,[26][27][28][29], particularly from the peptidomimetics class of compounds.
Weldon et al. recently designed, synthesised and evaluated a series of peptidomimetic pseudo-prolyl-homophenylalanyl ketones for their inhibition of the Pf cysteine proteases FP2 and FP3 [24]. One of these compounds showed nanomolar range activities against both enzymes (i.e., 80 nM against FP2 and 60 nM against FP3 [24]. These interesting results have been improved by the presence of the crystal structures of the FP3 apo structure co-crystallised with the inhibitor within the protein data bank [30,31]. These have constituted the foundation of this work, which involves the design of PEP2 peptidomimetic analogues with the goal of identifying even more potent candidates via quantitative structure-activity relationship (QSAR) with FP3 inhibition pharmacophores. This is intended to further orientate the design of more potent non-peptidic FP3 inhibitors.
In the present work, we have built and validated a Hansch-type "complexation"FP3 inhibition QSAR models based on the in vitro activities of twelve (12) selected PEP derivatives against FP3. As a starting point, we chosethe experimental (X-ray crystal) structure of the protein-ligand complex of the enzyme and the potent inhibitor K11017 (PDB ID:3BWK) to build each selected inhibitor by in situ modifying of the native ligand. [24,30]. This consisted of computing the Gibbs free energies for the formation of the ligand-receptor complexes ∆∆G based on Molecular Mechanics Poisson-Boltzmann (MM-PB) approach for the training set molecules, followed by the correlation with the experimentally tested biological activities pIC . The established QSAR equation was then used to predict the activities of newly designed analogues based on the initial compound scaffold. Additionally, a FP3 inhibition pharmacophore model (PH4) from the bound conformation of the training set of PEPs was used to screen the virtual library of proposed PEP analogues to identify the best candidates, which have predicted ADMET profiles within the acceptable range for 95% of known drugs.

Materials and Methods
Scheme 1 displays the workflow of different steps involved for the computer-aided drug design of the new peptidomimetics analogues. Scheme 1. Novel peptidomimetic analogues design methodology workflow.

Biological Activities of Compounds Included in the Training and Validation Sets
The biological activities (IC ) of the compounds included in training and validation sets of PEP PfFP3 inhibitors were found in the literature, covering a range of activities from 60 nM to about 47,230 nM [24]. Weldon and co-workers synthesised 22 molecules, but not all showed detected biological activities against FP3 to be included in a QSAR study (e.g., activities recorded as >50 μM would not be included in our study). Finally, 12 compounds (almost the threshold for an acceptable QSAR study) have been used in our study.

Molecular Modeling
3D models of the enzyme-inhibitor (E:I) complexes were built starting from the free enzyme (E) and the free inhibitors (I), both derived from a well -refined X-ray crystal structure (PDB ID: 3BWK) of the co-crystallised potent inhibitor K11017 (or Mu-Leu-Hph-VSPh, where VSPh: phenyl vinyl sulfone; Hph:homophenylalanyl;Mu:morpholino urea) retrieved from the PDB [32]. Chain A was employed in all computations and modellingwas carried out on the graphical user in-terface of Discovery Studio 2.5 [33], using a previously well described protocol [34][35][36][37][38][39][40][41][42][43][44][45][46][47][48]. Thus, the pH values were kept at 7.0, while all N-and C-terminal residues were kept neutral. All water molecules originally in the crystal structure were deleted. Finally, protonated and ionised amino acid residues were charged. Each inhibitor was built into the crystal reference structure by modifying the original K11017 inhibitor in situ. During this process, all rotatable bonds of the replacing residues underwent an exhaustive conformational search by a careful and gradual energy minimisation of each modified inhibitor within the active site residues of FP3 within 5Å of the inhibitor, leading to the identification of low-energy bound conformations of each modified inhibitor. The various low-energy structures of the E:I complexes were then carefully refined by energy minimisation of the whole complex.

Molecular Mechanics
The simulation of each inhibitor, FP3 and E:I complex were carried out by molecular mechanics (MM) as implemented in the CHARMm forcefield [49]. All MM calculations used a dielectric constant of 4 for representing dielectric shielding effects in the proteins. The optimisation (energy minimisation process) of the free enzymes E, free inhibitors I and enzyme-inhibitor complexes E:I were carried out by a gradual relaxation of the structures, beginning by adding H-atoms, then the residue side chain heavy atoms, and ending up with the relaxation of the protein backbone. A large number of the steepest descent, followed by conjugate gradient iterative cycles were employed. A convergence criterion for the average gradient was set to of 0.01kcal. mol Å in each geometry optimisation procedure.

Conformational Search
The conformation of each free inhibitor was obtained from its bound conformation in the E:I complex, which had been previously obtained by the gradual relaxation to the nearest local energy minimum (see Section 2.3). Next, low energy structures of the free inhibitors were found by the quenching dynamics protocol available in the module Forcite Plus of Accelrys Materials Studio 4.4 [50]. Quench molecular dynamics performs a standard molecular dynamics calculation with an additional geometry optimisation step, in which a geometry optimisation is performed on every frame in the trajectory file. Forcite Plus calculations were carried out using Compass forcefield [49], ultra-fine quality options and nonbond cut-off distance equal to 30 Å.For each free inhibitor, 5000 steps (time step is 1 fs, total simulation time equal 5 ps) are used to run dynamics simulation at 350 K. A quench, or geometry optimisation, is performed every 25 steps. On completion of the quench dynamics calculation, 200 unique conformations are generated per inhibitor. Finally, the lowest energy conformer of each inhibitor is selected and minimised again using CHARMm forcefield of Discovery Studio. During this minimisation, the inhibitor's dielectric constant was kept at ε = 4.

Solvation Gibbs Free Energy
The electrostatic component of solvation Gibbs free energy was computed using the DelPhi package in Discovery Studio [33]. This incorporates the effects of ionic strength by solving the nonlinear Poisson-Boltzmann equation [51,52]. This DelPhi program treats the solvent as a continuous medium of high dielectric constant ( 80) while the solute is treated as a cavity with low dielectric ( 4). Boundaries are linked to the solute's molecular surface, which enclose the solute's atomic charges. The molecular electrostatic potential and reaction field around the solute are solved by a finite difference method on a 235 235 235 cubic lattice grid for the complexes and free enzyme and 65 65 65 grid for the free inhibitors, implementing the full Coulombic boundary conditions. In both cases, two (02) subsequent focusing steps led to a similar final resolution of about 0.3 Å per grid unit at 70 % filling of the grid by the solute. Physiological ionic strength of 0.145 mol dm , atomic partial charges and radii defined in the CHARMm parameter set (Biovia DS) and a probe sphere radius of 1.4 Å were used. In this implementation, the electrostatic component of the solvation Gibbs free energy was calculated as the reaction field energy [53][54][55].

Calculation of the Entropic Term
During the simulation, the vibrational entropy change which occurs as the inhibitor binds to the enzyme was computed by normal mode analysis of the inhibitor vibrations, by using a simplified method previously described by Fischer and co-workers [56,57]. This approach involves the vibrational analyses of the inhibitor bound at the active site of a "frozen" enzyme, while the low-energy conformers of the free inhibitor were computed for fully minimised structures. This was carried out using the Discover module of Materials Studio 4.4 [50] and the formula: The method has previously provided a good approximation of the vibrational entropy change in the fully flexible system for small and relatively stiff ligands, i.e., including the degrees of freedom of the protein receptor [55,56]. In this calculation, the TS I term can explain vibrational motions of the free inhibitor and the conformational flexibility of the molecule, i.e., low frequency vibrations, which correspond to collective motions of atoms with larger amplitudes (conformational changes contribute mostly to this term). The relative values of ∆∆TS , with respect to the reference inhibitor, were used to compensate partially for the restricted flexibility of the E. In this respect, the entropy term is also recognised as an important factor for drug optimisation, even though an enthalpy contribution to binding affinity is known to be more essential [58].

Binding Affinity Calculations
It has been previously proven that the concentration of a competitive tight binding inhibitor that causes a 50% reduction in the rate of catalytic substrate conversion (IC ) of a reversible inhibitor depends on the enzyme-inhibition constant K as follows: where S and E are the substrate and enzyme concentrations, respectively, while K represents the Michaelis constant [59]. The IC value can thus be predicted from the standard Gibbs free energy (GFE) change during the enzyme:inhibitor complex formation: assuming that there is equilibrium in solution between the solvated protein (or enzyme) E , the solvated inhibitor I and the solvated protein-ligand complex E: I : the standard Gibbs free energy change in the above equilibrium (4) can be written as: in our calculations, the exact values of standard Gibbs free energies for larger systems such as enzyme: inhibitor complexes were approximated by the derived expressions from the works of Frecer and Miertus [36,60,61]: where E E: I stands for the molecular mechanics total energy of the complex (including bonding and non-bonding contributions), G E: I is the solvation GFE and TS E: I is the entropic term :  TS E: I  TS  E: I  TS E: I  TS E: I  (7) composed of a sum of contributions arising from translational, rotational and vibrational motions of the complex E:I. Assuming that the translational and rotational terms for the complex E:I and free enzyme E are approximately equal, we obtain: I  TS  I  TS I  TS E: I  TS E  TS I   ∆H  TS  I  TS I ∆TS ∆G with TS I and TS I describing the translational and rotational entropy terms of the free inhibitor, respectively, and ∆TS representing the vibrational entropy change during the formation of the enzyme-inhibitor complex ∆TS TS I TS I . By comparing between different inhibitors via relative changes in their respective complexation Gibbs free energies, with respect to a reference inhibitor, I (in this case PEP23), and by assuming the ideal gas behaviour for the rotational and translational motions of the inhibitors, it can be shown that: The advantage of such an approach is that the evaluation of relative changes is preferable, since it results in the partial cancellation of errors caused by the approximate nature of the MM method. Additionally, solvent and entropic effects are included in the description.

Interaction Energy Calculations
Interaction energy values were computed using Discovery Studio 2.5 [33]. The MM interaction energy (E ) protocol available in this program computes the (non-bonded) van der Waals and electrostatic interactions between enzyme residues and each inhibitor. The CHARMm force field [49] was used during the calculations, with a dielectric constant set at 4. The breakdown of E into the contributions by active site residues reveals the significance of individual interactions and permits us to carry out a comparative analysis. The approach leads to the identification of affinity values which would enhance the prediction of favourable and unfavourable PEP substitutions.

Pharmacophore (PH4) Modeling
By definition, a pharmacophore is often regarded as a set of features arranged in 3D space which are essential for a molecule to exert a certain biological activity. The perception of a pharmacophore is essential for understanding the interaction between a ligand and its receptor. The PH4 concept is based on the assumption that a set of structural features in a molecule is recognised at the receptor site and is responsible for the molecule's biological activity. Bound conformations of inhibitors taken from E:I complexes were used to construct 3D-QSAR pharmacophore models using the Catalyst HypoGen algorithm implemented in Discovery Studio 2.5 [33]. This consisted of building a top scoring pharmacophore hypothesis from the most active inhibitor. Three stages (construction, subtraction and optimisation) are involved, meanwhile the inactive ones were used to define the excluded volumes. A maximum number of five excluded volumes and five features were selected according to the PEP scaffold and substituents, i.e., hydrophobic aliphatic (HYd), hydrophobic aromatic (HYdAr), hydrogen-bond acceptor (HBA) hydrogen-bond donor (HBD) and ring aromatic (Ar). As per the adjustable parameters in the HypoGen protocol, all were kept by default except for the uncertainty on the activity and the minimum inter-feature distance, which were set to 1.1 Å and 2.5 Å (instead of 3), respectively. These parameters were carefully chosen in order to bring the uncertainty interval of experimental activity from a wide span IC 3 ⁄ , 3 IC to a relatively narrow one IC . This is important because the accuracy and homoge-neity of the measured inhibitory activities based on the fact that the experimental biological activities were derived from the same laboratory must be taken into account [24]. During the generation of 10 pharmacophores, 0 was set as the number of missing features and the best pharmacophore models were selected.

Generation of the Virtual Library
Molecular models of new analogues were generated using the Molecular Operating Environment (MOE) program [62]. This was carried out by attaching the R-groups (fragments, building blocks) onto the PEP scaffolds using the Quasar CombiDesign module of the MOE program. Chemical reagents considered in this study were taken from the directories of chemicals available from the commercial suppliers of chemicals [63]. Each analogue was built as a neutral molecule in the MOE program and its molecular geometry has been refined by the MM optimisation, implemented in the Discovery Studio 2.5 smart minimiser. Convergence criteria (energy difference of 10 kcal. mol , root-mean-square displacement (RMSD) of 10 −5 Å and a dielectric constant of 4 using the CHARMm force field were set, as described in Section 2.3.

In Silico Screening
The conformers with the best mapping on PH4 pharmacophores in each cluster was selected for virtual screening using the complexation QSAR model. For each E:I complex, the relative complexation Gibbs free energies ∆∆G was calculated. This was then used to compute the predicted activities (pIC of each of the newly designed analogues against FP3. The IC values were then calculated using the formula IC 10 .

Selection of Training and Validation (or Test) Data Sets
A data set of ten (10) FP3 inhibitors with a broad range of in vitro activities IC , obtained from the same laboratory, with a sufficiently broad range of activities (60-47,230nM) [24] were used to generate a 3D-QSAR model. This was divided into a training set of nine (9) inhibitors used to build the QSAR model and a validation set of three (3) inhibitors for evaluating the model (Table 1 and Figure 1).

Obtained QSAR Model
The relative Gibbs free energy of the non-covalent enzyme-inhibitor (E:I) complex formation from free enzyme (E) and free inhibitor (I), shown in the Experimental Section, were computed for each FP3-PEPx prepared complex. This was carried out by modifying in situ of the inhibitor K11017 within the binding site of FP3 of the refined crystal structure, with PDB ID: 3BWK [30,31]. Table 2 provides the computed values of complexes formation GFE ( ∆∆G ) and its components (see Experimental Section). Since the ∆∆G values were computed in an approximate way, the relevance of the binding model is evaluated by correlating it with the experimental activity data (IC ) using linear regression; see Table 3. For this training set, a plot of the linear correlation is shown in Figure 2 and the statistical data of the regression are provided in Table 3. For the correlation involving ∆∆G , the relatively high regression coefficient on the values, together with the statistical significance Fischer F-test, suggest that there is no chance correlation between the binding mode and the observed inhibitory potencies of the training set. [e] ∆∆G represents the relative GFE change related to the enzyme-inhibitor complex formation: (see Equation (9)). [f] IC [24] represents the inhibitor concentration that causes 50% decrease in the rate of substrate conversion by FP3 measured in the enzyme assay: pIC log .
[g] This is the ratio of the predicted activity on the experimental activity, pIC pIC . This ratio is close to 1, indicating the predictivity of the QSAR model.  The ratio of the predicted and observed inhibition constants pIC pIC for the validation set of three PEPs (not included in the training set) were closed to 1. This proves the predictive power of the QSAR model, suggesting that the regression equation (A) ( Table 3), and the computed ∆∆G quantities of the newly designed PEP analogues can be used to predict their inhibitory potencies (IC ) against FP3, on condition that the binding modes of the designed analogues and those of the training set compounds are the same relative to the receptor site. Such an approach could reduce the required number of molecules to be synthesised in a rational drug development project quite considerably. The above procedure has been previously applied by our team in several drug design projects [34,36,[38][39][40][41][42][43][44][45][46][47].

Inhibitor Binding Modes
The predicted binding mode of the best active PEP39 coming from the complexation model is illustrated in 3D depiction in Figure 3. The main interactions with the active site residues, namely the H-Bond with the catalytic residue Cys51, are in line with the docking study and WaterMap calculations [24] which, unfortunately, did not provide any statistical correlation between binding affinity and activity (results not shown). The bound conformation of PEP sheds light on the structural features for binding affinity, which are vital for the design of novel potent non-peptidic FP3 inhibitors by exploiting the S1' to S3 pockets (Figure 4) [30]. In order to verify whether other interesting interactions not displayed have to be taken into account in the description of PEP binding mode at the FP3 active site for the rational design of new analogues, the interaction energy (IE) between each active site residue and PEPx was computed. The breakdown of the interaction energy diagram into each S1′-S3 subsite residue contribution of FP3 for PEPs, displayed in Figure 3, indicates the highest interacting residues of the overall active site of FP3. Moreover, the predicted binding mode of PEP inhibitors highlights three main fa-vourable non-bond interactions ( Figure 4)   It was observed that the IE diagrams analysis could not significantly guide the choice of the R-groups in S1' and S2 subsites [30], when compared with the case for the design of thymine-like inhibitors of thymidine monophosphate kinase [41]. It would rather be suggested that a large and diverse combinatorial virtual library (VL) of PEPs be built and screened with our FP3 inhibition 3D-QSAR PH4, based on the complexation one descriptor QSAR model. A successful case study was the design of pyrrolidine carboxamide inhibitors of Mycobacterium tuberculosisInhA [42]. S1 S S2 S3

Ligand-Based 3D-QSAR PH4 Model of FP3 Inhibition
The 3D-QSAR PH4 pharmacophore generation process follows three main steps, namely the constructive, the subtractive and the optimisation steps [33]. The constructive phase of Hypo-Gen has automatically selected the most active compounds for which IC 1.1 60 nM as leads. Thus, only the most active compound PEP39 (IC = 60 nM) was used to generate the starting PH4 features. Only those features that matched this lead were retained. In the subtractive phase, which is normally used to remove pharmacophoric features present in poorly active molecules, none of the training set compounds were found to be inactive IC 60 10 . 189,736 nM). During the optimisation phase, the score of the pharmacophoric hypothesis is improved. Hypotheses are scored according to errors in activity estimates from regression and complexity via a simulated annealing approach. At the end, the top scoring 10 unique pharmacophoric hypotheses (Table 4) were kept, all displaying four features. The generated pharmacophore models were then assessed for their reliability based on the calculated cost parameters. The overall costs ranged from 24.13 (Hypo1) to 456.44 (Hypo10). The relatively small gap between the highest and lowest cost parameter corresponds well with the homogeneity of the generated hypotheses and the consistency of the training set. For this PH4 model, the fixed cost (21.24) is lower than the null cost (2317.26) by a difference ∆ 2296. This difference is a major quality indicator of the PH4 predictability (∆ > 70 corresponds to an excellent chance or a probability higher than 90% that the model represents a true correlation [33]). To be statistically significant, the hypotheses have to be as close as possible to the fixed cost and as far as possible from the null cost. The cost distance ∆ ≥ 1860 for the set of 10 hypotheses confirms the high quality of the pharmacophore model.  [c] overall cost parameter of the PH4 pharmacophore; [d] cost difference between null cost and hypothesis total cost; [e] lowest cost from 49 scrambled runs at a selected level of confidence of 98%. The fixed cost = 21.24, with RMSD = 0, the null cost = 2317.26, with RMSD = 22.65 and the configuration cost = 11.85.
The standard indicators such as the RMSDs between the hypotheses ranged from 0.79 to 9.83 and the squared correlation coefficient (R 2 ) falls to an interval from 0.90 to 0.99. The first PH4 hypothesis with the best RMSD and R 2 was retained for further analysis. The statistical data for the set of hypotheses (costs, RMSD, R 2 ) are listed in Table 4. The geometry of the Hypo1 pharmacophore of FP3 inhibition is displayed in Figure 5. Table 5 lists the regression equation (Table 3) for pIC vs. pIC estimated from Hy-po1 with related indicators such as R , R , Fisher F-test, and , while Figure 6 displays a plot of regression equation for pIC vs. pIC . To check the consistency of the generated pharmacophore model, we have computed the ratio of predicted and observed activities pIC / pIC ) for the validation set. The computed ratios are as follows: PEP26 (1.01), PEP28 (1.01), PEP36 (1.01) all of them relatively close to one, which documents the substantial predictive power of the regression for the best PH4 model.   The configuration cost (11.85 for all hypotheses) far below 17 confirms this pharmacophore as a reasonable one.
The evaluation of Hypo 1 was performed first through Fisher's randomisation cross-validation test. The Cat-Scramble program was used to randomise the experimental activities of the training set. At 98% confidence level, each of the 49 scramble runs created 10 valid hypotheses, using the same features and parameters as in the generation of the original 10 pharmacophore hypotheses.
Among them, the cost value of Hypo1 is the lowest compared with those of the 49 randomly generated hypotheses, as we can see in Table 4, where the lowest cost of the 49 random runs is listed for each original hypothesis, and none of them was as predictive as the original hypotheses generated shown in Table 4. Thus, there is a 98% probability that the best selected hypothesis Hypo1 represents a pharmacophore model for inhibitory activity of peptidomimetics with a similar level of predictive power as the complexation QSAR model, which relies on the PEPx active conformations from 3D structures of the FP3-PEPx complexes and computed GFE of enzyme-inhibitor binding ∆∆G . Another evaluation of Hypo 1 is the mapping of the best active training set PEP39 ( Figure 5) displaying the geometry of the Hypo1 pharmacophore of FP3 inhibition. The regression equation for pIC vs. pIC estimated from Hypo1: pIC 1.0002 pIC 0.0012 (n = 9, R = 0.99, R = 0.99, F-test = 5675.56, σ = 0.04, α > 98 %) is also plotted on Figure 5.

Library Design and ADME Focusing
In order to identify more potent PfFP3 peptidomimetics inhibitors, we have built a virtual library of new analogues inhibitors of PfFP3 based on substitutions at four positions (P1′, P1, P2 and P3) of a scaffold of a dipeptidic compound in order to better FP3 active site four pockets S1′, S1, S2 and S3 [30]. This virtual library was built from the side chains of 20 essential amino acids except the proline side chain. The 19 R-groups listed in Table 6 have been attached in positions R1 to R4 of the appropriate scaffold to provide a combinatorial library of the size: R R R R 19 130,321 PEPAs. It should be noted that one of the important criteria for the design of new antimalarials, is their oral bioavailability in the context of oral delivery. Structural information provided from these peptidomimetic computational studies can guide the design of novel antimalarial inhibitors of FP3 deliverable by injection. Therefore, no ADME-based focusing step in order to remove compounds with expected poor oral bioavailability and low drug likeness was performed for the enumerated combinatorial library [64].

Screening PEPs Virtual Library Using the Obtained in Silico Model
The library of PEP analogues has been further virtually screened for molecular structures matching to the 3D-QSAR PH4 pharmacophore model Hypo1 of FP3 inhibition. During this virtual screening, 1000 conformations were generated for each analogue using the BEST algorithm of Discovery Studio 2.5. Thus 130,321,000 conformations were screened to fit the 3D-QSAR PH4 pharmacophore model Hypo1 retained in this work. From the set of 130,321,000 analogues, few thousands of PEPAs mapped to at least two features, 242 of which mapped to four features of the pharmacophore.Out of then, only the 21 best-fitting analogues (PH4 hits) have been retained and submitted to screening with help of the complexation QSAR model. Their Gibbs free energy (GFE) upon complex formation with PfFP3 has been computed along with its component and their predicted half-maximal inhibitory concentration IC has been estimated with the correla-tion the equation on Table 3. The results obtained are given in Table 7. Of the 21 analogues whose inhibitory activities were predicted in Table 7, 13 showed better activities than the most active compound of the training set among them four showed even more activity: PEP-17-03-14-10 ( IC = 0.29 nM); PEP-08-15-18-19 ( IC = 0.19 nM); PEP-13-06-04-19 (IC = 0.10 nM); PEP-14-14-14-18 (IC = 0.07 nM). Table 7. Complexation Gibbs free energy and its components for the top 21 scoring virtually designed analogues. The analogue numbering concatenates the index of each substituent R1 to R4 numbered in Table 6. represents the relative Gibbs free energy change related to the enzyme-inhibitor FP3-PEP complex formation (see Equation (9)); [g] IC represents the predicted inhibition constant towards PfFP3 calculated from ∆∆G using correlation (Table 3).

Analysis of New Inhibitors
In order to identify the substituents that make the analogues predicted to be active, we have analysed the frequency of occurrence of certain substituents chosen from Table  7, on the predicted active analogues. From the four best analogues proposed (seen chemical structure in Figure 7), the following R-groups are present 3, 4, 6, 8, 10, 13, 14, 15, 17, 18 and 19. Additionally, Figure 8 displays the best virtual hit, analogue PEP-14-14-14-18 and the least active PEP32 mapped to a PH4. Figure 9   The FP3-PEPs' interaction energy breakdown to active site residues displayed on Figure 10 is classified according to the enzyme's four pockets S1, S2, S3 and S1′ for the top four analogues (PEP-Top4). The interaction energy between FP3 and the most active training set compound PEP39 (IC = 60 nM) breakdown is presented on Figure 3 where only residues with noticeable contributions are displayed. The sum of residues' contribution to E for the S1 pocket is almost the same for PEP39 and PEP-Top4 as is the case for the S1′ pocket, where a slight difference of 1 kcal mol in favour of PEP-Top4 is noticed. For pocket S2 PEP-Top4 the sum of energy is lower by about 5 kcal mol compared with PEP39. The same stabilising effect of 6 kcal mol in favour of PEP-Top4 is detectable for pocket S3.

Conclusions
Structural information from the crystal structure of the FP3-K11017 complex has been successfully used to establish a reliable QSAR model of non-covalent PfFP3 inhibition by peptidomimetic (PEP) inhibitors. This model correlates the unique descriptor, namely the computed Gibbs free energies (GFE) upon complex formation, with observed inhibitory potencies and is able to identify a few predicted low nanomolar range inhibitors of P. falciparum. As GFE is a combined descriptor involving the enthalpic gas phase, entropic contributions and solvation free energy, a precise insight into S1' and S3 pockets filling has been performed from the model by analysis of interactions between the enzyme active-site residues and the inhibitor. For this purpose, the breakdown of the interaction energy clearly indicated the residues involved in the affinity with the most active inhibitors. This information has helped to design an initial diversity virtual combinatorial library of new analogues to be screened by the pharmacophore models derived from the GFE QSAR. The screened library by mapping of the analogues to the PfFP3 inhibition PH4 pharmacophore permitted a library subset of 21 best virtual hits to be selected, which was further submitted to the computation of predicted PfFP3 inhibitory potencies by the formerly established complexation QSAR model. The best cross checked analogues showed predicted activities in the low nanomolar concentration range, with the most promising hits being PEP-08