Identification of Selective Novel Hits against Plasmodium falciparum Prolyl tRNA Synthetase Active Site and a Predicted Allosteric Site Using in silico Approaches

Abstract

Recently, there has been increased interest in aminoacyl tRNA synthetases (aaRSs) as potential malarial drug targets. These enzymes play a key role in protein translation by the addition of amino acids to their cognate tRNA. The aaRSs are present in all Plasmodium life cycle stages, and thus present an attractive malarial drug target. Prolyl tRNA synthetase is a class II aaRS that functions in charging tRNA with proline. Various inhibitors against Plasmodium falciparum ProRS (PfProRS) active site have been designed. However, none have gone through clinical trials as they have been found to be highly toxic to human cells. Recently, a possible allosteric site was reported in PfProRS with two possible allosteric modulators: glyburide and TCMDC-124506. In this study, we sought to identify novel selective inhibitors targeting PfProRS active site and possible novel allosteric modulators of this enzyme. To achieve this, virtual screening of South African natural compounds against PfProRS and the human homologue was carried out using AutoDock Vina. The modulation of protein motions by ligand binding was studied by molecular dynamics (MD) using the GROningen MAchine for Chemical Simulations (GROMACS) tool. To further analyse the protein global motions and energetic changes upon ligand binding, principal component analysis (PCA), and free energy landscape (FEL) calculations were performed. Further, to understand the effect of ligand binding on the protein communication, dynamic residue network (DRN) analysis of the MD trajectories was carried out using the MD-TASK tool. A total of ten potential natural hit compounds were identified with strong binding energy scores. Binding of ligands to the protein caused observable global and residue level changes. Dynamic residue network calculations showed increase in betweenness centrality (BC) metric of residues at the allosteric site implying these residues are important in protein communication. A loop region at the catalytic domain between residues 300 and 350 and the anticodon binding domain showed significant contributions to both PC1 and PC2. Large motions were observed at a loop in the Z-domain between residues 697 and 710 which was also in agreement with RMSF calculations that showed increase in flexibility of residues in this region. Residues in this loop region are implicated in ATP binding and thus a change in dynamics may affect ATP binding affinity. Free energy landscape (FEL) calculations showed that the holo protein (protein-ADN complex) and PfProRS-SANC184 complexes were stable, as shown by the low energy with very few intermediates and hardly distinguishable low energy barriers. In addition, FEL results agreed with backbone RMSD distribution plots where stable complexes showed a normal RMSD distribution while unstable complexes had multimodal RMSD distribution. The betweenness centrality metric showed a loss of functional importance of key ATP binding site residues upon allosteric ligand binding. The deep basins in average L observed at the allosteric region imply that there is high accessibility of residues at this region. To further analyse BC and average L metrics data, we calculated the ΔBC and ΔL values by taking each value in the holo protein BC or L matrix less the corresponding value in the ligand-bound complex BC or L matrix. Interestingly, in allosteric complexes, residues located in a loop region implicated in ATP binding had negative ΔL values while in orthosteric complexes these residues had positive ΔL values. An increase in contact frequency between residues Ser263, Thr267, Tyr285, and Leu707 at the allosteric site and residues Thr397, Pro398, Thr402, and Gln395 at the ATP binding TXE loop was observed. In summary, this study identified five potential orthosteric inhibitors and five allosteric modulators against PfProRS. Allosteric modulators changed ATP binding site dynamics, as shown by RMSF, PCA, and DRN calculations. Changes in dynamics of the ATP binding site and increased contact frequency between residues at the proposed allosteric site and the ATP binding site may explain how allosteric modulators distort the ATP binding site and thus might inhibit PfProRS. The scaffolds of the identified hits in the study can be used as a starting point for antimalarial inhibitor development with low human cytotoxicity.

Keywords: Aminoacyl tRNA synthetase; MD-TASK; allosteric modulators; dynamic residue network; free energy landscape; virtual screening.

Figure 1

2D representations of identified South African natural compound hits against PfProRS active and allosteric sites. (A) Identified orthosteric hit compounds against PfProRS. (B) Identified allosteric hit compounds against PfProRS allosteric site.

Figure 2

Binding energy (kcal/mol) heat maps and binding modes of the selected orthosteric and allosteric hit compounds on PfProRS and HsProRS. (A) Binding energies of orthosteric hit compounds. Selected orthosteric hits binding at the PfProRS active site are shown at the left side of the heat map. Binding energy increases from black to yellow. A black colour shows good binding energy. (B) Binding modes of identified orthosteric hit compounds on PfProRS. (C) The binding modes of the selected orthosteric hits in HsProRS. In HsProRS, selected orthosteric ligands do not bind the targeted active and allosteric sites. (D) Binding energies of allosteric hit compounds. Selected allosteric hits binding at the PfProRS allosteric site are shown at the left side of the heatmap. A black colour shows good binding energy. (E) Binding modes of identified allosteric hit compounds on PfProRS. (F) Binding poses of identified allosteric hits in HsProRS. Selected allosteric hits do not bind the targeted active and predicted allosteric site in HsProRS. Adenosine binding pose at the active site is shown in blue sticks. The active site is shown in red ellipse while the identified allosteric site is shown in green ellipse. Binding modes of identified orthosteric and allosteric hits on HsProRS are shown in purple and yellow ellipses, respectively.

Figure 3

Binding modes of selected (A) orthosteric and (B) allosteric hits. (A) 2D representation of (a) SANC152; (b) SANC235; (c) SANC236; (d) SANC244; (e) SANC318; (f) halofuginone (PDB ID: 4YDG) [37] and their binding modes with PfProRS active site. (B) 2D representation of (a) SANC184; (b) SANC257; (c) SANC264; (d) SANC456; (e) SANC622 (f) TCMDC-124506 (PDB ID: 4WI1); (g) glyburide (PDB ID: 5IFU) and their binding mode with PfProRS allosteric site.

Figure 4

PfProRS backbone RMSD analysis of all atom MD trajectories for the holo system and PfProRS-ligand complexes during the 200 ns simulation. (A) Holo system (PfProRS-ADN), (B) PfProRS-SANC152 complex, (C) PfProRS-SANC235 complex, (D) PfProRS-SANC236 complex, (E) PfProRS-SANC244 complex, (F) PfProRS-SANC318 complex, (G) PfProRS-SANC184 complex, (H) PfProRS-SANC257 complex, (I) PfProRS-SANC264 complex, (J) PfProRS-SANC456 complex and (K) PfProRS-SANC622 complex.

Figure 5

Backbone RMSD distribution plots for PfProRS ligand-bound complexes. Conformational flexibility can be assessed by comparing the backbone RMSD distribution of each ligand-bound complex and the holo protein. Comparison of the mean (µ) of the holo protein (black dashed line) to each ligand complex (green dashed line) demonstrates the shift in conformation distribution of the complexes during the 200 ns MD simulation. σ is the standard deviation. (A) Holo system (PfProRS-ADN), (B) PfProRS-SANC152 complex, (C) PfProRS-SANC235 complex, (D) PfProRS-SANC236 complex, (E) PfProRS-SANC244 complex, (F) PfProRS-SANC318 complex, (G) PfProRS-SANC184 complex, (H) PfProRS-SANC257 complex, (I) PfProRS-SANC264 complex, (J) PfProRS-SANC456 complex and (K) PfProRS-SANC622 complex.

Figure 6

3D representation of binding free energy landscape as a function of PC1 and PC2. Energy distribution is shown by the coloring pattern: Blue defines the conformational space with minimum energy (stable state) while red defines a conformational space with maximum energy (unstable state). Transient local energy states are defined by intermediate color patterns. PC1 and PC2 are displayed as a contour map at the bottom of each FEL plot with similar color pattern like the energy landscape. (A) Holo protein (PfProRS-ADN), (B) PfProRS-SANC152 complex, (C) PfProRS-SANC235 complex, (D) PfProRS-SANC236 complex, (E) PfProRS-SANC244 complex, (F) PfProRS-SANC318 complex, (G) PfProRS-SANC184 complex, (H) PfProRS-SANC257 complex, (I) PfProRS-SANC264 complex, (J) PfProRS-SANC456 complex and (K) PfProRS-SANC622 complex.

Figure 7

Bar graph representation of protein backbone root mean square deviation of FEL stable states compared to initial structures. The holo protein and protein-ligand complexes are shown on the X-axis.

Figure 8

Per residue root mean square fluctuation analysis of the holo system and PfProRS ligand-complexes during the 200 ns simulation. (A) PfProRS-SANC152 complex, (B) PfProRS-SANC235 complex, (C) PfProRS-SANC236 complex, (D) PfProRS-SANC244 complex, (E) PfProRS-SANC318 complex, (F) PfProRS-SANC184 complex, (G) PfProRS-SANC257 complex, (H) PfProRS-SANC264 complex, (I) PfProRS-SANC456 complex and (J) PfProRS-SANC622 complex.

Figure 9

Contribution of each residue to PC1 and PC2. PC1 is shown in black while PC2 is shown in red for all the ligand-bound complexes. (A) Holo protein, (B) PfProRS-SANC152 complex, (C) PfProRS-SANC235 complex, (D) PfProRS-SANC236 complex, (E) PfProRS-SANC244 complex, (F) PfProRS-SANC318 complex, (G) PfProRS-SANC184 complex, (H) PfProRS-SANC257 complex, (I) PfProRS-SANC264 complex, (J) PfProRS-SANC456 complex and (K) PfProRS-SANC622 complex.

Figure 10

Cartoon representation of PfProRS homology model showing structural mapping of residues with significant ΔBC (2 standard deviations). ΔBC values were calculated by taking the BC value of the holo protein less the ligand-bound complex. Structures were obtained at 200 ns of the MD simulations. Residues with negative ΔBC values are shown in blue while residues with positive ΔBC values are shown in red. The active site is shown in green ellipses, ATP binding site in yellow ellipses and the predicted allosteric site in orange ellipses. (A) PfProRS-SANC152 complex, (B) PfProRS-SANC235 complex, (C) PfProRS-SANC236 complex, (D) PfProRS-SANC244 complex, (E) PfProRS-SANC318 complex (F) PfProRS-SANC184 complex, (G) PfProRS-SANC257 complex, (H) PfProRS-SANC264 complex, (I) PfProRS-SANC456 complex and (J) PfProRS-SANC622 complex.

Figure 11

Cartoon representation of PfProRS homology model showing structural mapping of residues with significant ΔL (2 standard deviations). ΔL values were calculated by taking the L value of the holo protein less the ligand-bound complex. Structures were obtained at 200 ns of the MD simulations. Residues with negative ΔL values are shown in blue while residues with positive ΔL values are shown in red. The active site is shown in green ellipses, ATP binding site in yellow ellipses and the predicted allosteric site in orange ellipses. (A) PfProRS-SANC152 complex, (B) PfProRS-SANC235 complex, (C) PfProRS-SANC236 complex, (D) PfProRS-SANC244 complex, (E) PfProRS-SANC318 complex (F) PfProRS-SANC184 complex, (G) PfProRS-SANC257 complex, (H), PfProRS-SANC264 complex, (I) PfProRS-SANC456 complex and (J) PfProRS-SANC622 complex.

Figure 12

Frequency of contacts at the determined allosteric site in the holo (left sub-figure) and the ligand-bound complexes (right sub-figure). The allosteric site residue is shown in the middle of each contact map while residues at the ATP binding TXE loop are circled in red. (A) Contact frequency between residue Ser263 and Thr397 (B) Contact frequency between residue Thr267 and Pro398, (C) Contact frequency between residue Tyr285 and Thr402 and (D) Contact frequency between residue Leu707 and Gln395.

Figure 13

Frequency of contacts at the determined allosteric site in the holo (left sub-figure) and the ligand-bound complexes (right sub-figure). The allosteric site residue is shown in the middle of each contact map while residues at the Z-domain are circled in red. (A) Contact frequency between Gln395, Thr706 and Leu707, (B) Contact frequency between Arg472 and Ser708, (C) Contact frequency between Arg744 and Gly709 and (D) Contact frequency between Ser745 and Ser708.