Identification of Flavonoids as Inhibitors of Plasmodium falciparum Enoyl-ACP Reductase by Hierarchical Virtual Screening

Malaria, a parasitic infection caused by Plasmodium falciparum, is a serious global public health problem. The enzyme enoyl-acyl carrier protein (ACP) reductase present in a limiting step in the biosynthesis of fatty acids type II was used to search for novel molecules with potential inhibitor by virtual screening techniques. A flavonoid library obtained in ChEMBL database (n = 4,008) was filtered through physico-chemical similarity using the Euclidean distance as a criterion. The selected molecules were subjected to molecular docking using the scoring function GridScore in the DOCK 6.5 software. The top 30 ranked molecules by molecular docking were submitted to the AuPosSOM 2.1 software with the aim of clustering molecules by means of intermolecular interactions in order to understand the patterns of these interactions with enoyl-ACP reductase from Plasmodium falciparum (PfENR). We observed the importance of flavonoid moiety to hydrophobic interactions with the active site of enzyme, this confirmed the choice to study this metabolite class as PfENR inhibitors. Moreover, it was possible to recognize which intermolecular interactions contribute to the molecular recognition process.


Introduction
Malaria is one of the parasitic infections that causes most concern and can be considered a serious global public health problem, since it causes many cases of morbidity and mortality. 1,2 The Anopheles mosquitoes infected with Plasmodium protozoa are responsible for transmitting malaria in humans. Only five species promote this infection, P. falciparum being the most dangerous for causing cerebral malaria and in more severe cases leading to death. 3 Although antimalarial drugs play a central role in controlling and eliminating the disease, their activity against Plasmodium sp. showed reduced effectiveness due to resistance that parasites have developed. Plasmodium falciparum has a high capacity to adapt by mutation, causing resistance to antimalarials. Thereby, the development of new therapeutic antimalarial agents is a matter of great interest to public health. 4,5 The synthesis of fatty type II (FAS II) by P. falciparum is of great relevance to antimalarial research because its enzymes are absent in humans and it is part of the limiting step of biosynthesis in the FAS II route, being potential targets for development of new drugs. In addition, studies 6 have shown that inhibition in submicromolar concentrations of the antimicrobial agent triclosan has caused the death of the parasite.
Other reports 7 show that bioactive flavonoids extracted from medicinal plants from Africa have potential antimalarial activity. Therefore, medicinal plants are a potential source of new drugs against malaria, in view of the presence of a variety of metabolites with structural diversity and pharmacological activity. The promise of results in the use of plant materials in the treatment of malaria in African traditional medicine (ATM) has led to the search for natural products as a source of molecules with bioactive potential or as an inspiration for the development of new antimalarials. 8 The action mechanism of flavonoid with antimalarial activity is by the biosynthesis inhibition of FAS II of P. falciparum. 7 The first identified antimalarial natural product was a common glycosylated flavonoid, luteolin-7-O-glucoside, which acted as an inhibitor of enoyl-acyl carrier protein (ACP) reductase from Plasmodium falciparum (PfENR). 5 This underscores the importance of flavonoids as potential inhibitors to PfENR, which can be used as models for the research and development of new antimalarial drugs, through methods to adopt more direct approaches, made from the understanding of the molecular processes involved in disease. 9,10 The use of experimental techniques and a complementary manner cheminformatics provide an increase in the chances of success in the process of discovery of bioactive molecules, since that facilitates the identification of new targets, the elucidation of the mechanism and functions for the development of molecules with more selective properties. 10 The discovery of the activity of a molecule on a specific target and a careful analysis at a molecular level of the main interactions that may occur may be used to streamline the process of interaction between the compound and the receptor. 11 Thus, cheminformatics techniques, such as the virtual screening, were used to find new potential molecules for the treatment of malaria. Thus, flavonoids were evaluated as a promising antimalarial, because of their potential to inhibit the PfENR, present in a limiting step of FAS II of P. falciparum.

Ligand-based virtual screening method
Initially two libraries of chemical structures were constructed, the first with only inhibitors of PfENR (DB1), with half maximal inhibitory concentration (IC 50 ) < 1 μM (Table S1, Supplementary Information (SI) section). The second library was formed by flavonoids (DB2) obtained from ChEMBL database. 12 The data line entry system molecular simplified (SMILES) of the DB1 and DB2 were submitted to the ChemGPS-NP web server, 13,14 to generate the coordinates of each principal component (PC) to locate the molecules in physical-chemical space. The Euclidian distance (ED) was used to group DB1 molecules into clusters (ED < 2). Flavonoids with ED < 5 were selected using the DB clusters with spatial reference.

Structure-based virtual screening method
Nicotinamide adenine dinucleotide (NAD) is a cofactor that is central to metabolism. The molecular docking was performed with crystallographic structure 3LT0 (ligands: 4-(2,4-dichlorophenoxy)-3-hydroxybenzaldehyde (FT-1) and NAD501 cofactor). We used the UCSF Chimera (or simply Chimera) 1.9 software 15 to add polar hydrogens and Gasteiger charge (ff12SB method). For docking studies, the DOCK 6.5 software 16 was used. The choice of the optimal scoring function for the molecular docking was performed using the redocking and enrichment metrics.
The preparation of the binding sites involved the calculation of the solvent accessible surface of the 3LT0 (PDB ID), 17 without hydrogen atoms, using a probe radius of 1.4 Å with the Document Management Systems (DMS) software. 18 A negative image of the surface was generated using the program Sphgen (SPHGEN). 19 A subset of spheres near the crystallographic ligand (8 Å) was selected with the aid of the sphere selector software. A box around the binding site was constructed with the accessory in SHOWBOX software. 20 The box includes the selected spheres and a protocol dependent margin. The grid was computed with the accessory in GRID 20,21 software using a 0.4 Å grid spacing, a 9,999 Å distance cutoff, 6-9 van der Waals (VDW) exponents. The following scoring functions were used: Grid Score, Grid Score + Grid Hawkins GB / SA and Continuous Scoring.
The evaluation of docking parameters with performed using redocking metric, docking success was observed when the top scoring pose is within 2 Å heavy atom rootmean-square deviation (RMSD) of crystallographic ligand pose. 22 The RMSD value was calculated by DOCK 6.5. 16 For enrichment calculations the decoys database was obtained from a Database of Useful Decoys-Enhanced (DUD-E). 23 Receiver operating characteristic (ROC) curves were used to evaluate how well methods favorably rank known active molecules compared with decoys. 24 The analysis of scoring functions was performed from the recovery rate, calculating the area under the curve (AUC) which determines measure of global enrichment of the method. The structure-based virtual screening was performed using the scoring function selected; the results were organized by number.

Contact-based ligand-clustering analysis
The top 30 best results for molecular docking were submitted to an automatic analysis of poses using a selforganizing map with AuPosSOM 2.1 server 25 for the evaluation of docking results based on the clustering of compound by similarity of contacts.

The intermolecular interactions (hydrogen interactions, electrostatic interactions, lipophilic interactions or all)
were done by ROC curve analysis, therefore we used the same set of decoys and PfENR inhibitors used in the evaluation stage for virtual screening. The contact that was most efficient was selected for analysis. The results of ligand-clustering were analyzed by interpreting the Newick tree generated in the Dendroscope 3.0 software. 26 The molecules were separated into clusters based on the analysis of the intermolecular interactions in each subclusters.
The analysis of intermolecular interactions was performed in the online web server Poseview 1.97.0 27,28 for the construction of 2D interactions maps and Discovery Studio Visualizer 4.1 29 for the construction of 3D interactions maps.

Ligand-based virtual screening method
The enoyl-ACP reductase (ENR) inhibitors (n = 32) identified in the literature (Table S1, Supplementary Information (SI) section) were submitted to ChemGPS-NP web server. 13,14 Some inhibitors were located in the nearest regions in the chemical space; therefore, it was decided to group in 5 clusters with different chemotypes. The clusters were compared to a database of flavonoids obtained in ChEMBL (n = 4,008). A total of 3,057 molecules were selected using this chemical filter ( Figure 1).
The three dimensions of the ChemGPS-NP web map can be interpreted in such a way that the first dimension (principal component one, PC1) represents size, shape and polarizability; PC2 corresponds to aromatic and conjugation properties and PC3 describes lipophilicity, polarity, and H-bond capacity. 13,14 The screening by chemical similarity selected molecules with aromatic and polar nature characteristics, supporting both with the features of the inhibitors CHEMBL370595 (IC 50 = 0.18 μM), CHEMBL370448 (IC 50 = 0.81 μM), CHEMBL200492 (IC 50 = 0.56 μM), CHEMBL199081 (IC 50 = 0.59 μM) described in the literature for PfENR 30 which presented in the positive sense of PC2, indicative of larger aromaticity, as well as the chemical properties of flavonoids. [31][32][33][34] Evaluation of the docking parameters The redocking showed that the parameters entered in the DOCK 6.5 16 software were able to accurately predict the pose of crystallographic ligand (4-(2,4-dichlorophenoxy)-3-hydroxybenzaldehyde) ( Figure 2). The 0.54 Å RMSD of the initial pose of the crystal ligand is below the value considered as ideal. 22 The ROC analyses were used to assess the ability of each scoring function to recover the active molecules instead of inactive (false positive). The best recovery rate is one that can recover the greatest amount of true-positive early in the bank compounds and AUC analysis shows best performance scoring function to present value closer to 1, because an AUC = 1 indicates 100% differentiation of active compounds (true positives) over inactive. All score functions showed good value in AUC analysis: Grid Score (AUC = 0.86); Grid Score + Grid Hawkins GB / SA (AUC = 0.76); Continuous Scoring (AUC = 0.83). 35 The Grid Score and Continuous Scoring functions showed very similar behavior, difference between AUC = 0.04, however, the Grid Score recovery rate was higher, 80% of the molecules with the best score were true positives. It is, therefore, the function chosen for structure-based virtual screening with the selected molecules by ligandbased strategies. The Matthews correlation coefficient (MCC = 0.31) proposed by Matthews 36 is considered a  balanced measure that can be used even when the samples are very uneven in size. A correlation C = 1 indicates perfect agreement, C = 0 is expected for a prediction no better than random, and C = −1 indicates total disagreement between prediction and observation. 37

Structure-based virtual screening method
The molecules were ordered by Grid Score values (Table S2, SI section). Grid Score of 30 molecules selected ranged from −80.4 to −51.5 kcal mol −1 . When compared with triclosan (−51.5 kcal mol −1 ) and crystallographic ligand (−54.7 kcal mol −1 ), it was possible to detect all the selected molecules that have an affinity for the enzyme; such fact may increase the possibility as potential inhibitor.

Contact-based ligand-clustering analysis
The AuPosSOM 2.1 25 uses contact analysis approaches of molecules, the mode of interaction with the receptor, and in most cases, structurally similar ligands show high conservation of connection modes to the recipient thereof. AuPosSOM efficiency depends on the set of molecules and molecular target, it should evaluate the pattern of intermolecular interaction that best reproduces the system under study. 25 Through a ROC analysis it was possible to select the pattern of interaction that best fits the set of molecules to be organized by AuPosSOM 2.1. It can be observed that the hydrogen interactions standards (AUC = 0.85), as it evaluates all intermolecular interactions (hydrogen, electrostatic and lipophilic) (AUC = 0.88), describe satisfactorily contacts inhibitors PfENR with the catalytic site of this enzyme. However, the pattern of interaction chosen was the organizing molecules based on all the intermolecular interactions (hydrogen, electrostatic and lipophilic), as well as presenting the best AUC = 0.88 and also describing the types of interactions found in literature [30][31][32][33] that showed hydrogen interactions, electrostatic and hydrophobic between the catalytic site inhibitor and the PfENR.

Virtual testing by self-organizing maps
The top 30 ranked molecules by DOCK 6.5 were submitted to AuPosSOM 2.1 software for carrying out virtual testing in order to group the molecules by the interaction types. The interaction pattern was organized with molecules based on all intermolecular interactions (hydrogen, electrostatic and lipophilic) and taking in consideration the type of interaction between the inhibitors and the PfENR catalytic site in the literature. [30][31][32][33] The tree of Newick (Figure 3) has 10 leaves, and the group of molecules of each leaf has the same pattern of interaction with the receptor. Additionally, some leaves also display the same pattern of intermolecular interactions between them, so based on the proximity of the branches of Newick's tree, the 5 higher leaves were selected for analysis of their intermolecular interactions.
Cluster 1 has a higher number of conserved hydrophobic interactions, which indicates that the orientation of the ligands obtained by molecular docking is prioritizing interactions with residues highlighted in the literature [30][31][32][33] as important for molecular recognition. The presence of aromatic amino acids such as tyrosine (Tyr) and phenylalanine (Phe) in the active site of the enzyme can contribute to the formation of aromatic interactions between inhibitors and PfENR.
Clusters 2 and 5 were the ones who presented hydrogen bonds with Tyr277, this residue performs this same pattern of coordination with triclosan inhibitor. Hydrogen interaction is fairly weak, but it is of utmost importance in biological systems. 38 The hydrogen bonds, hydrophobic and electrostatic highlighted on the 2D map ( Figure 4) illustrate those already listed in Table 1. However, the visualization of 3D interaction maps results in more thorough analysis of the geometry of the interactions and allows you to analyze the interaction of molecule of NAD501 cofactor. The analysis of interactions with the NAD501 is interesting, since these interactions are of great importance for the catalytic activity of PfENR. Several studies [30][31][32][33] describe the ternary complex interactions PfENR-NAD501-triclosan. Thus, a representative molecule of each cluster was selected in order to illustrate the interaction pattern within these clusters.

Intermolecular interactions representative of cluster 1 (CHEMBL128043)
The aromatic region of the CHEMBL128043 forms aromatic-aromatic interaction (face-to-face) with the nicotinamide ring from NAD501. Furthermore, the secondary amide hydrogen attached to the aromatic ring Table 1. Conserved interactions of ligands with PfENR active site residues within clusters and sub-clusters (SC), and the ligand CHEMBL codes. Amino acids in bold refer to amino acids already described in the literature. 30 of the ligand forms a hydrogen bond donor interaction with acetamide of NAD501 nicotinamide ring. Another type of aromatic interaction of great importance for the stability of the complex is the T-shaped aromatic-aromatic interaction with the amino acids Tyr277 and Tyr267. It is important to note that T-shaped aromatic-aromatic interactions are considered important because they favor greater stability in the interaction between ligand and receptor. 39 The CHEMBL128043 flavonoid moiety has hydrophobic interaction PfENR, wherein ring B interacts with Arg318 and Trp131 amino acids, and the region of the rings A and C interact with Val222. Aliphatic regions of the molecule also perform hydrophobic interactions with the amino acids present in the hydrophobic region of PfENR as Asn218, Tyr267, Tyr277, Met281, Pro314, Ala319, Ile323 and Ala372. π-Cation interaction was observed between the rings A and C of CHEMBL128043 with side chain of Asn218 and the secondary amide of Ala219 (Figure 4a).

Intermolecular interactions representative of cluster 2 (CHEMBL1095399)
This compound does not form interactions with the NAD501 cofactor, this can be justified because CHEMBL1095399 has only an aliphatic chain substituted in flavonoid moiety and the interactions described with NAD501 are mostly between an aromatic region of the molecule and NAD501 nicotinamide ring. The B ring region of the CHEMBL1095399 flavonoid moiety performs hydrophobic interactions with residues Arg318 and Trp131 (Table 1). Furthermore, the aliphatic chain of this compound also performs this type of interactions with the amino acids Ala217, Asn218, Ala319, Ala320 and Ile323. In the regions of hydrophobic interactions that characterize the catalytic site of the enzyme PfENR have been described in the literature 32 with a presence of hydrogen bonds with Tyr277, which can be observed between the hydroxyl side chain of Tyr277 and the carbonyl group of the aliphatic region of CHEMBL1095399. The C ring of CHEMBL1095399 molecule shows π-donor aromatic interaction with the primary amide side chain of Asn218 (Table 1).
In the three-dimensional interactions map, π-donor interactions between the C ring of the flavonoid moiety and the CHEMBL1095399 molecule with the primary amide side chain of Asn218 can be observed (Figure 4b).

Intermolecular interactions representative of cluster 3 (CHEMBL1484629)
The aromatic ring of this molecule forms face-to-face aromatic-aromatic interaction with the nicotinamide ring of NAD501. T-shaped aromatic-aromatic interaction occurs between the CHEMBL1484629 aromatic ring and Tyr267 side chain (Table 1). Furthermore, it was also possible to observe the hydrophobic interactions with the flavonoid moiety, wherein ring B of CHEMBL1484629 forms hydrophobic interactions with Trp131 and Ala319, whereas Asn218 and Met281 interact with the rings A and C of CHEMBL1484629. The C ring of flavonoid moiety interacts with Asn218. The 3D map shows π-donor interaction between secondary amide of the Asn218 and the C ring flavonoid moiety (Figure 4c).

Intermolecular interactions representative of cluster 4 (CHEMBL2042175)
The aromatic region of the molecule carries face-to-face aromatic-aromatic interaction with NAD501 nicotinamide ring. T-shaped aromatic-aromatic interaction was observed with Tyr267. The B ring of the CHEMBL2042175 forms hydrophobic interaction with the amino acids Val134, Trp131 and Arg318. The region of the rings A and C interact with Val222 and Asn218. The CHEMBL2042175 also formed hydrophobic interactions with Tyr277, Ala319 and Ile323 (Table 1). The three-dimensional map also displays possible interactions of π-donor type between the rings A and C of flavonoid moiety of the CHEMBL2042175 molecule with primary amide side chain of Asn218 (Figure 4d).

Intermolecular interactions representative of cluster 5 (CHEMBL2042175)
The aromatic substituent of this flavonoid molecule carries aromatic interaction with the nicotinamide ring NAD501. Moreover, there is a T-shaped aromaticaromatic interaction between the aromatic substituent of CHEMBL1668048 and the aromatic ring of Tyr267 side chain. It is also possible to observe hydrophobic interaction with the amino acids Arg318, Trp131, Asn218 and Val222 (Table 1). The 3D map shows the interaction of π-donor type between rings A and C of flavonoid moiety of CHEMBL1668048 molecule with primary amide side chain of Asn218 (Figure 4e).

Analysis of the pattern of interactions that favor the stability of the ligand-receptor complex
The detailed analysis of intermolecular interactions of the five molecules representing each cluster and PfENR enzyme showed that there is an essential interactions pattern favoring the ligand-receptor complex stability and consequently result in better energy affinity values.
The analyzed molecules presented electrostatic interactions with cofactor NAD+, with the Tyr267 or Tyr277. These forms of electrostatic interactions of the molecules with NAD501 and PfENR have been described as a triplet of interactions between an aromatic group in PfENR inhibitors, the nicotinamide ring of NAD501 and residues Tyr267 and Tyr277 are apparently important for the affinity of the compound to the site active enzyme, in addition, it is believed that this type of interaction contributes significantly to the stability of the ternary complex PfENR-NAD501-compound inhibitor. 40 The 3D maps also allowed viewing π-donor interactions between the rings of the flavonoid moiety of the molecules and the primary amide side chain of Asn218 and in some cases it is also noted that interact with Ala219. Moreover, it was possible to detect hydrophobic interactions between Arg318, Trp131, Asn218, Val222, Val134, Met281 and Ala319 amino acids with flavonoid moiety, and in most of the molecule's interactions with Arg318, Trp131, Asn218 and Val222 residues remained. It is important to note that hydrophobic interactions with Asn218 and Val222 residues have been described in the literature. [30][31][32][33] Figure 5 shows the hydrophobic pocket that perform flavonoid moiety interactions with the core. This analysis shows the importance of flavonoid moiety for inhibiting PfENR, confirming the choice of this metabolite class to study and development of new antimalarial candidates. Furthermore, it is possible to infer the importance of a candidate molecule inhibitor of the enzyme aromatic PfENR in having regions that allow electrostatic interactions with the NAD501 and Tyr277 or Tyr267 or regions that allow the hydrogen interaction with these residues, especially Tyr277. These interactions are apparently important for the affinity of the compound to the active site of the enzyme and can assist in the stability of the complex. Studies 38 report that Tyr277 have been proposed as being important for the catalytic mechanism of the enzyme. It is likely that the interaction between the compounds and enzyme Tyr277 avoids this residue connection with its natural substrate. Thus, it can be proposed that a successful PfENR inhibitor needs to interact with the Tyr277, and thus occupy the space of interaction with the natural substrate in order to block its role in catalysis.
Aromatic regions, such as the flavonid nucleus, that perform hydrophobic interactions, in addition to interactions π-donor type with Asn218 are also important for a PfENR inhibiting compound. Thus, it is possible to determine that PfENR inhibition occurs through the sum of non-covalent interactions, such as electrostatic, hydrogen and hydrophobic. Therefore, as inhibitory molecules, they must present in their structure the ability to perform this type of interactions.
The specificity and affinity of the enzyme and its inhibitor will depend on the types of interactions between the molecules. 41,42 The analysis of intermolecular interactions made it possible to identify the presence of interactions described in the literature and are considered essential for the catalytic mechanism of the enzyme. All 30 selected molecules are potential inhibitors to PfENR based on the comparison of the affinity energy of these molecules with the inhibitors already described in the literature, as well as the standard identified intermolecular interactions.

Conclusions
The intermolecular interactions analyzed between 30 synthetic flavonoid molecules with a receptor showed that the affinity energy is related to the sum of hydrogen interactions, electrostatic and hydrophobic, a fact reinforced by the literature available on the inhibition of PfENR.
In addition, we observed the importance of flavonoid moiety to hydrophobic interactions with the active site of PfENR, a fact that confirms the choice of this metabolite class to study PfENR inhibitors.

Supplementary Information
Supplementary data (molecules and structures) are available free of charge at http://jbcs.sbq.org.br as PDF file.