Accepted Manuscript Title: Insights into the functions of M-T hook structure in HIV fusion inhibitor using molecular modeling Author: Jianjun Tan Hongling Yuan Chunhua Li Xiaoyi Zhang Cunxin Wang PII: DOI: Reference:
S1476-9271(16)30043-3 http://dx.doi.org/doi:10.1016/j.compbiolchem.2016.01.006 CBAC 6500
To appear in:
Computational Biology and Chemistry
Received date: Revised date: Accepted date:
19-1-2015 15-1-2016 21-1-2016
Please cite this article as: Tan, Jianjun, Yuan, Hongling, Li, Chunhua, Zhang, Xiaoyi, Wang, Cunxin, Insights into the functions of M-T hook structure in HIV fusion inhibitor using molecular modeling.Computational Biology and Chemistry http://dx.doi.org/10.1016/j.compbiolchem.2016.01.006 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Insights into the Functions of M-T Hook Structure in HIV Fusion Inhibitor Using Molecular Modeling Jianjun Tan* [email protected], Hongling Yuan, Chunhua Li* [email protected], Xiaoyi Zhang, Cunxin Wang College of life science and Bioengineering, Beijing University of Technology, Beijing, 100124, China * Corresponding authors: Tel: +86-10-67392724
Graphical Abstract
Highlights We explored the mechanism why M-T could improve the activity of inhibitors. The system with the M-T hook structure formed more stable Hydrogen bonds than without M-T hook structure. The M-T hook could further fortify the interaction by “hooking” the pocket residues tightly; in contrast, this structure do not appear in the system without the M-T hook. Binding free energy of the ligand with M-T was higher than the other without M-T.
Abstract HIV-1 membrane fusion plays an important role in the process that HIV-1 entries host cells. As a treatment strategy targeting HIV-1 entry process, fusion inhibitors have been proposed. Nevertheless, development of a short peptide possessing high anti-HIV potency is considered a daunting challenge. He et al. found that two residues, Met626 and Thr627, located the upstream of the C-terminal heptad repeat of the gp41, formed a unique hook-like structure (M-T hook) that can dramatically improve the binding stability and anti-HIV activity of the inhibitors. In this work, we explored the molecular mechanism why M-T hook structure could improve the anti-HIV activity of inhibitors. Firstly, molecular dynamic simulation was used to obtain information on the time evolution between gp41 and ligands. Secondly, based on the simulations, molecular mechanics Poisson-Boltzmann surface area (MM-PBSA) and molecular mechanics Generalized Born surface area (MM-GBSA) methods were used to calculate the binding free energies. The binding free energy of the ligand with M-T hook was considerably higher than the other without M-T. Further studies showed that the hydrophobic interactions made the dominant contribution to the binding free energy. The numbers of Hydrogen bonds between gp41 and the ligand with M-T hook structure were more than the other. These findings should provide insights into the inhibition mechanism of the short peptide fusion inhibitors and be useful for the rational design of novel fusion inhibitors in the future. Keywords: Molecular dynamic simulation; MM-PB/GBSA; Binding free energy; HIV-1; Fusion inhibitors
INTRODUCTION Acquired immune deficiency syndrome (AIDS), which is induced by human immunodeficiency virus (HIV) infection, has become one of the major medical and humanitarian challenges(Liu, et al. 2011). So far, the US Food and Drug Administration (FDA) approves four classes of anti-HIV drugs: (i) reverse transcriptase inhibitors (RTIs), (ii) protease inhibitors (PIs), (iii) fusion and entry inhibitors, and (iv) integrase inhibitors (Maguire, et al. 2014). In particular, HIV-1 fusion/entry inhibitors can inhibit early steps of the HIV replication cycle, thereby inhibiting the spread of infection in the first time (Lin and Xiong 2013). Fusion of HIV virus membrane with host cell membrane is an essential event for virus infection and proliferation. This process is mediated by HIV envelope glycoprotein 120 (gp120) (Xu, et al. 2013) and glycoprotein 41 (gp41) (AIhemaiti·abudureyimu, et al. 2013), thus, it is served as an attractive target for designing a novel HIV-1fusion/entry inhibitors(Wang, et al. 2013c). T20 (Enfuvirtide, Fuzeon) is the first HIV fusion inhibitor approved by U.S. FDA. Although T20 is effective (Liang, et al. 2012; Tan, et al. 2013), T20 has some shortcomings that limit its clinical usefulness. Firstly, owing to its nature characteristic, T20 is marketed in an injectable form. Patients are subcutaneously injected twice in every day. Because of the chronic nature of the disease therapy, this dosage form may be one of major problems for the patient's adherence to this drug regimen. Secondly, T20 clinical efficacy is also negatively affected by gp41 resistance mutations that arise owing to continued use. To overcome these disadvantages, a number of strategies have been applied to discover the next generation inhibitors with improved potency and pharmaceutical properties(Wang, et al. 2013b). However, the resulting peptides usually inherit a longer sequence(Gao, et al. 2013), such as SC35EK (35-mer)(Xie, et al. 2012b), SFT (sifuvirtide, 36-mer) (Xie, et al. 2012a), and T2635 (39-mer)(Zhang, et al. 2011). Undoubtedly, the length of peptide closely relates to difficulty and cost of its production, but short peptides (Amino acid residue numbers < 14) lead to a dramatic loss of antiviral potency and binding stability(Wang, et al. 2013a). Therefore, the development of peptide’s fusion inhibitor with highly active antiretroviral efficacy and reasonable residue numbers faces great challenges(Panthong, et al. 2013). Promisingly, He Y.’s group discoved a novel strategy to solve above problem. By studying the inhibit's mechanism of gp41, they found that two key residues,
Met626 and Thr627, which forms a hook-like structure (M-T hook), located upstream of the pocket-binding domain (PBD) of the C-terminal heptad repeat (CHR). The M-T hook can dramatically improve the binding stability of the inhibitors(Yang, et al. 2011a; Yang, et al. 2011b); and the addition of the M-T hook structure to the N terminus of C-peptide could markedly enhance its antiviral activity(Liu, et al. 2013; Yang, et al. 2011a). Based on their former research, Chong H. et al proposed a prototype that the short C-peptide fusion inhibitors could be developed by adding an M-T hook structure. They introduced the methionine and threonine to the N terminus of a panel of C peptides (24-mer or less), and found one pepide, named MT-SC22EK, was highly active against diverse HIV-1 variants, including those T20- and SC29EK-resistant variants. Based on the crystal structures of MT-SC22EK, their paper reported that MT-SC22EK primarily targets on the hydrophobic pocket region on N-heptad repeat (NHR) trimer with its PBD and M-T hook regions, whereas the sites conferring T20-resistance are not bound by the inhibitor. Therefore, the M-T hook structure provides a novel strategy for designing short peptide fusion inhibitors that can overcome the problem of drug resistance. The major problem affected activity of HIV-1 fusion inhibitors is molecular resistance. Owing to point mutations in the GIV motif of NHR (L33 to L45), the inhibit activity of T20 is on the decrease (Liu, et al. 2005). Some studies reported that mutations on NHR of gp41 would induce the compensatory mutation on CHR, which may balance those damages and remedy the function of gp41(Xu, et al. 2005). In our former work, we studied the two kinds of mutants systems of gp41, N43D single mutation and compensatory mutation S138A (yielding N43D/S138A), respectively, and the analyses of structural alteration indicated that the inhibitor C34 displays different binding modes in wild type and two mutants(AIhemaiti·abudureyimu, et al. 2013). In N43D single mutation system, it was founded that the N43D mutation introduced a negative charge on receptor, thus led to repulsion between receptor and ligand, and induced the conformation changes of complex, and led to greater losses in hydrophobic interactions and hydrogen bond interaction. In N43D/S138A double mutation system, we founded CHR have more contact with hydrophobic residues of receptor, induced more van der Waals contributions. This great hydrophobic interaction of double mutant could partially offset the energy loss arising from N43D mutation. The change of binding mode on N43D/S138A could explain why compensatory mutation S138A could partially restores infectivity, lead to increased
infectivity and greater resistance to inhibitors than mutation N43D. In order to explore the role of M-T hook structure in improving the activity of fusion inhibitors, we used molecular dynamics (MD) simulations and molecular mechanics Poisson-Boltzmann/ Generalized Born surface area (MM-PB/GBSA) approach to analyze structure-function relationships between gp41 and two ligands (SC22EK and MT-SC22EK) (figure. 1). The only difference of the two ligands is that MT-SC22EK includes M-T hook, the other without. These results may be useful for designing novel HIV fusion inhibitors (Tan, et al. 2016), or optimizing the lead compounds.
MATERIALS AND METHODS Preparation of protein system The initial structures were downloaded from the RSCB protein data bank (www.rcsb.org) (PDB ID: 3VU5, 3VU6). The amino acid sequences of the two ligands are shown in Figure 1. The receptor was one N-terminal heptad repeat (NHR) chain of gp41. We named the system including NHR and MT-SC22EK as complex 1, the other as complex 2.
Molecular dynamic simulation The MD simulations were carried out for each complex system using NAMD2.9 software and amber ff03 force-field. Firstly, the each system was solvated with TIP3P (Jorgensen, et al. 1983) water box, 10 Å of space between the protein and the simulation box boundaries on each direction. After model generation, each system was neutralized separately with Cl- ions. Then, two successive minimization steps were performed: (i) the system was minimized by 20000 steps using the conjugate gradient with 2.0 kcal/mol restraints on the backbone of the protein; (ii) the system was minimized by 20000 steps using the conjugate gradient without constrained atoms. Next, the systems were gradually heated from 10 K to 300 K using Langevin dynamics in a slow, stepwise fashion (20 K every 100 ps). During heating, no restraint was utilized in the complex. Finally, the MD simulations were performed in the NPT ensemble (a constant temperature of 300K and a pressure of 1 atm) for 50 ns (the time step is 2fs) with all restraints removed. Throughout the MD simulations, a cut-off of 12 Å was used for nonbonded interactions while the SHAKE and Particle Mesh
Ewald (PME) were used with periodic boundary conditions(Ryckaert, et al. 1977).
Binding free energy calculation and energy decomposition The binding free energies were calculated by subtracting the free energies of the unbound receptor and ligand from the free energy of the bound complex, shown in eq (1):(Kollman, et al. 2000; Kuhn and Kollman 2000)
G Gcomplex Greceptor Gligand where
Gcomplex ,
Greceptor
, Gligand
(1)
are the average values of the Gibbs free energy for
the complex, receptor, and ligand, respectively. Binding free energies were calculated using the MMPBSA.py approaches as implemented in AMBER12 and AmberTools13. MMPBSA(Zhang, et al. 2012) computed the binding free energy by using a thermodynamic cycle that combined the molecular mechanical energies with the solvated approaches.
MM MM sol sol GMMPB (GB ) SA Eele EvdW Gnonpol Gpol
(2)
sol sol Where the Gnonpol and G pol are nonpolar and polar terms in solution, respectively.
MM MM Eele and EvdW are the electrostatic and van der Waals contributions to the
molecular mechanics energy, respectively. The polar part expresses the electrostatic contribution to solvation, which is obtained by solving the linear Poisson Boltzmann equation in a continuum model of solvent. On the other hand, the other part accounted for the nonpolar contribution is related linearly to the solvent accessible surface area. The solvation free energies are calculated using an implicit solvent model. Several implicit solvent models are available to calculate the solvation free energies, including Generalized-Bron (GB), Poisson-Boltzmann (PB), and Reference Interaction Site Model (RISM) etc. In this experiment, we calculated the solvation free energies with the PB and GB methods. sol Using the eq. (3) calculated the nonpolar solvation energy, Gnonpol : sol Gnonpol A
(3)
where A is the solvent accessible surface area (SASA), using the MSMS program (Sanner, et al. 1996) with a 1.4 Å radius probe; γ is the surface tension, and β is the offset. The value of constants γ and β are 0.0052 kcal mol−1 Å−2 and 0.92 kcal mol−1,
respectively. The nonpolar solvation term in the MM(GB)PBSA method combines an implicit estimate of the entropic changes associated with the insertion of a solute into the solvent. However, it is obvious that the entropic impact of changes in the configurational freedom of the protein and ligand upon complex formation in vacuo. As a general rule, in the process of protein−ligand binding, the restrictions to the number of conformations causes a reduction in entropy; this contribution is known as the configurational (or conformational) entropy. For some systems, some researchers have found that the calculated binding affinity values is more colse experimental binding affinity values by incorporating a normal mode estimation of the entropic component of the binding free energy alongside MMPB(GB)SA(Men and Guo 2010; Tu 2011). Including this contribution, the final binding affinity estimate ΔGtheor is given by
Gtheor PB (GB ) GMMPB (GB ) SA T S
(4)
where GMMPB (GB ) SA is the MMPBSA or MMGBSA binding free energy, T is the thermodynamic temperature, and −TΔS is the normal mode estimate of the configurational entropy penalty of the binding process of ligand. As implemented within the AMBER 12 package employed in this study, S is calculted as eq. (5): S Stra Srot Svib
(5)
where ΔStra, ΔSrot, and ΔSvib are the contributions to changes in translational, rotational, and vibrational freedoms, respectively.To explore the impact of the binding residue for binding energies, the intermolecular interactions between each residue of the NHR of gp41 and inhibitors were calculated in terms of the pair interaction decomposition of free energy(Li, et al. 2013). The calculations were also performed over MD trajectories using the decomposition energy module in AMBER12 and AmberTools13. The steps of the experiment are shown in Figure 2.
RESULTS AND DISCUSSIONS The stability of the simulation systems The root-mean-square deviation (RMSD) of the protein backbone residues in regard
to the initial structure versus simulation time is shown in Figure 3. The results indicate that the RMSD of the complex1 fluctuates reaching a value of about 3 Å at 10 ns and keep stable throughout the course of the trajectory. The RMSD of the complex 2 fluctuates reaching a value of 4 Å at 10 ns and remains stable from 10 ns to 30 ns. From 30 ns to 35 ns, an increasing trend in RMSD behavior is observed. Afterwards complex 2 seems to relatively stable for the next 15 ns about 5 Å. Overall, the fluctuation of RMSD of the complex 2 is larger than the complex 1, which means that the complex 2 is relatively unstable (figure 3). To know the extent of
variation of individual residues of inhibitor (MT-SC22EK and SC22EK) and receptor respectively, in complex systems, the per-residue root-mean-square fluctuations (RMSFs) of backbone atoms were computed. In the case of receptors, we observed the NO.571-590 residue of the complex 1 shows higher RMSF than that found in complex 2 (Fig. 4). On the other hand, RMSF of the two inhibitors is same. The RMSD dissimilarity matrix about the structure of the complex along the trajectory (from the beginning of the optimization to 10 ns MD simulations).The blue color indicates the correlation; the red color represents the uncorrelation (figure 5). To see this can provide some implications: (1) RMSD of complex 1 reaches 3Å (Figure 5a), and complex 2 reaches 4Å (Figure 5b). This also explains that compared with complex 2, complex 1 is more stable. (2) The similarity of the initial structure and the MD stable structure is weak, which means that the final structures obtained from the stable trajectory that had no relation to the starting structure.
Binding free energy calculation Contributions of the binding free energies are summarized in Table 1. The binding free energy of the complex 1 and 2 is -46.70 kcal/mol -38.17 kcal/mol, respectively. The binding free energy of the complex 1 is considerably higher than the complex 2. In accordance with the energy components of the binding free energy (Table 1), we found that the electrostatic and van der Waals energies contributed to strong favorable sol binding of the receptor with the ligand. The non-polar solvation Energies ( Gnonpol )
have similar small favorable contribution to the binding affinity in the two complexes. The entropy (-TS ) for the two complexes is rather close because the structures of the two ligands are very similar (Figure 6, Table 1).
Ranking aggregation result of the per-residue decomposition list Many key residues contribute to binding free energies that ensure the stability of the protein-ligand systems. The binding free energy was decomposed into per-residue to create a per-residue spectrum (figure 7a, figure 7b). The approach of residue decomposition is tremendously useful to clarify the inhibition mechanism at a molecular level and helpful to locate contribution of individual residue to the ligand-receptor interactions, too. The decomposition of Gtheor values on a MM per-residue basis into contributions from EvdW , the sum of electrostatic energy in MM sol gas phase and polar solvation energy ( Eele + G pol ), and non-polar solvation energy sol ( Gnonpol ). It can be seen that the van der Waals energy contributes mostly to the
binding free energies. To further get insight into the difference contribution between the complex 1 and the complex 2 in the van der Waals energies, we calculated each residue’s contribution to the van der Waals (figure 8). We can see: Ile573, Leu576, Gln577 (in the receptor) and Trp623, Trp626, Tyr633 (in the ligand) have a positive contribution to the van der Waals energy in complex 1; In contrast, the contribution to the binding energy of the residues significantly decrease in complex 2. Overall, the van der Waals interaction of each residues in the two systems has a significantly difference.
The hydrogen bond analysis The hydrogen bonds (H-bonds) that occupied more than 20% are listed in Table 2. It shows the system of complex 1 has more H-bonds than the other system. In the system of complex 1, there are two inter-hydrogen bonds (Trp626 and Thr567, Trp623 and Ile571) and two intra-hydrogens (Thr634 and Ile630, Ser644 and Leu640). However, in the other system, there only two intra-hydrogen bonds (Thr632 and Ile628, Ser642 and Leu638). The M-T structure is conducive to form more H-bonds between receptor and the ligand (Table 2). This also explains the reason why the
Trp623 and Trp624 have a positive contribution to the van der Waals energy in complex 1 (AIhemaiti·abudureyimu, et al. 2013).
The analysis of the contact residues We analyzed the contact residues near the hydrophobic pocket to further explore the differences of the van der Waals contribution in the two systems (figure 9). The pocket is formed by a cluster of hydrophobic residues, including Leu565, Leu566, Leu568, Thr569, Val570, Trp571, Gly572, Ile573, Lys574, Leu576, and Gln577(Barlow, et al. 2012). According to the experiment reported, when the distance between molecules is greater than 8.5Å, the van der Waals forces between the molecules is negligible(Menendez-Arias and Alvarez 2014). We calculated the number of contact residues between the hydrophobic pocket of gp41 and the ligands (Table 3). By comparison, we found that there are more contact residues in the complex 1 than the other (Table 3). From the binding model of the gp41 with the ligand (figure 8), we find that the MT-SC22EK inserts the groove of gp41 completely. However, the C-terminal of the SC22EK tend to be away from the combination (figure 8b). The analysis of the contact residues reveals that the M-T hook can further fortify the interaction by “hooking” the pocket residues tightly; in contrast, this structure do not appear between SC22EK and the NHR of gp41. To a certain extent, it can explain why the van der Waals interactions between MT-SC22EK (with M-T) and NHR of gp41 are a few stronger than that of the other.
CONCLUSION The primary objective of this work focused on exploring the molecular mechanism that M-T hook structure in fusion inhibitors improved the anti-HIV activity of inhibitors. MM-PBSA (GBSA) was used in conjunction with normal mode analysis to successfully account for experimental binding free energies between the gp41 protein and MT-SC22EK (or SC22EK). Average energies were obtained from samples of 100 snapshots extracted from MD simulations. Free energy decomposition on a residual basis indicates differences between MT-SC22EK and SC22EK. We found the M-T hook structure was contributed to the van der Waals energy, which was favorable to the binding energy. The analysis of the H-bond indicated that the system with the M-T
hook structure formed more stable H-bonds than the system without M-T hook structure, and these hydrophobic interactions also played a key role in the binding process. The analysis of the contact residues revealed that the M-T hook could further fortify the interaction by “hooking” the pocket residues tightly; in contrast, this structure do not appear in the system without the M-T hook. This work has provided the detailed molecular mechanism that M-T hook structure in fusion inhibitors. In order to improve lead compounds activity, during drug design, these interactions should be considered carefully. Our research illustrates that the interactions between gp41 and short peptide ligands may lead to possible anti-HIV drug design in the future.
Acknowledgements The authors thank Chinese Natural Science Foundation project (No. 21173014, and 11474013).
References AIhemaiti·abudureyimu, et al., 2013. Clinical Curative Effect of Traditional Chinese Medicine on Immune Function of 20 Cases HIV /AIDS Patients. Chinese Archives of Traditional of Chinese Medicine, 31, 491-493. Barlow, D. J., et al., 2012. In-silico studies in Chinese herbal medicines' research: evaluation of in-silico methodologies and phytochemical data sources, and a review of research to date. Journal of Ethnopharmacology, 140, 526-34. Gao, Guojian, et al., 2013. Research of Chinese herbal compound on promoting HIV immune reconstitution China Journal of Chinese Materia Medica, 38, 2523-2526. Jorgensen, W. L., et al., 1983. Comparison of Simple Potential Functions for Simulating Liquid Water. Journal of Chemical Physics, 79, 926-935. Kollman, P. A., et al., 2000. Calculating structures and free energies of complex molecules: Combining molecular mechanics and continuum models. Accounts of Chemical Research, 33, 889-897. Kuhn, B., and P. A. Kollman, 2000. Binding of a diverse set of ligands to avidin and streptavidin: An accurate quantitative prediction of their relative affinities by a combination of molecular mechanics and continuum solvent models. Journal of Medicinal Chemistry, 43, 3786-3791. Li, R., et al., 2013. Imbalanced network biomarkers for traditional Chinese medicine Syndrome in gastritis patients. Sci Rep, 3, 1543. Liang, Feili, et al., 2012. Effect of traditional Chinese medicine combined with acyclovir on herpes zoster in elderly with HIV/AIDS. Modern Journal of Integrated Traditional Chinese and Western Medicine, 21, 917-918. Lin, Jiang, and Yu Xiong, 2013. Current clinical research and prospect on Chinese medicine for AIDS. China Journal of TCM and Pharmacy, 28, 2070-2073. Liu, S. W., et al., 2005. Different from the HIV fusion inhibitor C34, the anti-HIV drug fuzeon (T-20) inhibits HIV-1 entry by targeting multiple sites in gp41 and gp120. Journal of Biological Chemistry, 280, 11259-11273. Liu, Ying, et al., 2011. Retrospective Analysis of Clinical Literature on Chinese Medicine for AIDS. Journal of Traditional Chinese Medicine
54, 1686.
Liu, Zhen, et al., 2013. Effect of Immune 2 with highly active antiretroviral treatment on immune function of HIV/AIDS patients with poor immune reconstitution China Journal of Chinese Materia Medica, 38, 2458-2462. Maguire, J. J., et al., 2014. The CCR5 chemokine receptor mediates vasoconstriction and stimulates intimal hyperplasia in human vessels in vitro. Cardiovascular Research. Men, Jiuzhang, and Lei Guo 2010 A general introduction to traditional Chinese medicine: Science Press & CRC Press. Menendez-Arias, L., and M. Alvarez, 2014. Antiretroviral therapy and drug resistance in human immunodeficiency virus type 2 infection. Antiviral Research, 102, 70-86. Panthong, K., et al., 2013. Benzene, coumarin and quinolinone derivatives from roots of Citrus hystrix. Phytochemistry, 88, 79-84. Ryckaert, J. P., G. Ciccotti, and H. J. C. Berendsen, 1977. Numerical-Integration of Cartesian Equations of Motion of a System with Constraints - Molecular-Dynamics of N-Alkanes. Journal of Computational Physics, 23, 327-341.
Sanner, M. F., A. J. Olson, and J. C. Spehner, 1996. Reduced surface: An efficient way to compute molecular surfaces. Biopolymers, 38, 305-320. Tan, J. J., et al., 2013. The Current Status and Challenges in the Development of Fusion Inhibitors as Therapeutics for HIV-1 Infection. Current Pharmaceutical Design, 19, 1810-1817. Tan, J., et al., 2016. Design, synthesis and activity evaluation of novel peptide fusion inhibitors targeting HIV-1 gp41. Bioorg Med Chem, 24, 201-206. Tu, Y., 2011. The discovery of artemisinin (qinghaosu) and gifts from Chinese medicine. Natural Medicines, 17, 1217-20. Wang, J., W. Zou, and Y. Liu, 2013a. [Analysis on efficacy action points of HIV/AIDS treatment with traditional Chinese medicine and antiretroviral drugs]. Zhongguo Zhong Yao Za Zhi (China Journal of Chinese Materia Medica), 38, 2504-6. Wang, Jian, et al., 2013b. A sixty-month study of the change of CD4T－lymphocyte cell count treated with traditional Chinese medicine in 565HIV / AIDS patients. Chinese Journal of AIDS & STD, 19, 468-470. Wang, Jie, et al., 2013c. Study of drug intervention in diversity changes of TCRVβ gene in AIDS patients with incomplete immune reconstitution. China Journal of Chinese Materia Medica, 38, 2424-2428. Xie, Zheng, et al., 2012a. The Clinical Research of TCM on HIV/AIDS Depression Traditional Chinese Medicine, 6, 172-173. Xie, Zheng, Li. Pengyu, and Ziqiang Jiang, 2012b. Clinical Observation on 26 Cases of HIV -pruritic Popular Eruption. Chinese Medicine Modern Distance Education of China, 10, 10-12. Xu, L., et al., 2005. Emergence and evolution of enfuvirtide resistance following long-term therapy involves heptad repeat 2 mutations within gp41. Antimicrobial Agents and Chemotherapy, 49, 1113-1119. Xu, Liran, et al., 2013. Study on quality of life of asymptomatic HIV infected persons with traditional Chinese medicine. China Journal of Chinese Materia Medica, 38, 2480-2483. Yang, Liumeng, et al., 2011a. In vitro anti-HIV-1 activities of Qishile,a Chinese medicine effective fraction formula. Chinese Pharmacological Bulletin, 27, 566-571. Yang, Shaochun, et al., 2011b. A Clinical Study on the Treatment of Methadone Maintenance HIVinfected Persons with Fuzheng Ant-i Drug Capsules. Yunnan Journal of Traditional ChineseMedicine and Materia Medica, 32, 18-20. Zhang, Aihua, et al., 2012. Future perspectives of personalized medicine in traditional Chinese medicine: A systems biology approach. Complementary Therapies in Medicine, 20, 93-99. Zhang, Ying, et al., 2011. Chinese Medicine Treatment of HIV Infection. Chinese Archives of Traditional Chinese Medicine, 29, 2628-2630.
Figure Captions Figure 1. Positively charged lysine residues are in blue; negatively charged glutamic acid residues are in red. Black solid lines indicate salt bridges observed in the structure. Green solid line between Thr and Glu indicates hydrogen bond. The sequences of the SC22EK and MT-SC22EK have only one difference is the two residues in the front of MT-SC22EK
Figure 2. General workflow for performing steady state calculated with MMPBSA.py. tleap is a program in Amber used to create topology files for dynamics. “Dry” topologies and ensembles are systems without explicit solvent that is subsequently treated using an implicit solvent model (Li, et al. 2013)
Figure 3. Backbone RMSD are shown as a function of time for complex 1 and complex 2 using the ptraj module in AMBER10
Figure 4. Comparisons of backbone atoms RMSF in inhibitor SC22EK and MT-SC22EK and gp41 NHR
Figure 5. RMSD dissimilarity plot for complex 1 (a) and complex 2 (b). The horizontal and vertical coordinate of the figure indicate the label of simulated trajectories
Figure 6. Energy components (kcal/mol) for the binding of complex 1 and complex 2
Figure 7. The decomposition of G on a per-residue basis into contributions from the MM van der Waals energy ( EvdW ), the sum of electrostatic interactions and polar
MM sol sol solvation energy ( Eele + G pol ), and nonpolar solvation energy ( Gnonpol ) for the
residues of G
MM Figure 8. Energy shown as the contribution from the van der Waals energy ( EvdW )
Figure 9. The binding model of the gp41 with the ligand. (a) The model of complex 1; (b) The model of complex 2. The structural composites from the average structure of the 20-50 ns molecular dynamic simulation trajectories.
Tables Table1. Binding free energies for two complexes by using the MM-PB (GB) SA method (800 snapshots from 20 ns-50 ns)
Component
Complex 1
Complex 2
MM Eele
-56.15 (0.09)*
MM EvdW
-62.47 (0.48)*
-48.09 (0.08)*
sol G pol ( GB )
85.19 (0.47)*
171.87 (0.56)*
sol G pol ( PB )
79.16 (0.47)*
165.08 (0.57)*
sol Gnonpol ( GB )
-8.35 (0.01)*
-7.33 (0.01)*
sol Gnonpol ( PB )
-56.15 (0.09)*
-TS
-149.39
-149.39
(0.62)*
(0.62)*
33.26 (0.27)*
32.97 (0.38)*
Gtheor GB
-41.74
(0.16)*
-32.99 (0.25)*
Gtheor PB
-46.70
(0.16)*
-38.17 (0.24)*
Gexpt(Liu, et al. 2013)
-11.71
All values are given in kcal/mol * standard error of mean values
-9.53
Table2. Hydrogen-bond analysis from the results of MD simulation System
Donor
Accepter
Thr639-HG1-OG1
Ile635-O
Ser649-HGI-OG
Distance
(Å)
Angle (Degrees)
Occupied（%）
2.72 ( 0.11)*
164.61 (8.63)*
97.12
Leu645-O
2.74 ( 0.11)*
162.42 (10.32)*
74.75
Trp631-HE1-NE1
Thr569-O
2.77 ( 0.17)*
160.00 (10.64)*
43.83
Trp628-HE1-NE1
Ile573-O
2.87 ( 0.20)*
153.85(12.13)*
21.83
Thr639-HG1-OG1
Ile635-O
2.73 ( 0.11)*
164.27 (8.85)*
99.42
Ser649-HG-OG
Leu645-O
2.74 ( 0.11)*
160.98 (9.94)*
87.83
Complex 1
Complex 2
The H-bonds are determined by the donor…acceptor atom distance of ≤3.5 Å and acceptor…H donor angle of ≥120° * standard error
Table 3. The hydrophobic interactions between the hydrophobic pocket and the ligand
Residue
MT-SC22EK
SC22EK
Leu565
8*
6*
Leu566
14*
13*
Leu568
6*
3*
Thr569
13*
7*
Val570
13*
10*
Trp571
4*
5*
Gly572
8*
3*
Ile573
11*
9*
Lys574
5*
6*
Leu576
10*
2*
Gln577
8*
7*
The hydrophobic interactions with distance of less than 8.5Å. * The number of hydrophobic contact
S1476-9271(16)30043-3 http://dx.doi.org/doi:10.1016/j.compbiolchem.2016.01.006 CBAC 6500
To appear in:
Computational Biology and Chemistry
Received date: Revised date: Accepted date:
19-1-2015 15-1-2016 21-1-2016
Please cite this article as: Tan, Jianjun, Yuan, Hongling, Li, Chunhua, Zhang, Xiaoyi, Wang, Cunxin, Insights into the functions of M-T hook structure in HIV fusion inhibitor using molecular modeling.Computational Biology and Chemistry http://dx.doi.org/10.1016/j.compbiolchem.2016.01.006 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Insights into the Functions of M-T Hook Structure in HIV Fusion Inhibitor Using Molecular Modeling Jianjun Tan* [email protected], Hongling Yuan, Chunhua Li* [email protected], Xiaoyi Zhang, Cunxin Wang College of life science and Bioengineering, Beijing University of Technology, Beijing, 100124, China * Corresponding authors: Tel: +86-10-67392724
Graphical Abstract
Highlights We explored the mechanism why M-T could improve the activity of inhibitors. The system with the M-T hook structure formed more stable Hydrogen bonds than without M-T hook structure. The M-T hook could further fortify the interaction by “hooking” the pocket residues tightly; in contrast, this structure do not appear in the system without the M-T hook. Binding free energy of the ligand with M-T was higher than the other without M-T.
Abstract HIV-1 membrane fusion plays an important role in the process that HIV-1 entries host cells. As a treatment strategy targeting HIV-1 entry process, fusion inhibitors have been proposed. Nevertheless, development of a short peptide possessing high anti-HIV potency is considered a daunting challenge. He et al. found that two residues, Met626 and Thr627, located the upstream of the C-terminal heptad repeat of the gp41, formed a unique hook-like structure (M-T hook) that can dramatically improve the binding stability and anti-HIV activity of the inhibitors. In this work, we explored the molecular mechanism why M-T hook structure could improve the anti-HIV activity of inhibitors. Firstly, molecular dynamic simulation was used to obtain information on the time evolution between gp41 and ligands. Secondly, based on the simulations, molecular mechanics Poisson-Boltzmann surface area (MM-PBSA) and molecular mechanics Generalized Born surface area (MM-GBSA) methods were used to calculate the binding free energies. The binding free energy of the ligand with M-T hook was considerably higher than the other without M-T. Further studies showed that the hydrophobic interactions made the dominant contribution to the binding free energy. The numbers of Hydrogen bonds between gp41 and the ligand with M-T hook structure were more than the other. These findings should provide insights into the inhibition mechanism of the short peptide fusion inhibitors and be useful for the rational design of novel fusion inhibitors in the future. Keywords: Molecular dynamic simulation; MM-PB/GBSA; Binding free energy; HIV-1; Fusion inhibitors
INTRODUCTION Acquired immune deficiency syndrome (AIDS), which is induced by human immunodeficiency virus (HIV) infection, has become one of the major medical and humanitarian challenges(Liu, et al. 2011). So far, the US Food and Drug Administration (FDA) approves four classes of anti-HIV drugs: (i) reverse transcriptase inhibitors (RTIs), (ii) protease inhibitors (PIs), (iii) fusion and entry inhibitors, and (iv) integrase inhibitors (Maguire, et al. 2014). In particular, HIV-1 fusion/entry inhibitors can inhibit early steps of the HIV replication cycle, thereby inhibiting the spread of infection in the first time (Lin and Xiong 2013). Fusion of HIV virus membrane with host cell membrane is an essential event for virus infection and proliferation. This process is mediated by HIV envelope glycoprotein 120 (gp120) (Xu, et al. 2013) and glycoprotein 41 (gp41) (AIhemaiti·abudureyimu, et al. 2013), thus, it is served as an attractive target for designing a novel HIV-1fusion/entry inhibitors(Wang, et al. 2013c). T20 (Enfuvirtide, Fuzeon) is the first HIV fusion inhibitor approved by U.S. FDA. Although T20 is effective (Liang, et al. 2012; Tan, et al. 2013), T20 has some shortcomings that limit its clinical usefulness. Firstly, owing to its nature characteristic, T20 is marketed in an injectable form. Patients are subcutaneously injected twice in every day. Because of the chronic nature of the disease therapy, this dosage form may be one of major problems for the patient's adherence to this drug regimen. Secondly, T20 clinical efficacy is also negatively affected by gp41 resistance mutations that arise owing to continued use. To overcome these disadvantages, a number of strategies have been applied to discover the next generation inhibitors with improved potency and pharmaceutical properties(Wang, et al. 2013b). However, the resulting peptides usually inherit a longer sequence(Gao, et al. 2013), such as SC35EK (35-mer)(Xie, et al. 2012b), SFT (sifuvirtide, 36-mer) (Xie, et al. 2012a), and T2635 (39-mer)(Zhang, et al. 2011). Undoubtedly, the length of peptide closely relates to difficulty and cost of its production, but short peptides (Amino acid residue numbers < 14) lead to a dramatic loss of antiviral potency and binding stability(Wang, et al. 2013a). Therefore, the development of peptide’s fusion inhibitor with highly active antiretroviral efficacy and reasonable residue numbers faces great challenges(Panthong, et al. 2013). Promisingly, He Y.’s group discoved a novel strategy to solve above problem. By studying the inhibit's mechanism of gp41, they found that two key residues,
Met626 and Thr627, which forms a hook-like structure (M-T hook), located upstream of the pocket-binding domain (PBD) of the C-terminal heptad repeat (CHR). The M-T hook can dramatically improve the binding stability of the inhibitors(Yang, et al. 2011a; Yang, et al. 2011b); and the addition of the M-T hook structure to the N terminus of C-peptide could markedly enhance its antiviral activity(Liu, et al. 2013; Yang, et al. 2011a). Based on their former research, Chong H. et al proposed a prototype that the short C-peptide fusion inhibitors could be developed by adding an M-T hook structure. They introduced the methionine and threonine to the N terminus of a panel of C peptides (24-mer or less), and found one pepide, named MT-SC22EK, was highly active against diverse HIV-1 variants, including those T20- and SC29EK-resistant variants. Based on the crystal structures of MT-SC22EK, their paper reported that MT-SC22EK primarily targets on the hydrophobic pocket region on N-heptad repeat (NHR) trimer with its PBD and M-T hook regions, whereas the sites conferring T20-resistance are not bound by the inhibitor. Therefore, the M-T hook structure provides a novel strategy for designing short peptide fusion inhibitors that can overcome the problem of drug resistance. The major problem affected activity of HIV-1 fusion inhibitors is molecular resistance. Owing to point mutations in the GIV motif of NHR (L33 to L45), the inhibit activity of T20 is on the decrease (Liu, et al. 2005). Some studies reported that mutations on NHR of gp41 would induce the compensatory mutation on CHR, which may balance those damages and remedy the function of gp41(Xu, et al. 2005). In our former work, we studied the two kinds of mutants systems of gp41, N43D single mutation and compensatory mutation S138A (yielding N43D/S138A), respectively, and the analyses of structural alteration indicated that the inhibitor C34 displays different binding modes in wild type and two mutants(AIhemaiti·abudureyimu, et al. 2013). In N43D single mutation system, it was founded that the N43D mutation introduced a negative charge on receptor, thus led to repulsion between receptor and ligand, and induced the conformation changes of complex, and led to greater losses in hydrophobic interactions and hydrogen bond interaction. In N43D/S138A double mutation system, we founded CHR have more contact with hydrophobic residues of receptor, induced more van der Waals contributions. This great hydrophobic interaction of double mutant could partially offset the energy loss arising from N43D mutation. The change of binding mode on N43D/S138A could explain why compensatory mutation S138A could partially restores infectivity, lead to increased
infectivity and greater resistance to inhibitors than mutation N43D. In order to explore the role of M-T hook structure in improving the activity of fusion inhibitors, we used molecular dynamics (MD) simulations and molecular mechanics Poisson-Boltzmann/ Generalized Born surface area (MM-PB/GBSA) approach to analyze structure-function relationships between gp41 and two ligands (SC22EK and MT-SC22EK) (figure. 1). The only difference of the two ligands is that MT-SC22EK includes M-T hook, the other without. These results may be useful for designing novel HIV fusion inhibitors (Tan, et al. 2016), or optimizing the lead compounds.
MATERIALS AND METHODS Preparation of protein system The initial structures were downloaded from the RSCB protein data bank (www.rcsb.org) (PDB ID: 3VU5, 3VU6). The amino acid sequences of the two ligands are shown in Figure 1. The receptor was one N-terminal heptad repeat (NHR) chain of gp41. We named the system including NHR and MT-SC22EK as complex 1, the other as complex 2.
Molecular dynamic simulation The MD simulations were carried out for each complex system using NAMD2.9 software and amber ff03 force-field. Firstly, the each system was solvated with TIP3P (Jorgensen, et al. 1983) water box, 10 Å of space between the protein and the simulation box boundaries on each direction. After model generation, each system was neutralized separately with Cl- ions. Then, two successive minimization steps were performed: (i) the system was minimized by 20000 steps using the conjugate gradient with 2.0 kcal/mol restraints on the backbone of the protein; (ii) the system was minimized by 20000 steps using the conjugate gradient without constrained atoms. Next, the systems were gradually heated from 10 K to 300 K using Langevin dynamics in a slow, stepwise fashion (20 K every 100 ps). During heating, no restraint was utilized in the complex. Finally, the MD simulations were performed in the NPT ensemble (a constant temperature of 300K and a pressure of 1 atm) for 50 ns (the time step is 2fs) with all restraints removed. Throughout the MD simulations, a cut-off of 12 Å was used for nonbonded interactions while the SHAKE and Particle Mesh
Ewald (PME) were used with periodic boundary conditions(Ryckaert, et al. 1977).
Binding free energy calculation and energy decomposition The binding free energies were calculated by subtracting the free energies of the unbound receptor and ligand from the free energy of the bound complex, shown in eq (1):(Kollman, et al. 2000; Kuhn and Kollman 2000)
G Gcomplex Greceptor Gligand where
Gcomplex ,
Greceptor
, Gligand
(1)
are the average values of the Gibbs free energy for
the complex, receptor, and ligand, respectively. Binding free energies were calculated using the MMPBSA.py approaches as implemented in AMBER12 and AmberTools13. MMPBSA(Zhang, et al. 2012) computed the binding free energy by using a thermodynamic cycle that combined the molecular mechanical energies with the solvated approaches.
MM MM sol sol GMMPB (GB ) SA Eele EvdW Gnonpol Gpol
(2)
sol sol Where the Gnonpol and G pol are nonpolar and polar terms in solution, respectively.
MM MM Eele and EvdW are the electrostatic and van der Waals contributions to the
molecular mechanics energy, respectively. The polar part expresses the electrostatic contribution to solvation, which is obtained by solving the linear Poisson Boltzmann equation in a continuum model of solvent. On the other hand, the other part accounted for the nonpolar contribution is related linearly to the solvent accessible surface area. The solvation free energies are calculated using an implicit solvent model. Several implicit solvent models are available to calculate the solvation free energies, including Generalized-Bron (GB), Poisson-Boltzmann (PB), and Reference Interaction Site Model (RISM) etc. In this experiment, we calculated the solvation free energies with the PB and GB methods. sol Using the eq. (3) calculated the nonpolar solvation energy, Gnonpol : sol Gnonpol A
(3)
where A is the solvent accessible surface area (SASA), using the MSMS program (Sanner, et al. 1996) with a 1.4 Å radius probe; γ is the surface tension, and β is the offset. The value of constants γ and β are 0.0052 kcal mol−1 Å−2 and 0.92 kcal mol−1,
respectively. The nonpolar solvation term in the MM(GB)PBSA method combines an implicit estimate of the entropic changes associated with the insertion of a solute into the solvent. However, it is obvious that the entropic impact of changes in the configurational freedom of the protein and ligand upon complex formation in vacuo. As a general rule, in the process of protein−ligand binding, the restrictions to the number of conformations causes a reduction in entropy; this contribution is known as the configurational (or conformational) entropy. For some systems, some researchers have found that the calculated binding affinity values is more colse experimental binding affinity values by incorporating a normal mode estimation of the entropic component of the binding free energy alongside MMPB(GB)SA(Men and Guo 2010; Tu 2011). Including this contribution, the final binding affinity estimate ΔGtheor is given by
Gtheor PB (GB ) GMMPB (GB ) SA T S
(4)
where GMMPB (GB ) SA is the MMPBSA or MMGBSA binding free energy, T is the thermodynamic temperature, and −TΔS is the normal mode estimate of the configurational entropy penalty of the binding process of ligand. As implemented within the AMBER 12 package employed in this study, S is calculted as eq. (5): S Stra Srot Svib
(5)
where ΔStra, ΔSrot, and ΔSvib are the contributions to changes in translational, rotational, and vibrational freedoms, respectively.To explore the impact of the binding residue for binding energies, the intermolecular interactions between each residue of the NHR of gp41 and inhibitors were calculated in terms of the pair interaction decomposition of free energy(Li, et al. 2013). The calculations were also performed over MD trajectories using the decomposition energy module in AMBER12 and AmberTools13. The steps of the experiment are shown in Figure 2.
RESULTS AND DISCUSSIONS The stability of the simulation systems The root-mean-square deviation (RMSD) of the protein backbone residues in regard
to the initial structure versus simulation time is shown in Figure 3. The results indicate that the RMSD of the complex1 fluctuates reaching a value of about 3 Å at 10 ns and keep stable throughout the course of the trajectory. The RMSD of the complex 2 fluctuates reaching a value of 4 Å at 10 ns and remains stable from 10 ns to 30 ns. From 30 ns to 35 ns, an increasing trend in RMSD behavior is observed. Afterwards complex 2 seems to relatively stable for the next 15 ns about 5 Å. Overall, the fluctuation of RMSD of the complex 2 is larger than the complex 1, which means that the complex 2 is relatively unstable (figure 3). To know the extent of
variation of individual residues of inhibitor (MT-SC22EK and SC22EK) and receptor respectively, in complex systems, the per-residue root-mean-square fluctuations (RMSFs) of backbone atoms were computed. In the case of receptors, we observed the NO.571-590 residue of the complex 1 shows higher RMSF than that found in complex 2 (Fig. 4). On the other hand, RMSF of the two inhibitors is same. The RMSD dissimilarity matrix about the structure of the complex along the trajectory (from the beginning of the optimization to 10 ns MD simulations).The blue color indicates the correlation; the red color represents the uncorrelation (figure 5). To see this can provide some implications: (1) RMSD of complex 1 reaches 3Å (Figure 5a), and complex 2 reaches 4Å (Figure 5b). This also explains that compared with complex 2, complex 1 is more stable. (2) The similarity of the initial structure and the MD stable structure is weak, which means that the final structures obtained from the stable trajectory that had no relation to the starting structure.
Binding free energy calculation Contributions of the binding free energies are summarized in Table 1. The binding free energy of the complex 1 and 2 is -46.70 kcal/mol -38.17 kcal/mol, respectively. The binding free energy of the complex 1 is considerably higher than the complex 2. In accordance with the energy components of the binding free energy (Table 1), we found that the electrostatic and van der Waals energies contributed to strong favorable sol binding of the receptor with the ligand. The non-polar solvation Energies ( Gnonpol )
have similar small favorable contribution to the binding affinity in the two complexes. The entropy (-TS ) for the two complexes is rather close because the structures of the two ligands are very similar (Figure 6, Table 1).
Ranking aggregation result of the per-residue decomposition list Many key residues contribute to binding free energies that ensure the stability of the protein-ligand systems. The binding free energy was decomposed into per-residue to create a per-residue spectrum (figure 7a, figure 7b). The approach of residue decomposition is tremendously useful to clarify the inhibition mechanism at a molecular level and helpful to locate contribution of individual residue to the ligand-receptor interactions, too. The decomposition of Gtheor values on a MM per-residue basis into contributions from EvdW , the sum of electrostatic energy in MM sol gas phase and polar solvation energy ( Eele + G pol ), and non-polar solvation energy sol ( Gnonpol ). It can be seen that the van der Waals energy contributes mostly to the
binding free energies. To further get insight into the difference contribution between the complex 1 and the complex 2 in the van der Waals energies, we calculated each residue’s contribution to the van der Waals (figure 8). We can see: Ile573, Leu576, Gln577 (in the receptor) and Trp623, Trp626, Tyr633 (in the ligand) have a positive contribution to the van der Waals energy in complex 1; In contrast, the contribution to the binding energy of the residues significantly decrease in complex 2. Overall, the van der Waals interaction of each residues in the two systems has a significantly difference.
The hydrogen bond analysis The hydrogen bonds (H-bonds) that occupied more than 20% are listed in Table 2. It shows the system of complex 1 has more H-bonds than the other system. In the system of complex 1, there are two inter-hydrogen bonds (Trp626 and Thr567, Trp623 and Ile571) and two intra-hydrogens (Thr634 and Ile630, Ser644 and Leu640). However, in the other system, there only two intra-hydrogen bonds (Thr632 and Ile628, Ser642 and Leu638). The M-T structure is conducive to form more H-bonds between receptor and the ligand (Table 2). This also explains the reason why the
Trp623 and Trp624 have a positive contribution to the van der Waals energy in complex 1 (AIhemaiti·abudureyimu, et al. 2013).
The analysis of the contact residues We analyzed the contact residues near the hydrophobic pocket to further explore the differences of the van der Waals contribution in the two systems (figure 9). The pocket is formed by a cluster of hydrophobic residues, including Leu565, Leu566, Leu568, Thr569, Val570, Trp571, Gly572, Ile573, Lys574, Leu576, and Gln577(Barlow, et al. 2012). According to the experiment reported, when the distance between molecules is greater than 8.5Å, the van der Waals forces between the molecules is negligible(Menendez-Arias and Alvarez 2014). We calculated the number of contact residues between the hydrophobic pocket of gp41 and the ligands (Table 3). By comparison, we found that there are more contact residues in the complex 1 than the other (Table 3). From the binding model of the gp41 with the ligand (figure 8), we find that the MT-SC22EK inserts the groove of gp41 completely. However, the C-terminal of the SC22EK tend to be away from the combination (figure 8b). The analysis of the contact residues reveals that the M-T hook can further fortify the interaction by “hooking” the pocket residues tightly; in contrast, this structure do not appear between SC22EK and the NHR of gp41. To a certain extent, it can explain why the van der Waals interactions between MT-SC22EK (with M-T) and NHR of gp41 are a few stronger than that of the other.
CONCLUSION The primary objective of this work focused on exploring the molecular mechanism that M-T hook structure in fusion inhibitors improved the anti-HIV activity of inhibitors. MM-PBSA (GBSA) was used in conjunction with normal mode analysis to successfully account for experimental binding free energies between the gp41 protein and MT-SC22EK (or SC22EK). Average energies were obtained from samples of 100 snapshots extracted from MD simulations. Free energy decomposition on a residual basis indicates differences between MT-SC22EK and SC22EK. We found the M-T hook structure was contributed to the van der Waals energy, which was favorable to the binding energy. The analysis of the H-bond indicated that the system with the M-T
hook structure formed more stable H-bonds than the system without M-T hook structure, and these hydrophobic interactions also played a key role in the binding process. The analysis of the contact residues revealed that the M-T hook could further fortify the interaction by “hooking” the pocket residues tightly; in contrast, this structure do not appear in the system without the M-T hook. This work has provided the detailed molecular mechanism that M-T hook structure in fusion inhibitors. In order to improve lead compounds activity, during drug design, these interactions should be considered carefully. Our research illustrates that the interactions between gp41 and short peptide ligands may lead to possible anti-HIV drug design in the future.
Acknowledgements The authors thank Chinese Natural Science Foundation project (No. 21173014, and 11474013).
References AIhemaiti·abudureyimu, et al., 2013. Clinical Curative Effect of Traditional Chinese Medicine on Immune Function of 20 Cases HIV /AIDS Patients. Chinese Archives of Traditional of Chinese Medicine, 31, 491-493. Barlow, D. J., et al., 2012. In-silico studies in Chinese herbal medicines' research: evaluation of in-silico methodologies and phytochemical data sources, and a review of research to date. Journal of Ethnopharmacology, 140, 526-34. Gao, Guojian, et al., 2013. Research of Chinese herbal compound on promoting HIV immune reconstitution China Journal of Chinese Materia Medica, 38, 2523-2526. Jorgensen, W. L., et al., 1983. Comparison of Simple Potential Functions for Simulating Liquid Water. Journal of Chemical Physics, 79, 926-935. Kollman, P. A., et al., 2000. Calculating structures and free energies of complex molecules: Combining molecular mechanics and continuum models. Accounts of Chemical Research, 33, 889-897. Kuhn, B., and P. A. Kollman, 2000. Binding of a diverse set of ligands to avidin and streptavidin: An accurate quantitative prediction of their relative affinities by a combination of molecular mechanics and continuum solvent models. Journal of Medicinal Chemistry, 43, 3786-3791. Li, R., et al., 2013. Imbalanced network biomarkers for traditional Chinese medicine Syndrome in gastritis patients. Sci Rep, 3, 1543. Liang, Feili, et al., 2012. Effect of traditional Chinese medicine combined with acyclovir on herpes zoster in elderly with HIV/AIDS. Modern Journal of Integrated Traditional Chinese and Western Medicine, 21, 917-918. Lin, Jiang, and Yu Xiong, 2013. Current clinical research and prospect on Chinese medicine for AIDS. China Journal of TCM and Pharmacy, 28, 2070-2073. Liu, S. W., et al., 2005. Different from the HIV fusion inhibitor C34, the anti-HIV drug fuzeon (T-20) inhibits HIV-1 entry by targeting multiple sites in gp41 and gp120. Journal of Biological Chemistry, 280, 11259-11273. Liu, Ying, et al., 2011. Retrospective Analysis of Clinical Literature on Chinese Medicine for AIDS. Journal of Traditional Chinese Medicine
54, 1686.
Liu, Zhen, et al., 2013. Effect of Immune 2 with highly active antiretroviral treatment on immune function of HIV/AIDS patients with poor immune reconstitution China Journal of Chinese Materia Medica, 38, 2458-2462. Maguire, J. J., et al., 2014. The CCR5 chemokine receptor mediates vasoconstriction and stimulates intimal hyperplasia in human vessels in vitro. Cardiovascular Research. Men, Jiuzhang, and Lei Guo 2010 A general introduction to traditional Chinese medicine: Science Press & CRC Press. Menendez-Arias, L., and M. Alvarez, 2014. Antiretroviral therapy and drug resistance in human immunodeficiency virus type 2 infection. Antiviral Research, 102, 70-86. Panthong, K., et al., 2013. Benzene, coumarin and quinolinone derivatives from roots of Citrus hystrix. Phytochemistry, 88, 79-84. Ryckaert, J. P., G. Ciccotti, and H. J. C. Berendsen, 1977. Numerical-Integration of Cartesian Equations of Motion of a System with Constraints - Molecular-Dynamics of N-Alkanes. Journal of Computational Physics, 23, 327-341.
Sanner, M. F., A. J. Olson, and J. C. Spehner, 1996. Reduced surface: An efficient way to compute molecular surfaces. Biopolymers, 38, 305-320. Tan, J. J., et al., 2013. The Current Status and Challenges in the Development of Fusion Inhibitors as Therapeutics for HIV-1 Infection. Current Pharmaceutical Design, 19, 1810-1817. Tan, J., et al., 2016. Design, synthesis and activity evaluation of novel peptide fusion inhibitors targeting HIV-1 gp41. Bioorg Med Chem, 24, 201-206. Tu, Y., 2011. The discovery of artemisinin (qinghaosu) and gifts from Chinese medicine. Natural Medicines, 17, 1217-20. Wang, J., W. Zou, and Y. Liu, 2013a. [Analysis on efficacy action points of HIV/AIDS treatment with traditional Chinese medicine and antiretroviral drugs]. Zhongguo Zhong Yao Za Zhi (China Journal of Chinese Materia Medica), 38, 2504-6. Wang, Jian, et al., 2013b. A sixty-month study of the change of CD4T－lymphocyte cell count treated with traditional Chinese medicine in 565HIV / AIDS patients. Chinese Journal of AIDS & STD, 19, 468-470. Wang, Jie, et al., 2013c. Study of drug intervention in diversity changes of TCRVβ gene in AIDS patients with incomplete immune reconstitution. China Journal of Chinese Materia Medica, 38, 2424-2428. Xie, Zheng, et al., 2012a. The Clinical Research of TCM on HIV/AIDS Depression Traditional Chinese Medicine, 6, 172-173. Xie, Zheng, Li. Pengyu, and Ziqiang Jiang, 2012b. Clinical Observation on 26 Cases of HIV -pruritic Popular Eruption. Chinese Medicine Modern Distance Education of China, 10, 10-12. Xu, L., et al., 2005. Emergence and evolution of enfuvirtide resistance following long-term therapy involves heptad repeat 2 mutations within gp41. Antimicrobial Agents and Chemotherapy, 49, 1113-1119. Xu, Liran, et al., 2013. Study on quality of life of asymptomatic HIV infected persons with traditional Chinese medicine. China Journal of Chinese Materia Medica, 38, 2480-2483. Yang, Liumeng, et al., 2011a. In vitro anti-HIV-1 activities of Qishile,a Chinese medicine effective fraction formula. Chinese Pharmacological Bulletin, 27, 566-571. Yang, Shaochun, et al., 2011b. A Clinical Study on the Treatment of Methadone Maintenance HIVinfected Persons with Fuzheng Ant-i Drug Capsules. Yunnan Journal of Traditional ChineseMedicine and Materia Medica, 32, 18-20. Zhang, Aihua, et al., 2012. Future perspectives of personalized medicine in traditional Chinese medicine: A systems biology approach. Complementary Therapies in Medicine, 20, 93-99. Zhang, Ying, et al., 2011. Chinese Medicine Treatment of HIV Infection. Chinese Archives of Traditional Chinese Medicine, 29, 2628-2630.
Figure Captions Figure 1. Positively charged lysine residues are in blue; negatively charged glutamic acid residues are in red. Black solid lines indicate salt bridges observed in the structure. Green solid line between Thr and Glu indicates hydrogen bond. The sequences of the SC22EK and MT-SC22EK have only one difference is the two residues in the front of MT-SC22EK
Figure 2. General workflow for performing steady state calculated with MMPBSA.py. tleap is a program in Amber used to create topology files for dynamics. “Dry” topologies and ensembles are systems without explicit solvent that is subsequently treated using an implicit solvent model (Li, et al. 2013)
Figure 3. Backbone RMSD are shown as a function of time for complex 1 and complex 2 using the ptraj module in AMBER10
Figure 4. Comparisons of backbone atoms RMSF in inhibitor SC22EK and MT-SC22EK and gp41 NHR
Figure 5. RMSD dissimilarity plot for complex 1 (a) and complex 2 (b). The horizontal and vertical coordinate of the figure indicate the label of simulated trajectories
Figure 6. Energy components (kcal/mol) for the binding of complex 1 and complex 2
Figure 7. The decomposition of G on a per-residue basis into contributions from the MM van der Waals energy ( EvdW ), the sum of electrostatic interactions and polar
MM sol sol solvation energy ( Eele + G pol ), and nonpolar solvation energy ( Gnonpol ) for the
residues of G
MM Figure 8. Energy shown as the contribution from the van der Waals energy ( EvdW )
Figure 9. The binding model of the gp41 with the ligand. (a) The model of complex 1; (b) The model of complex 2. The structural composites from the average structure of the 20-50 ns molecular dynamic simulation trajectories.
Tables Table1. Binding free energies for two complexes by using the MM-PB (GB) SA method (800 snapshots from 20 ns-50 ns)
Component
Complex 1
Complex 2
MM Eele
-56.15 (0.09)*
MM EvdW
-62.47 (0.48)*
-48.09 (0.08)*
sol G pol ( GB )
85.19 (0.47)*
171.87 (0.56)*
sol G pol ( PB )
79.16 (0.47)*
165.08 (0.57)*
sol Gnonpol ( GB )
-8.35 (0.01)*
-7.33 (0.01)*
sol Gnonpol ( PB )
-56.15 (0.09)*
-TS
-149.39
-149.39
(0.62)*
(0.62)*
33.26 (0.27)*
32.97 (0.38)*
Gtheor GB
-41.74
(0.16)*
-32.99 (0.25)*
Gtheor PB
-46.70
(0.16)*
-38.17 (0.24)*
Gexpt(Liu, et al. 2013)
-11.71
All values are given in kcal/mol * standard error of mean values
-9.53
Table2. Hydrogen-bond analysis from the results of MD simulation System
Donor
Accepter
Thr639-HG1-OG1
Ile635-O
Ser649-HGI-OG
Distance
(Å)
Angle (Degrees)
Occupied（%）
2.72 ( 0.11)*
164.61 (8.63)*
97.12
Leu645-O
2.74 ( 0.11)*
162.42 (10.32)*
74.75
Trp631-HE1-NE1
Thr569-O
2.77 ( 0.17)*
160.00 (10.64)*
43.83
Trp628-HE1-NE1
Ile573-O
2.87 ( 0.20)*
153.85(12.13)*
21.83
Thr639-HG1-OG1
Ile635-O
2.73 ( 0.11)*
164.27 (8.85)*
99.42
Ser649-HG-OG
Leu645-O
2.74 ( 0.11)*
160.98 (9.94)*
87.83
Complex 1
Complex 2
The H-bonds are determined by the donor…acceptor atom distance of ≤3.5 Å and acceptor…H donor angle of ≥120° * standard error
Table 3. The hydrophobic interactions between the hydrophobic pocket and the ligand
Residue
MT-SC22EK
SC22EK
Leu565
8*
6*
Leu566
14*
13*
Leu568
6*
3*
Thr569
13*
7*
Val570
13*
10*
Trp571
4*
5*
Gly572
8*
3*
Ile573
11*
9*
Lys574
5*
6*
Leu576
10*
2*
Gln577
8*
7*
The hydrophobic interactions with distance of less than 8.5Å. * The number of hydrophobic contact
