BioMed Research International
Volume 2017 (2017), Article ID 5270940, 15 pages
https://doi.org/10.1155/2017/5270940

Research Article

Computational Exploration for Lead Compounds That Can Reverse the Nuclear Morphology in Progeria

Shailima Rampogu,¹ Ayoung Baek,¹ Minky Son,¹ Amir Zeb,¹ Chanin Park,¹ Raj Kumar,¹ Gihwan Lee,¹ Donghwan Kim,¹ Yeonuk Choi,¹ Yeongrae Cho,¹ Yohan Park,² Seok Ju Park,³ and Keun Woo Lee¹

¹Division of Applied Life Science (BK21 Plus), Plant Molecular Biology and Biotechnology Research Center (PMBBRC), Systems and Synthetic Agrobiotech Center (SSAC), Research Institute of Natural Science (RINS), Gyeongsang National University (GNU), 501 Jinju-daero, Jinju 52828, Republic of Korea
²College of Pharmacy, Inje University, 197 Inje-ro, Gimhae, Gyeongnam 50834, Republic of Korea
³Department of Internal Medicine, College of Medicine, Busan Paik Hospital, Inje University, Gyeongnam, Republic of Korea

Correspondence should be addressed to Seok Ju Park; psjgenesis@hanmail.net and Keun Woo Lee; kwlee@gnu.ac.kr

Received 25 July 2017; Accepted 24 September 2017; Published 26 October 2017

Academic Editor: Rituraj Purohit

Copyright © 2017 Shailima Rampogu et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Progeria is a rare genetic disorder characterized by premature aging that eventually leads to death and is noticed globally. Despite alarming conditions, this disease lacks effective medications; however, the farnesyltransferase inhibitors (FTIs) are a hope in the dark. Therefore, the objective of the present article is to identify new compounds from the databases employing pharmacophore based virtual screening. Utilizing nine training set compounds along with lonafarnib, a common feature pharmacophore was constructed consisting of four features. The validated Hypo1 was subsequently allowed to screen Maybridge, Chembridge, and Asinex databases to retrieve the novel lead candidates, which were then subjected to Lipinski’s rule of 5 and ADMET for drug-like assessment. The obtained 3,372 compounds were forwarded to docking simulations and were manually examined for the key interactions with the crucial residues. Two compounds that have demonstrated a higher dock score than the reference compounds and showed interactions with the crucial residues were subjected to MD simulations and binding free energy calculations to assess the stability of docked conformation and to investigate the binding interactions in detail. Furthermore, this study suggests that the Hits may be more effective against progeria and further the DFT studies were executed to understand their orbital energies.

1. Introduction

Hutchinson-Gilford Progeria Syndrome, Progeria, is a rare genetic disorder seen in children and is manifested by premature aging [1]. This fatal disorder was studies by two scientists Jonathan Hutchison in the year 1886 and Hastings Gilford [2] in 1897 and hence it was named Hutchinson-Gilford Progeria Syndrome (HGPS) [3, 4]. Progeria is originated from the Greek terminology “progeros” that refers to “prematurely old” [5, 6] and affects 1 in 4–8 million [7, 8]. This syndrome can be observed in both sexes with 2 : 1 male/female ratio and across different countries showing no geographic and ethnic bias [8, 9]. The general life expectancy is approximately an average of 13 years and the affected may die due to several reasons [10]. Conversely, only a single case of a patient who lived 45 years of age exists [11]. The affected demonstrates a characteristic features by displaying delayed growth, osteoporosis, cardiovascular ailments, alopecia, pinched nose, and sclerodermatous skin [12, 13]. However, they display no change in the mental ability [4]. This is because the brain largely synthesizes lamin C and very little prelamin A [14].

Progeria is defined as sporadic autosomal dominant mutation [15], whose progression begins in utero [16]. Though the circumoral pallor symptom was found associated with the child at the time of birth, the delay in the phenomenal representation of the disease is due to the low levels of progerin during the undifferentiated embryonic cells and is demonstrated after the levels are elevated [17, 18].

The nuclear morphology is imperative in demonstrating progeria [19]. The abnormality of the nucleus is due to the mutations of two genes lamin (LMNA) and ZMPSTE 24 [20]. More specifically, the point mutations that occur in the lamin A/C are vital in demonstrating the disease [19, 21]. Dominant negative form of lamin A protein is produced as a consequence of the mutation G608G (GGC to GGT) within the exon 11 of lamin [22–24]. This mutation results in the formation of cryptic splice site subsequently causing the cleavage of 50 amino acid residues in the C-terminus of lamin A [19, 25, 26] and thus forms a protein named progerin with a distorted nucleus. This results in the deletion of site ZMPSTE24, demonstrating a permanent farnesylated protein leading to abnormal nucleus. Such cells with abnormal nucleus are prone to develop several diseases which are collectively referred to as laminopathies [27, 28] such as Emery-Dreifuss muscular dystrophy [29, 30], Dunnigan-type familial partial lipodystrophy [31, 32], and mandibuloacral dysplasia [33, 34]. It is therefore evident that the defective lamin A influences the instability of the protein thereby developing the diseases. Additionally, the accumulated progerin renders abnormalities in the behaviour of chromosome segregation and the reassembly of the nuclear envelop [16, 35]. Additionally, it dislocates the centromere protein-F (CENP-F) from kinetochores [16]. Subsequently, the genetic instability is elevated, thus favouring premature aging.

Despite alarming condition, promising therapeutic treatments are still under trails. Under such circumstances, drugs that were originally developed for certain diseases have proven to be effective against progeria [4]. Pravastatin, originally developed against cardiovascular diseases [36–38], zoledronic acid, a bisphosphate employed for treating osteoporosis [36–38], and farnesyltransferase inhibitors (FTIs) [36–38], used to treat cancers, have improved the condition of the progeroid children. Among them, the FTIs have ameliorated the diagnostic conditions in the affected by reverting the abnormalities of the nucleus [39–42]. Additionally, they have effectively improved the nuclear blebbing in the mouse models [43–46]. Mechanistically, FTIs operate by inhibiting the conversion of prelamin A to mature lamin A [13, 16, 47–49] and further improve the cardiovascular and skeletal pathologies besides gaining weight [4, 49]. Lonafarnib, one of the FTIs that has gained increasing popularity for treating progeria, has reached the clinical trials [16, 50, 51]. Owing to the beneficial effects of FTIs, it is mightily essential to identify new drugs which can perform with similar strategy. Therefore, in the current study, we focused on screening new chemical compounds that might be able to treat progeria using ligand-based pharmacophore method. In order to sieve the potential candidate compounds, the chemical features of lonafarnib have also been considered.

2. Materials and Methods

2.1. Generation of the Pharmacophore Model

For the generation of the pharmacophore nine inhibitors, Figure 1, from various literatures have been considered along with lonafarnib with the known inhibitory activities [52, 53]. These inhibitors have displayed different IC₅₀ values and diverse structure. The common feature pharmacophore generation protocol available on the Discovery Studio (DS) v4.5 (Accelrys, San Diego, CA) was used maintaining a minimum interfeature distance of 2.00 with fast flexible conformation generation. Common feature pharmacophore is generated using the HipHop algorithm that determines the common chemical features associated with the 3D spatial arrangements in a given training set. Additionally, the identification of the configurations is an exhaustive process that begins with a small set of features and proceeds on larger configurations. Subsequently, the pharmacophores are ranked as they are built, a method employed to explain the ability to map onto the pharmacophores and further their rarity. From the generated pharmacophores the best model was selected based upon the rank scores [54].

Figure 1: Training set compounds employed to build the pharmacophore. The chemical features of lonafarnib are also exploited in its construction. IC₅₀ values are indicated in parentheses.

2.2. Validation of the Pharmacophore Model

The generated pharmacophore was validated adapting the Receiver Operating Characteristic (ROC) and the decoy set methods to assess its robustness in identifying the lead molecules. The ROC was initiated simultaneously during the generation of the pharmacophore. The nine inhibitors that were used to generate the pharmacophore model were taken as known actives and six other compounds were recruited as known inactives. On the other hand, the decoy set method was evaluated constituting a database of 107 compounds comprised of 9 active molecules. The quality of the pharmacophore was evaluated by computing the goodness of fit (GF) [55] score employing the formulae and .

2.3. Database Screening for Retrieving the Virtual Candidates and Assessment of Drug-Like Properties

To identify the candidate drug compounds with high or similar potentiality, the validated pharmacophore has been used as a 3D query to search the chemical databases. Three databases have been used, namely, the Chembridge, Maybridge, and Asinex consisting of 50,000, 59,652, and 21,3262 compounds, respectively. Employing the Ligand Pharmacophore Mapping protocol implemented on the DS, the compounds were retrieved and were further culled setting the fit value as 3. Subsequently, obtained compounds were filtered on the bases of Lipinski’s rule of five [56] and ADMET properties [57]. Lipinski’s rule of five logically indicates that the well-absorbed compounds exhibits a of less than 5, less than 5 hydrogen bond donors, and less than 500 molecular weight, respectively. Additionally, they also have the hydrogen bond accepting capability of less than 10. Furthermore, the compounds were checked by the ADMET descriptors to evaluate if they can cross the bold-brain barrier (BBB), possess low toxicity, and have good solubility and human intestinal absorption. Among them the key descriptors are the BBB and displaying no hepatotoxicity. The absorption, solubility, and the BBB were fixed at 3, 3, and 0, respectively [58, 59]. The resultant compounds were escalated to molecular docking along with the 9 compounds used for the generation of the pharmacophore.

2.4. Molecular Docking Mechanism

Molecular docking is one of the superior methods employed as a sampling method to identify the most accurate conformation [60]. Furthermore this technique also identifies the compounds that can fit into the active site of the target molecule and their corresponding interactions with the residues. For the current study, CDOCKER, implemented on the DS, was employed to understand the binding affinities of the protein and the ligand. The results were evaluated based upon the CDOCKER energy, while the CDOCKER interaction energy is used as a rescore. Highest CDOCKER interaction energy implies greater favourable binding [61, 62]. The protein target for the present research is the farnesyltransferase of high resolution (1.8 Å) with the PBD code 1TN6 bound with the cocrystal imported form the protein data bank [63]. The heteroatoms were removed performing the clean protein available in the DS. The protein was further minimized employing the CHARMm force field. Furthermore the binding site was evaluated 15 Å around the cocrystal and the histidine tautomers were protonated to ND1H state as was observed in the crystal structure. The key residues located at the active site have been identified as His748, Arg791, Lys794, and Tyr800 [63]. Each ligand was allowed to generate 100 conformers and depending upon the scoring functions and molecular interactions between the protein and the ligand, the ideal dock pose was chosen.

2.5. Molecular Dynamics Simulation

To further affirm the potentiality of the selected compounds and to evaluate the dynamic behaviour of the prospective drug molecules in the binding site pocket, they were subjected to MD simulations along with the reference compounds using GROMACS 4.5.7, employing CHARMm27 force field [64]. MD simulations were executed to examine the binding stability of the identified lead compounds in comparison with the reference compound. Ligand topologies were generated using SwissParam [65]. Corresponding counter ions were added to neutralize the solvated TIP3P water model present in the dodecahedron box. Unwanted contacts from the initial structure were dislodged by performing the energy minimization, adapting the steepest descent algorithm which was followed by the NVT and NPT equilibration steps. During this process, the solvent molecules along with the counter ions were allowed to move restraining the protein backbone. Both the processes were executed for 100 ps at 300 K and a pressure of 1 bar, respectively. Parrinello-Rahman barostat was employed to maintain the pressure of the system [66]. The geometry of the water molecules and the bonds involving hydrogen atoms were constrained employing the SETTLE [67] and LINCS [68], respectively. Particle Mesh Ewald (PME) [69] was used to calculate long-range electrostatic interactions. A cut-off distance of 12 Å was attributed for Coulombic and van der Waals interactions. The equilibrated structures were then subjected to production run which was conducted for 20 ns and the results were evaluated using VMD [70] and DS.

2.6. Binding Free Energy Calculations

To delineate further on the protein-ligand complex, time-dependent binding free energy calculations were performed. Molecular Mechanics/Poisson-Boltzmann Surface Area (MM/PBSA) [71, 72] was used for its accomplishment and was performed after the MD simulations. The obtained should take into account the protein fluctuations and the ligand conformations, which therefore ensures proper positioning of the ligand within the binding pocket.

The binding free energy protein-ligand complex in solvent system is stated as

Herein, refers to the total free energy of the complex and and indicate the separated protein and ligand in the solvent. Their free energies can be computed by

where can be a protein, ligand, or its complex. represents the average molecular mechanics potential energy in vacuum, while interprets the free energy present in the solvation.

Additionally, molecular mechanics potential energy in vacuum can be evaluated by adapting the equation

represents the bonded interactions, while the nonbonded interactions are denoted by . is generally regarded as zero [73].

The solvation free energy () is expressed by the sum of electrostatic solvation free energy () and a polar solvation free energy () and is given as follows:

is computed recruiting the Poisson-Boltzmann (PB) equation [74] while is computed from the solvent-accessible surface area (SASA) and can be written as follows:

Here, γ is the coefficient of the surface tension of the solvent, whereas is its fitting parameter, whose values are 0.02267 kJ/mol/Å² or 0.0054 kcal/mol/Å² and 3.849 kJ/mol or 0.916 kcal/mol, respectively.

2.7. Density Functional Theory

The MD optimized conformations were further examined by density functional theory (DFT) analysis implemented with DS. With DFT we can calculate the Kohn-Sham orbital energies by means of the Highest Occupies Molecular Orbital and (HUMO) and Lowest Unoccupied Molecular Orbital (LUMO) [75, 76]. Moreover, the ionization potential (electron donor capacity) was demonstrated by the HOMO, while the electron affinity (electron acceptor) was represented by LUMO. In order to understand the energy transfer and further the stability of the small molecules within the binding site [58, 77], the MD optimized docked conformations were subjected to DFT. A lower energy gap between the molecules demonstrates that the molecules are highly reactive while the higher energy gap implies low reactivity [78–80]. The Dmol3 and Becke, three-parameter, Lee-Yang-Parr (B3LYP) [81, 82], using the DND basis set with self-consistent field (SCF) density convergence of , available on the DS, was employed for computing the DFT. Additionally, the DFT studies were carried out to evaluate the electronic properties of the obtained Hits and the reference compound.

3. Results and Discussion

3.1. Generation of the Pharmacophore Model

Common feature pharmacophore module that utilizes the HipHop [83] algorithm implemented on the DS was employed to generate the pharmacophore with minimum features and maximum features being 1 and 10, respectively. Delineating on the pharmacophore features, it was observed that the HHDA feature was present in all the pharmacophores; however it was absent in one model. Furthermore, the ring aromatic (RA) feature was noticed in six models. It could therefore be understood that these features are the key features to be possessed by the drug molecules. Accordingly, care was taken in selecting an ideal pharmacophore that contains these features and hence Hypo 1 was selected as the best model displaying a rank score of 26.307, Table 1. Subsequently, a four-featured pharmacophore model was generated comprised of two-hydrophobic, one-ring aromatic, and one-hydrogen bond acceptor, Figure 2.

Table 1: Common feature pharmacophore results as generated by the HipHop.

Hypo number
Features
Rank
Direct Hit
Partial Hit
Max fit

1
RHHDA
26.307
1111
0000
4

2
RHHDA
25.624
1111
0000
4

3
ZHHDA
25.574
1111
0000
4

4
RHHDA
25.352
1111
0000
4

5
RHHDA
25.352
1111
0000
4

6
HHHDA
24.569
1111
0000
4

7
RHHDA
24.450
1111
0000
4

8
HHHDA
24.362
1111
0000
4

9
ZHHDA
24.290
1111
0000
4

10
RHAA
24.069
1111
0000
4

R: ring aromatic, H: hydrophobic, HD: hydrogen bond donor, A: hydrogen bond acceptor, Z: zinc binder.

Figure 2: The best pharmacophore model (Hypo 1) consisting of four pharmacophoric features with its geometry.

3.2. Validation of the Pharmacophore Model

The generated pharmacophore was validated utilizing two different methods, the ROC curve and the decoy set method. The ROC validation was performed during the pharmacophore generation. In this, the nine training set ligands are labeled as known actives and other nonspecific ligands are referred to as known inactives. Ideally, a good pharmacophore should recognize the known actives form the inactives and thus evaluates the quality of the pharmacophore. The generated pharmacophore was successful in identifying the known actives and displayed an excellent quality, Figure 3.

Figure 3: ROC curve. The quality of the pharmacophore was observed to be good.

In the decoy set method, an external database of 107 compounds has been instituted consisting of nine actives (). Ligand pharmacophore mapping available on DS was then launched. Pharmacophore has mapped with nine compounds (Ht) in which eight compounds were the active compounds (Ha) conferring 88.8% yield of actives. Furthermore, the (Goodness of Hit Score) GH and the (Enrichment Factor) EF scores have been computed to be 0.71 and 9.4, respectively. Generally, the GH score lies between 0 and 1, where 0 represents null model while 1 indicates the best model. Furthermore, good model should have a GH score above 0.7 [84, 85]. Accordingly, the generated Hypo1 is considered good as it demonstrated a score of 0.71. Different calculations computed as a part of decoy set validation are tabulated in Table 2.

Table 2: Different parameters computed by decoy set method. The GF score confirms the predictive ability of the pharmacophore.

S. number
Parameters
Values

Total number of molecules in database ()
107

Total number of actives in database ()
9

Total number of Hit molecules from the database (Ht)
9

Total number of active molecules in Hit list (Ha)
8

% yield of active (Ha/Ht)
88.8

% ratio of actives [(Ha/) × 100]
88

Enrichment factor (EF)
9.4

False negatives
1

False positives
1

Goodness of fit score (GF)
0.71

GH = and .

3.3. Database Screening for Retrieving the Virtual Candidates and Assessment of Drug-Like Properties

Virtual screening is a process of determining every small molecule provided in the databases so as to judge their capability of binding with the target protein. Pharmacophore based virtual screening largely depends upon the chemical features present on the pharmacophore in critically assessing and identifying the candidate compounds from the databases. In the current study, the pharmacophore has mapped with 24037, 27513, and 41385 compounds of Chembridge, Maybridge, and Asinex databases, respectively. Subsequently, compounds were filtered based upon the fit value greater than five. As a result, the obtained 5682, 6955, and 10296 compounds were examined for their drug-like properties using the Lipinski’s rule of 5 and ADMET. Consequently, a total of 3372 compounds were obtained, Figure 4, and were forwarded to molecular docking along with the training set compounds.

Figure 4: Different steps involved in screening the candidate compounds. 3,372 compounds are obtained from screening and 9 compounds belong to the training set.

3.4. Molecular Docking Mechanism

The compounds obtained from the previous steps along with the nine compounds were subjected to molecular docking. Hereinafter, the compound with the lowest IC₅₀ value from the nine training set compounds is referred to as the reference compound. In order to assess the suitability and accuracy of the CDOCKER prior to the commencement of the protocol, the cocrystal ligand was docked into the crystal structure of the target. The resultant pose has generated an acceptable RMSD of 1.4 Å, Figure 5. For assessing the dock results, the reference dock scores and the scores of lonafarnib were considered. The reference has generated an interaction energy of 23.5208 kcal/mol, while the lonafarnib has displayed 50.6141 kcal/mol, Table 3. In pursuit of identifying the most biologically active compounds, the compounds that have demonstrated higher CDOCKER interaction energy and higher CDOCKER energy than the reference and lonafarnib have been considered for further studies. Consequently, a total of six compounds have been retrieved from Asinex, Chembridge, and Maybridge databases. These compounds were further probed for their interactions with the residues that reside at the active site groove. Among the six Hit compounds listed in Table 3, only two compounds (CHEM, AXN_4) were observed to show interactions with the key residues and mapped with all the features exhibited by Hypo 1, Figure 6, and therefore these compounds were escalated to the MD simulations.

Table 3: Docking results according to CDOCKER interaction energy of the potential candidates. Dock scores higher than reference and lonafarnib were considered.

S. number
Name of the compound
-CDOCKER energy
(kcal/mol)
-CDOCKER interaction energy (kcal/mol)

AXN_1
23.97
56.60

AXN_2
20.47
51.07

AXN_3
20.51
51.67

AXN_4
24.13
52.01

CHEM
28.53
51.33

MAY
40.37
52.38

Lonafarnib
20.16
50.61

Reference
22.57
23.52

Figure 5: Validating the docking parameters using the cocrystal. Pink refers to the cocrystal and cyan refers to the docked pose.

Figure 6: Mapping of the Hits onto the pharmacophore. The Hits are found to map with all the features of the pharmacophore.

3.5. Molecular Dynamics Simulations

To further investigate the structural stability and to evaluate the dynamic behaviour of the Hit compounds in the binding site of farnesyltransferase, the MD simulation was executed for the reference, lonafarnib, and the Hits. Accordingly, 20 ns run was initiated to understand the conformational variations and the nature of the ligand within the active site. For this assessment, the docked poses were determined as the initial conformations. Subsequently, the root mean square deviation (RMSD) of the protein backbone atoms and the radius of gyration (Rg) were evaluated. The RMSD of all the complexes was found to be within 0.27 nm; however, Hit 1 has displayed slightly higher RMSD of 0.25 nm, while the others were stable at 0.2 nm, Figure 7. Furthermore, it was noticed that towards the last 3 ns the systems were converged. Additionally, the radius of gyration that implies the compactness of the proteins, revealed that the four systems are finely folded with no major aberrations after 8000 ps and are represented between 2.13 and 2.15 nm, Figure 8.

Figure 7: Structural stability analysis through RMSD of the four complexes.

Figure 8: Radius of gyrations of the four systems to understand their compactness.

To further analyse the binding modes of each compound, the representative structures of the last 3 nm were obtained as it was observed to be well converged, Figure 7. Upon superimposition of the corresponding structures, it was established that the binding fashion of the Hits was in agreement with lonafarnib and the reference compound, Figure 9. The Hits were thereafter analysed for the intermolecular interactions with the crucial residues. Reference compound has formed two hydrogen bonds with Arg 791 and Lys794 with a distance of 2.2 Å and 1.9 Å, respectively. Lonafarnib has demonstrated two hydrogen bonds through Cys754 and Tyr800 represented by a length of 2.2 Å and 2.1 Å. On the contrary, the Hits have demonstrated three hydrogen bonds each. Hit 1 has produced two hydrogen bonds with Arg791 and Lys794, represented by a bond length of 2.4 Å, 2.0 Å, and 2.0 Å, respectively. Similarly, Hit 2 also has displayed three hydrogen bonds, one with Lys794 and two with Arg791, portraying a length of 2.6, 1.8, and 1.9, correspondingly, Figure 10. Additionally, the benzene ring of the reference has formed the π interaction while Tyr751 and Trp803 have been involved in π-π T-shaped interactions with the other benzene ring. Lonafarnib has formed a π-anion bond with Asp797. Additionally, the Br group of the ligand has interacted with Alkyl with Leu795. Furthermore, the benzene ring has formed π-alkyl bond with Arg702 and Trp803, respectively. The two extreme benzene rings Hit 1 have displayed π-π stacked bonds with Tyr800 and Trp803 residues. Additionally, the Arg702, Tyr751, Cys799, and His862 have formed the alkyl and π-alkyl bonds. Hit 2 on the contrary has projected π-π stacked bonds with Tyr800. Upon closer look at the interactions, it is evident that Hit 1 has showed greater number of bonds which may imply its higher efficiency. The details of the interactions are tabulated in Table 4 and Supplementary 1 in Supplementary Material available online at https://doi.org/10.1155/2017/5270940. Moreover, to delineate on the hydrogen bonds and to gain insight and understand the nature of the ligand in the active site, the hydrogen bond interactions were monitored throughout the MD run. The Hits have demonstrated greater hydrogen bonds as compared to the reference and lonafarnib, displaying average hydrogen bonds of 0.4 and 1.6, respectively. The reference has shown 0.07 and lonafarnib has projected 0.2 hydrogen bonds at an average. Furthermore, the reference has rendered marginal hydrogen bonds between 12000 to 15000 ps. On the contrary, the Hits have demonstrated regular hydrogen bonds during the MD run, Figure 11. Additionally, their 2D structures are depicted in Figure 12, Hit 1 was retrieved from Chembridge database, and Hit 2 was obtained from Asinex database. Furthermore, probing the PubChem online search tool, it was affirmed that these compounds have not been explored against progeria.

Table 4: Different interactions rendered between the protein and the ligands.

S. number
Name of the compound
Hydrogen bonds < 3 Å
π-bonds
Van der Waals interactions

Residue
atom
Ligand atom
Bond length

Reference
Lys794
Arg791
HZ2
HH21
O21
O21
2.2
1.9
His748
Tyr751
Trp803
Leu795, Asp797,
Lys853, Gly790,
Tyr861, Gly750,
Arg702, Tyr800

Lonafarnib
Cys754
Tyr800
HG
HH
O5
N8
2.2
2.1
Arg702, Asp797, Leu795,
Trp 803
Trp602, Tyr654,
Tyr705, Cys706,
His748, Gly750,
Phe753, Arg791,
Lys794, Val796,
Asp852, Lys853

Hit 1
Arg791
Arg791
Lys794
HE
HH21
HZ2
O24
O24
O24
2.4
2.0
2.0
Arg702, Tyr751,
Cys754,
Tyr800
His862
Trp602, Ala651,
Tyr705, Cys706,
Tyr 705, His748,
Gly790, Leu795,
Asp797, Cys799,
Asp852, Lys856,
Tyr861

Hit 2
Lys794
Arg791
Arg791
HZ2
HH21
HE
O10
O10
O10
2.6
1.8
1.9
Tyr800
Trp606, Tyr751,
Gly750, Cys754,
Leu795, Val796,
Asp797, Trp803,
Asp852, Lys853,
Trp803, Tyr861,
Tyr865

Figure 9: Binding mode analysis of the reference the Hits. Picture on the left represents the superimposed structure and on the right is its enlarged form.

Figure 10: Hydrogen bond interactions between the protein and the ligands. Purple dashed lines denote the hydrogen bonds.

Figure 11: The number of intermolecular hydrogen bonds between protein and the ligands during 20 ns MD simulations.

Figure 12: 2D structures of the Hit compounds.

3.6. Binding Free Energy Calculations

MM/PBSA method was employed to compute the binding free energies for a given set of protein-ligand complex. In order to carry out this study, 20 snapshots [86] were evenly extracted for the four systems. These systems have displayed a between −20 kJ/mol~−110 kJ/mol; however, slight variations were noticed because the conformational space was not sampled enough. Additionally, the reference and lonafarnib have demonstrated a −32.78 kJ/mol and −60.39 kJ/mol, respectively. The Hits on the other hand were conferred with −64.38 kJ/mol (Hit 1) and −65.74 kJ/mol (Hit 2) demonstrating much lower binding energies than the reference and the lonafarnib, Figure 13. The generated binding energies are the sum of ligand conformations and the protein fluctuations so as to affirm proper positioning of the ligand in the active site cavity. Accordingly, the ligands (Hits) are seated in the charged active site of the farnesyltransferase through the hydrogen bond and the van der Waals interactions. Additionally, various hydrophobic interactions such as π-π, π-π T-stalked and π-alkyl bonds hold the ligand firmly within the active site. Furthermore, the Hits have displayed higher CDOCKER interaction energy and lower binding free energies, thus making themselves valuable in treating progeria. Delineating on the perresidue energy contribution, it was evident that in Hit 1 and Hit 2 Lys794 and Tyr800 have contributed majorly to the total energy as was seen with lonafarnib. Furthermore, His748 in the reference compound was a major contributor; however, in the other ligands the same impact was not displayed. Across all the ligands, Tyr800 was the largest contributor of the respective energy terms. Therefore, it can be deduced that Try800 might be imperative in inducing the inhibitory activity, Figure 14.

Figure 13: Estimation of binding free energy using MM/PBSA approach conducted during 20 ns.

Figure 14: Perresidue decomposition of the binding energies of the key residues.

3.7. Density Functional Theory

The molecular orbital energies were calculated in order to assess the HOMO and LUMO that are responsible for the charge transfers in a given chemical reaction and further describes a molecule to be encountered by the electrophiles and nucleophiles, respectively. Additionally, the band gap generated between the HOMO and LUMO demonstrates the reactivity of the molecules corresponding to smaller gap being more reactive and wider gap implies less reactive, Table 5, and therefore, the molecules with smaller band gap are considered. For the current study, the MD optimized lead conformations of the ligands were forwarded to the HOMO and LUMO analysis. The Hit molecules were chosen based upon the least energy gap as compared to the reference compound. Additionally, the electrostatic potential maps were computed to probe into the structural aspects of a molecule as it significantly plays a key role in the receptor-ligand interactions and are computed for a set of points in the molecule [87]. The electrostatic potential maps of the corresponding compounds have been obtained by the superimposition of the electrostatic potentials upon the electron density surfaces of the compound [88]. Furthermore, the electron density is plotted by the intensity of the colour that reflects the characteristic feature of a molecule. Subsequently, the red colour refers to the high negatively charged region and thus corresponds to high charge accumulation and the blue represents the charge depletion [89] and is positively charged region. The intermediate colours, such as orange, yellow, and green demonstrate the charges midway between both the extremes [90]. The order of the colour magnitude can be demonstrated as

 (Highly negative) red < orange < yellow < green < blue (highly positive).

Table 5: HOMO, LUMO, and band gap of the Hits and the reference compounds computed employing the DFT approach.

Name
HOMO
(eV)
LUMO
(eV)
Band gap
(eV)

Reference
−0.163
−0.0547
0.108

Hit 1
−0.154
−0.095
0.058

Hit 2
−0.179
−0.1149
0.064

Lonafarnib
−0.182
−0.0766
0.105

Furthermore, the electrostatic potentials have been plotted to evaluate the electrostatic interactions that exists between the protein and the ligand. In the current study, the electron density isovalue was taken as 0.03 with a default colour scheme that ranged between 0.05 and 0.1. The positive phase of the molecular orbitals utilizes an isovalue of 0.01 and is demonstrated in blue, while the negative phase uses an isovalue of −0.01 and is depicted in red. The half-way potentials between them are represented by the other aforementioned colours. Furthermore, the molecular features that are resulted from the self-consistent field (SCF) corresponds to the electrostatic potential atomic-centered charges depicted by the molecular electrostatic potential (MEP) Figure 15 [91–93]. Additionally, MEP was calculated by a host of points determined on the 3D structures. Additionally, the electrostatic potential charges can be used to compute the intermolecular properties rather than the intramolecular properties. However, the buried atoms that are not on the outer region of the molecular van der Waals (vdW) surface are not identified as MEP points and are read by the lower relative root mean square (RRMS) values [94, 95]. The lower the RRMS, the higher the accuracy of the fit demonstrated by the MEP that in turn is calculated by QM produced by the fitted charges for the individual atoms [94, 95], Supplementary 2. As compared to the reference compound, the Hits have generated lower RRMS fit values and are relatively equal to the lonafarnib, Supplementary 2. Furthermore, the Mulliken atomic charges [96] were computed with spin 0 for the MD optimized ligands, specifically to those atoms that are involved in the hydrogen bonds and are tabulated in Table 6. These charges are related to the overlap matrix of the atomic orbitals [97]. Mulliken atomic charges of oxygen atoms have ranged between −0.44 and −0.52 au, while the nitrogen atom has displayed relatively lower charge of −0.303 au.

Table 6: Mulliken atomic charges of the atoms involved in hydrogen bonds.

Compound name
Ligand atoms
Mulliken atomic charges (au)

Reference
O2
−0.449

Lonafarnib
O5
−0.527

N8
−0.303

Hit 1
O24
−0.448

Hit 2
O10
−0.437

Figure 15: DFT studies for the MD optimized conformations. The HOMO, LUMO, and the molecular electrostatic profiles of the four systems.

These results together with the HOMO and LUMO significantly portray the high electronegativity of the atoms. Furthermore, it can be observed that the Hits have demonstrated higher CDOCKER interaction energy and lower binding free energies and lower band gap as compared to the reference and lonafarnib. It can therefore be stated that the identified Hits have a similar efficacy or better electronic properties than the reference compounds in treating progeria.

4. Conclusion

Progeria is one of the rare genetic disorders manifested by premature aging in children leading to death. However, the currently available drugs are limited to clinical trials and therefore there is an essentiality for the discovery of new lead compounds. In the present article, we have successfully evaluated the potentiality of the novel lead candidates from the pharmacophore modelling coupled with molecular modelling studies. Furthermore, the putative binding modes of Hits have been examined upon comparison with the reference compound and the lonafarnib. Overall the results obtained from molecular docking, MD simulations, binding free energy, and DFT calculations indicate that the Hits have shown similar binding patterns with the reference and lonafarnib and have rendered strong molecular interactions with the active site residues of farnesyltransferase, thus qualifying as substantial compounds of interest in treating the Hutchinson-Gilford Progeria Syndrome.

Conflicts of Interest

The authors declare that there are no conflicts of interest.

Acknowledgments

This research was supported by Next-Generation BioGreen 21 Program from Rural Development Administration, Republic of Korea (Grant no. PJ01106202) and was also supported by the Pioneer Research Center Program through the National Research Foundation (NRF) funded by the Ministry of Science, ICT and Future Planning (NRF-2015M3C1A3023028).

References

N. J. Ullrich and L. B. Gordon, “Hutchinson–Gilford progeria syndrome,” in Neurocutaneous Syndromes, vol. 132 of Handbook of Clinical Neurology, pp. 249–264, Elsevier, 2015.

JAMA, “Progeria, a form of senilism,” Journal of the American Medical Association, vol. XLIII, no. 13, p. 894, 1904.

J. Hutchinson, “Congenital Absence of Hair and Mammary Glands with Atrophic Condition of the Skin and its Appendages in a Boy Whose Mother Had Been Almost Wholly Bald from Alopecia Areata from the Age of Six,” Journal of the Royal Society of Medicine, vol. MCT-69, no. 1, pp. 473–477, 1886.

J. K. Sinha, S. Ghosh, and M. Raghunath, “Progeria: A rare genetic premature ageing disorder,” Indian Journal of Medical Research, vol. 139, no. MAY, pp. 667–674, 2014.

S. L. Chandravanshi, A. K. Rawat, P. C. Dwivedi, and P. Choudhary, “Ocular manifestations in the Hutchinson-Gilford progeria syndrome,” Indian Journal of Ophthalmology, vol. 59, no. 6, pp. 509–512, 2011.

G. I. Stables and W. N. Morley, “Hutchinson-Gilford syndrome,” Journal of the Royal Society of Medicine, vol. 87, pp. 243-244, 1994.

R. L. Pollex and R. A. Hegele, “Hutchinson-Gilford progeria syndrome,” Clinical Genetics, vol. 66, no. 5, pp. 375–381, 2004.

F. Coppedè, “The epidemiology of premature aging and associated comorbidities,” Clinical Interventions in Aging, vol. 8, pp. 1023–1032, 2013.

S. Beauregard and BA. Gilchrest, “Gilchrest BA. Syndromes of premature aging,” Dermatologic Clinics, vol. 5, pp. 109–121, 1987.

S. Pereira, P. Bourgeois, C. Navarro et al., “HGPS and related premature aging disorders: From genomic identification to the first therapeutic approaches,” Mechanisms of Ageing and Development, vol. 129, no. 7-8, pp. 449–459, 2008.

K. Fukuchi, “LMNA mutation in a 45 year old Japanese subject with Hutchinson-Gilford progeria syndrome,” Journal of Medical Genetics, vol. 41, no. 5, pp. e67–e67, 2004.

M. A. Merideth, L. B. Gordon, S. Clauss et al., “Phenotype and Course of Hutchinson–Gilford Progeria Syndrome,” The New England Journal of Medicine, vol. 358, no. 6, pp. 592–604, 2008.

S. H. Yang, M. Meta, X. Qiao et al., “A farnesyltransferase inhibitor improves disease phenotypes in mice with a Hutchinson-Gilford progeria syndrome mutation,” The Journal of Clinical Investigation, vol. 116, no. 8, pp. 2115–2121, 2006.

H.-J. Jung, C. Coffinier, Y. Choe et al., “Regulation of prelamin A but not lamin C by miR-9, a brain-specific microRNA,” Proceedings of the National Acadamy of Sciences of the United States of America, vol. 109, no. 7, pp. E423–E431, 2012.

C.-E. Yu, J. Oshima, Y.-H. Fu et al., “Positional cloning of the Werner's syndrome gene,” Science, vol. 272, no. 5259, pp. 258–262, 1996.

S. Gonzalo, R. Kreienkamp, and P. Askjaer, “Hutchinson-Gilford Progeria Syndrome: A premature aging disease caused by LMNA gene mutations,” Ageing Research Reviews, vol. 33, pp. 18–29, 2017.

D. Constantinescu, H. L. Gray, P. J. Sammak, G. P. Schatten, and A. B. Csoka, “Erratum: Lamin A/C expression is a marker of mouse and human embryonic stem cell differentiation (Stem Cells (2006) 24 (177-185) doi: 10.1634/stemcells. 2004-0159),” Stem Cells, vol. 24, no. 2, p. 474, 2006.

J. Zhang, Q. Lian, G. Zhu et al., “A human iPSC model of hutchinson gilford progeria reveals vascular smooth muscle and mesenchymal stem cell defects,” Cell Stem Cell, vol. 8, no. 1, pp. 31–45, 2011.

M. Eriksson, W. T. Brown, L. B. Gordon et al., “Recurrent de novo point mutations in lamin A cause Hutchinson-Gilford progeria syndrome,” Nature, vol. 423, no. 6937, pp. 293–298, 2003.

K. H. Schreiber and B. K. Kennedy, “When lamins go bad: Nuclear structure and disease,” Cell, vol. 152, no. 6, pp. 1365–1375, 2013.

H. Cao and R. A. Hegele, “LMNA is mutated in Hutchinson-Gilford progeria (MIM 176670) but not in Wiedemann-Rautenstrauch progeroid syndrome (MIM 264090),” Journal of Human Genetics, vol. 48, no. 5, pp. 271–274, 2003.

D. McClintock, D. Ratner, M. Lokuge et al., “The mutant form of Lamin A that causes Hutchinson-Gilford progeria is a biomarker of cellular aging in human skin,” PLoS ONE, vol. 2, no. 12, Article ID e1269, 2007.

H. K. Kim, J. Y. Lee, E. J. Bae et al., “Hutchinson-Gilford progeria syndrome with G608G LMNA mutation,” Journal of Korean Medical Science, vol. 26, no. 12, pp. 1642–1645, 2011.

J. Mazereeuw-Hautier, L. C. Wilson, S. Mohammed et al., “Hutchinson-Gilford progeria syndrome: Clinical findings in three patients carrying the G608G mutation in LMNA and review of the literature,” British Journal of Dermatology, vol. 156, no. 6, pp. 1308–1314, 2007.

A. De Sandre-Giovannoli, R. Bernard, P. Cau et al., “Lamin A truncation in Hutchinson-Gilford progeria,” Science, vol. 300, no. 5628, p. 2055, 2003.

S. G. Young, S. H. Yang, B. S. Davies, H. Jung, and L. G. Fong, “Targeting Protein Prenylation in Progeria,” Science Translational Medicine, vol. 5, no. 171, pp. 171ps3–171ps3, 2013.

B. C. Capell and F. S. Collins, “Human laminopathies: nuclei gone genetically awry,” Nature Reviews Genetics, vol. 7, no. 12, pp. 940–952, 2006.

S. H. Yang, M. O. Bergo, J. I. Toth et al., “Blocking protein farnesyltransferase improves nuclear blebbing in mouse fibroblasts with a targeted Hutchinson-Gilford progeria syndrome mutation,” Proceedings of the National Acadamy of Sciences of the United States of America, vol. 102, no. 29, pp. 10291–10296, 2005.

D.-A. M. Pillers and N. H. Von Bergen, “Emery–dreifuss muscular dystrophy: A test case for precision medicine,” The Application of Clinical Genetics, vol. 9, pp. 27–32, 2016.

G. Bonne and S. Quijano-Roy, “Emery-Dreifuss muscular dystrophy, laminopathies, and other nuclear envelopathies,” Handbook of Clinical Neurology, vol. 113, pp. 1367–1376, 2013.

H. Imachi, K. Murao, S. Ohtsuka et al., “A case of Dunnigan-type familial partial lipodystrophy (FPLD) due to lamin A/C (LMNA) mutations complicated by end-stage renal disease,” Endocrine Journal, vol. 35, no. 1, pp. 18–21, 2009.

S. P. M. Belo, Â. C. Magalhães, P. Freitas, and D. M. Carvalho, “Familial partial lipodystrophy, Dunnigan variety - Challenges for patient care during pregnancy: A case report,” BMC Research Notes, vol. 8, no. 1, article no. 140, 2015.

G. Novelli, A. Muchir, F. Sangiuolo et al., “Mandibuloacral dysplasia is caused by a mutation in LMNA-encoding lamin A/C,” American Journal of Human Genetics, vol. 71, no. 2, pp. 426–431, 2002.

D.-Q. Luo, X.-Z. Wang, Y. Meng et al., “Mandibuloacral dysplasia type A-associated progeria caused by homozygous LMNA mutation in a family from Southern China,” BMC Pediatrics, vol. 14, no. 1, article no. 256, 2014.

U. Kutay and M. W. Hetzer, “Reorganization of the nuclear envelope during open mitosis,” Current Opinion in Cell Biology, vol. 20, no. 6, pp. 669–677, 2008.

L. B. Gordon, M. E. Kleinman, J. Massaro et al., “Clinical Trial of the Protein Farnesylation Inhibitors Lonafarnib, Pravastatin, and Zoledronic Acid in Children with Hutchinson-Gilford Progeria Syndrome,” Circulation, vol. 134, no. 2, pp. 114–125, 2016.

L. B. Gordon, J. Massaro, R. B. D'Agostino Sr. et al., “Impact of farnesylation inhibitors on survival in Hutchinson-Gilford progeria syndrome,” Circulation, vol. 130, no. 1, pp. 27–34, 2014.

F. S. Collins, “Seeking a Cure for One of the Rarest Diseases: Progeria,” in Circulation, vol. 134 of 129, p. 126, 2016.

M. W. Glynn and T. W. Glover, “Incomplete processing of mutant lamin A in Hutchinson-Gilford progeria leads to nuclear abnormalities, which are reversed by farnesyltransferase inhibition,” Human Molecular Genetics, vol. 14, no. 20, pp. 2959–2969, 2005.

Y. Wang, C. Östlund, and H. J. Worman, “Blocking protein farnesylation improves nuclear shape abnormalities in keratinocytes of mice expressing the prelamin A variant in Hutchinson-Gilford progeria syndrome,” Nucleus, vol. 1, no. 5, Article ID 12972, 2010.

S. Blondel, A.-L. Egesipe, P. Picardi et al., “Drug screening on Hutchinson Gilford progeria pluripotent stem cells reveals aminopyrimidines as new modulators of farnesylation,” Cell Death & Disease, vol. 7, no. 2, Article ID e2105, 2016.

H. T. Kang, J. T. Park, K. Choi et al., “Chemical screening identifies ROCK as a target for recovering mitochondrial function in Hutchinson-Gilford progeria syndrome,” Aging Cell, vol. 16, no. 3, pp. 541–550, 2017.

B. C. Capell, M. R. Erdos, J. P. Madigan et al., “Inhibiting farnesylation of progerin prevents the characteristic nuclear blebbing of Hutchinson-Gilford progeria syndrome,” Proceedings of the National Acadamy of Sciences of the United States of America, vol. 102, no. 36, pp. 12879–12884, 2005.

J. I. Toth, S. H. Yang, X. Qiao et al., “Blocking protein farnesyltransferase improves nuclear shape in fibroblasts from humans with progeroid syndromes,” Proceedings of the National Acadamy of Sciences of the United States of America, vol. 102, no. 36, pp. 12873–12878, 2005.

M. Bifulco, A. D'Alessandro, S. Paladino et al., “N6-isopentenyladenosine improves nuclear shape in fibroblasts from humans with progeroid syndromes by inhibiting the farnesylation of prelamin A,” FEBS Journal, vol. 280, no. 23, pp. 6223–6232, 2013.

I. S. Mehta, C. H. Eskiw, H. D. Arican, I. R. Kill, and J. M. Bridger, “Farnesyltransferase inhibitor treatment restores chromosome territory positions and active chromosome dynamics in Hutchinson-Gilford progeria syndrome cells,” Genome Biology, vol. 12, no. 8, article no. R74, 2011.

L. G. Fong, D. Frost, M. Meta et al., “A protein farnesyltransferase inhibitor ameliorates disease in a mouse model of progeria,” Science, vol. 311, no. 5767, pp. 1621–1623, 2006.

I. Varela, S. Pereira, A. P. Ugalde et al., “Combined treatment with statins and aminobisphosphonates extends longevity in a mouse model of human premature aging,” Nature Medicine, vol. 14, no. 7, pp. 767–772, 2008.

L. B. Gordon, M. E. Kleinman, D. T. Miller et al., “Clinical trial of a farnesyltransferase inhibitor in children with Hutchinson-Gilford progeria syndrome,” Proceedings of the National Acadamy of Sciences of the United States of America, vol. 109, no. 41, pp. 16666–16671, 2012.

M. W. Kieran, L. B. Gordon, and M. E. Kleinman, “The role of the farnesyltransferase inhibitor lonafarnib in the treatment of Progeria,” Expert Opinion on Orphan Drugs, vol. 2, no. 1, pp. 95–105, 2014.

N. J. Ullrich, M. W. Kieran, D. T. Miller, L. B. Gordon, Y.-J. Cho, V. M. Silvera et al., “Neurologic features of Hutchinson-Gilford progeria syndrome after lonafarnib treatment,” in Neurology, M. D. Hagerstown, Ed., vol. 81, pp. 427–430, Lippincott Williams & Wilkins, 2013.

S. B. Cammerer, C. Jimenez, S. Jones et al., “Quinuclidine derivatives as potential antiparasitics,” Antimicrobial Agents and Chemotherapy, vol. 51, no. 11, pp. 4049–4061, 2007.

T. Liu, Y. Lin, X. Wen, R. N. Jorissen, and M. K. Gilson, “BindingDB: a web-accessible database of experimentally determined protein-ligand binding affinities,” Nucleic Acids Research, vol. 35, supplement 1, pp. D198–D201, 2007.

Q. A. Al-Balas, H. A. Amawi, M. A. Hassan, A. M. Qandil, A. M. Almaaytah, and N. M. Mhaidat, “Virtual lead identification of farnesyltransferase inhibitors based on ligand and structure-based pharmacophore techniques,” Pharmaceuticals, vol. 6, no. 6, pp. 700–715, 2013.

J. Fei, L. Zhou, T. Liu, and X.-Y. Tang, “Pharmacophore modeling, virtual screening, and mo-lecular docking studies for discovery of novel Akt2 inhibitors,” International Journal of Medical Sciences, vol. 10, no. 3, pp. 265–275, 2013.

P. Leeson, “Drug discovery: Chemical beauty contest,” Nature, vol. 481, no. 7382, pp. 455-456, 2012.

M. P. Gleeson, A. Hersey, D. Montanari, and J. Overington, “Probing the links between in vitro potency, ADMET and physicochemical parameters,” Nature Reviews Drug Discovery, vol. 10, no. 3, pp. 197–208, 2011.

S. Sakkiah and K. W. Lee, “Pharmacophore-based virtual screening and density functional theory approach to identifying novel butyrylcholinesterase inhibitors,” Acta Pharmacologica Sinica, vol. 33, no. 7, pp. 964–978, 2012.

X. Zhou, S. Yu, J. Su, and L. Sun, “Computational Study on New Natural Compound Inhibitors of Pyruvate Dehydrogenase Kinases,” International Journal of Molecular Sciences, vol. 17, no. 3, p. 340, 2016.

S. Sakkiah, S. Thangapandian, S. John, Y. J. Kwon, and K. W. Lee, “3D QSAR pharmacophore based virtual screening and molecular docking for identification of potential HSP90 inhibitors,” European Journal of Medicinal Chemistry, vol. 45, no. 6, pp. 2132–2140, 2010.

S. Rampogu and M. Rampogu Lemuel, “Network Based Approach in the Establishment of the Relationship between Type 2 Diabetes Mellitus and Its Complications at the Molecular Level Coupled with Molecular Docking Mechanism,” BioMed Research International, vol. 2016, Article ID 6068437, 2016.

M. K. Abdel-Hamid and A. McCluskey, “In Silico docking, molecular dynamics and binding energy insights into the bolinaquinone-clathrin terminal domain binding site,” Molecules, vol. 19, no. 5, pp. 6609–6622, 2014.

T. S. Reid, K. L. Terry, P. J. Casey, and L. S. Beese, “Crystallographic analysis of CaaX prenyltransferases complexed with substrates defines rules of protein substrate selectivity,” Journal of Molecular Biology, vol. 343, no. 2, pp. 417–433, 2004.

D. van der Spoel, E. Lindahl, B. Hess, G. Groenhof, A. E. Mark, and H. J. C. Berendsen, “GROMACS: fast, flexible, and free,” Journal of Computational Chemistry, vol. 26, no. 16, pp. 1701–1718, 2005.

V. Zoete, M. A. Cuendet, A. Grosdidier, and O. Michielin, “SwissParam: a fast force field generation tool for small organic molecules,” Journal of Computational Chemistry, vol. 32, no. 11, pp. 2359–2368, 2011.

M. Parrinello and A. Rahman, “Polymorphic transitions in single crystals: a new molecular dynamics method,” Journal of Applied Physics, vol. 52, no. 12, pp. 7182–7190, 1981.

S. Miyamoto and P. A. Kollman, “Settle: an analytical version of the SHAKE and RATTLE algorithm for rigid water models,” Journal of Computational Chemistry, vol. 13, no. 8, pp. 952–962, 1992.

B. Hess, H. Bekker, H. J. C. Berendsen, and J. G. E. M. Fraaije, “LINCS: A linear constraint solver for molecular simulations,” Journal of Computational Chemistry, vol. 18, pp. 1463–1472, 1997.

T. Darden, D. York, and L. Pedersen, “Particle mesh Ewald: an N·log(N) method for Ewald sums in large systems,” The Journal of Chemical Physics, vol. 98, no. 12, pp. 10089–10092, 1993.

W. Humphrey, A. Dalke, and K. Schulten, “VMD: visual molecular dynamics,” Journal of Molecular Graphics, vol. 14, no. 1, pp. 33–38, 1996.

R. Kumari, R. Kumar, and A. Lynn, “g_mmpbsa—A GROMACS tool for high-throughput MM-PBSA calculations,” Journal of Chemical Information and Modeling, vol. 54, no. 7, pp. 1951–1962, 2014.

N. A. Baker, D. Sept, S. Joseph, M. J. Holst, and J. A. McCammon, “Electrostatics of nanosystems: application to microtubules and the ribosome,” Proceedings of the National Acadamy of Sciences of the United States of America, vol. 98, no. 18, pp. 10037–10041, 2001.

N. Homeyer and H. Gohlke, “Free energy calculations by the Molecular Mechanics Poisson-Boltzmann Surface Area method,” Molecular Informatics, vol. 31, no. 2, pp. 114–122, 2012.

P. A. Kollman, I. Massova, C. Reyes et al., “Calculating structures and free energies of complex molecules: combining molecular mechanics and continuum models,” Accounts of Chemical Research, vol. 33, no. 12, pp. 889–897, 2000.

A. Khandelwal, V. Lukacova, D. Comez, D. M. Kroll, S. Raha, and S. Balaz, “A combination of docking, QM/MM methods, and MD simulation for binding affinity estimation of metalloprotein ligands,” Journal of Medicinal Chemistry, vol. 48, no. 17, pp. 5437–5447, 2005.

R. Balajee, V. Srinivasadesikan, M. Sakthivadivel, and P. Gunasekaran, “In silico screening, alanine mutation, and DFT approaches for identification of NS2B/NS3 protease inhibitors,” Biochemistry Research International, vol. 2016, Article ID 7264080, 2016.

U. Abdulfatai, A. Uzairu, and S. Uba, “Quantitative structure-activity relationship and molecular docking studies of a series of quinazolinonyl analogues as inhibitors of gamma amino butyric acid aminotransferase,” Journal of Advanced Research, vol. 8, no. 1, pp. 33–43, 2017.

R. Rakhi and C. H. Suresh, “A DFT study on dihydropyrazine annulated linear polyacenes: aromaticity, stability and HOMO–LUMO energy modulation,” Phys. Chem. Chem. Phys., vol. 18, no. 35, pp. 24631–24641, 2016.

M. Prabhaharan, A. R. Prabakaran, S. Gunasekaran, and S. Srinivasan, “DFT studies on vibrational spectra, HOMO-LUMO, NBO and thermodynamic function analysis of cyanuric fluoride,” Spectrochimica Acta Part A: Molecular and Biomolecular Spectroscopy, vol. 136, pp. 494–503, 2015.

M. Bourass, A. T. Benjelloun, M. Benzakour et al., “DFT and TD-DFT calculation of new thienopyrazine-based small molecules for organic solar cells,” Chemistry Central Journal, vol. 10, no. 1, article no. 67, 2016.

C. Lee, W. Yang, and R. G. Parr, “Development of the Colle-Salvetti correlation-energy formula into a functional of the electron density,” Physical Review B: Condensed Matter and Materials Physics, vol. 37, no. 2, pp. 785–789, 1988.

A. D. Becke, “A multicenter numerical integration scheme for polyatomic molecules,” The Journal of Chemical Physics, vol. 88, no. 4, pp. 2547–2553, 1988.

B. R. Kumar and X. M. Suresh, “Pharmacophore mapping based inhibitor selection and molecular interaction studies for identification of potential drugs on calcium activated potassium channel blockers, tamulotoxin,” Pharmacognosy Magazine, vol. 9, no. 34, pp. 89–95, 2013.

F. Ntie-Kang, C. V. Simoben, B. Karaman et al., “Pharmacophore modeling and in silico toxicity assessment of potential anticancer agents from African medicinal plants,” Drug Design, Development and Therapy, vol. 10, pp. 2137–2154, 2016.

S.-H. Lu, J. W. Wu, H.-L. Liu et al., “The discovery of potential acetylcholinesterase inhibitors: a combination of pharmacophore modeling, virtual screening, and molecular docking studies,” Journal of Biomedical Science, vol. 18, no. 1, article 8, 2011.

L. L. Duan, G. Q. Feng, and Q. G. Zhang, “Large-scale molecular dynamics simulation: Effect of polarization on thrombin-ligand binding energy,” Scientific Reports, vol. 6, no. 1, 2016.

X.-B. Qiao, B. Jiang, T.-J. Hou, and X.-J. Xu, “Representation of Molecular Electrostatic Potentials of Biopolymer by Self-organizing Feature Map,” Chinese Journal of Chemistry, vol. 19, no. 12, pp. 1172–1178, 2001.

M. Arooj, S. Sakkiah, S. Kim, V. Arulalapperumal, and K. W. Lee, “A Combination of Receptor-Based Pharmacophore Modeling & QM Techniques for Identification of Human Chymase Inhibitors,” vol. 8, 2013.

S. Sakkiah, M. Arooj, M. R. Kumar, S. H. Eom, and K. W. Lee, “Identification of Inhibitor Binding Site in Human Sirtuin 2 Using Molecular Docking and Dynamics Simulations,” PLoS ONE, vol. 8, no. 1, Article ID e51429, 2013.

R. Haunschild, A. Barth, and W. Marx, “Evolution of DFT studies in view of a scientometric perspective,” Journal of Cheminformatics, vol. 8, no. 1, 2016.

F. A. Momany, “Determination of partial atomic charges from ab initio molecular electrostatic potentials. Application to formamide, methanol, and formic acid,” The Journal of Physical Chemistry C, vol. 82, no. 5, pp. 592–601, 1978.

L. E. Chirlian and M. M. Francl, “Atomic charges derived from electrostatic potentials: A detailed study,” Journal of Computational Chemistry, vol. 8, no. 6, pp. 894–905, 1987.

T. A. Manz and D. S. Sholl, “Chemically meaningful atomic charges that reproduce the electrostatic potential in periodic and nonperiodic materials,” Journal of Chemical Theory and Computation, vol. 6, no. 8, pp. 2455–2468, 2010.

D. E. Williams and A. Abraha, “Site charge models for molecular electrostatic potentials of cycloalkanes and tetrahedrane,” Journal of Computational Chemistry, vol. 20, pp. 579–585, 1999.

Cox S. R. and D. E. Williams, “Representation of the molecular electrostatic potential by a net atomic charge model,” Journal of Computational Chemistry, vol. 2, pp. 304–323, 1981.

R. S. Mulliken, “Electronic population analysis on LCAO-MO molecular wave functions. IV. bonding and antibonding in LCAO and valence-bond theories,” The Journal of Chemical Physics, vol. 23, no. 12, pp. 2343–2346, 1955.

J. A. Montgomery Jr., M. J. Frisch, J. W. Ochterski, and G. A. Petersson, “A complete basis set model chemistry. VII. Use of the minimum population localization method,” The Journal of Chemical Physics, vol. 112, no. 15, pp. 6532–6542, 2000.