INTERACTION MECHANISMS OF NOVEL ANTIVIRAL PEPTIDES AGAINST N1 SUBTYPE NEURAMINIDASE: A COMPUTATIONAL EXPLORATION
F. WUa, b, Z. YANGa,b*, S. YINc, C. SONGd, X. YUANe
a School of Basic Medical Sciences, Jiamusi University, Jiamusi 154007, China
bEngineering Research Center of Forest Bio-preparation, Ministry of Education, Northeast Forestry University, Harbin 150040, China
cThe Central Hospital of Jiamusi City, Jiamusi 154002, China
dHeilongjiang Nursing College, Harbin 150086, China
eInstitute of Biomedicine, Jinan University, Guangzhou 510632, China
In light of the current pandemic threat and the emergence of drug resistances, the development of next-generation anti-influenza agents must be a high priority, with many efforts to bind with the 150-cavity of N1 subtype neuraminidase (NA). In this work, docking and explicit solvent molecular dynamics simulations were combined to study the interactions between N1 subtype NA and novel inhibitory peptides EISYIHAEAYRRGELK, YIHAEAYRRG and HAEAYR. It was found that 16-mer peptide EISYIHAEAYRRGELK has the binding specificity to N1 subtype NA and the capability to lock the 150-cavity. The binding-site residues Asp198, Glu227 and Glu229 play important roles during the binding process. Compared with the 16-mer peptide, YIHAEAYRRG and HAEAYR even tightly coordinated to the 150-cavity, especially of the latter. The total interaction energies (Etotal) of them with N1 subtype NA were calculated to be -578.21, -442.21 and -358.80 kcal mol-1, respectively. Further energetic and geometric analyses revealed that the portion of Tyr-Arg in peptides is more helpful for maintaining favorable contacts with residues around 150-cavity and is of great importance in the scaffold modification. We hope that the results will be helpful for designing novel anti-influenza drugs.
(Received May 27, 2014; Accepted July 18, 2014)
Keywords: Neuraminidase; Antiviral peptides; Interaction mechanisms; Docking;
Molecular dynamics
1. Introduction
The influenza virus, a major cause of respiratory disease, constantly introduces possible public health threat [1]. In the early years of this century, highly pathogenic A/H5N1 and A/H1N1 (2009) influenza viruses have emerged and spread quickly around the world [2]. Beyond the annual production of vaccines, the mainstay of strategies to prevent and treat the virus infections is stockpiling of antiviral drugs, most commonly the neuraminidase (NA) inhibitors zanamivir (Relenza), oseltamivir (Tamiflu) and peramivir [3]. The frequent emergences of resistant strains to these agents have spurred the development of novel potential therapies [4-5].
Influenza NA can specifically identify sialic acid (SA) species and assist the cleavage of sialic acid residues from infected cell surfaces [1]. Blocking its function will interrupt the release of nascent virions, reducing the spread of the infection to neighboring cells [6-7]. NA has thereby been proposed as an attractive target for the design of anti-influenza drugs [7-8]. Currently, there are two phylogenetically distinct groups: Group-1 NA comprises N1, N4, N5 and N8, whereas the N2, N3, N6, N7 and N9 subtypes belong to Group-2 [7]. Group-1 differs from Group-2 by having the “150-cavity”, consisting of residues 147-152 (the 150-loop), adjacent to the binding site which
* Corresponding authors: [email protected]
generally composed of largely conserved residues (Fig. 1) [9-10].
Fig. 1. A) Ribbon schematic representation of N1 protein in complex with oseltamivir carboxylate (balls and sticks). B) Close up view of the highly flexible and charged binding-site residues. Protein coordinates from PDB entry 2HU0. The colors of the ribbons distinguish between helices (red), β-sheets (cyan), hydrogen-bonded turns (green) and random coils (white), and the colors of atoms between oxygen (red), nitrogen (blue) and carbon (green).
Most known NA inhibitors are mainly developed by referring to the binding site domain from Group-2 subtypes N2 and N9, such as zanamivir and oseltamivir [6-8]. In the case of oseltamivir binding (Fig. 1B), a triad of arginines (Arg118, Arg292 and Arg371) binds to the carboxylate anion; a negatively charged region (Glu119, Asp151 and Glu227) interacts with the amino group; residue Arg152 forms a hydrogen bond with the acetylamino group; residue Glu276 has hydrophobic stacking with the ethylpropoxy group. However, several studies have revealed the 150-cavity (open state) appears to be closed by the movement of 150-loop in response to ligand binding (closed state) [9, 11]. This motion may not be restricted to Group-1, and oseltamivir can induce partially opening of the N2 150-loop [12-14]. It indicates that the flexibility of 150-loop should be considered in the drug design, to improve binding affinity and specificity [12, 14]. For instance, numerous efforts have been devoted to the development of NA inhibitors that projected into the 150-cavity [10, 15-18].
The peptides against influenza viruses have shown potential as therapeutic agents [19-22].
Recently, a 16-mer peptide derived from Acetes Chinensis (EISYIHAEAYRRGELK, N- to C-terminus) has been reported to stand against influenza A/H1N1 virus, with an IC50 of 96.1 μmol/L [22]. Albeit its binding with NA has been indicated by evidences of UV spectrum analysis and docking/scoring method [22], however, the interaction mechanisms remain ambiguous, which is adverse for the design of improved inhibitors [3, 6]. Molecular dynamics (MD) simulations have been extensively used to expound the features of interactions between inhibitors and the NA binding site [3, 23-27]. Here, docking and explicit solvent MD simulations were combined to analyze the interactions involving 16-mer peptide with N1 subtype NA. Besides, the bindings of 10-mer (YIHAEAYRRG) and 6-mer (HAEAYR) peptides to N1 subtype NA were also studied to reveal the importance of individual functional groups of peptides to the binding affinities. We hope that the results can aid to understand the interaction profiles, and will be of value in the rational design of improved inhibitors across all NA subtypes.
150-cavit y
A B
2. Materials and methods 2.1. System setup
Because of high-resolution and well-defined binding site, the coordinates of H5N1 NA complexed with oseltamivir carboxylate (entry code 2HU0 [9]) were retrieved from Protein Data Bank. For convenience, it is named as N1 throughout this work. The calcium ion (Ca2+) near the binding site was retained. Missing hydrogen atoms of the protein were added, using Biopolymer module (InsightII 2005) [28]. Protonation states for the residues were determined by manual verification according to the interacting partners [29-31]. The protein structure was then neutralized with Na+ ions [29-31], and optimized with the conjugated gradient algorithm (Discover 3.0 module), using the consistent-valence force-field (CVFF). The convergence criterion was set to 0.01 kcal mol-1 Å-1.
The structures of peptides EISYIHAEAYRRGELK, YIHAEAYRRG and HAEAYR were generated with Builder module. The geometries and partial atomic charges of peptides were handled by the BFGS algorithm [32], with a convergence criterion of 0.01 kcal mol-1 Å-1.
2.2. Docking
The docking simulations were performed by the general and popular protocol in InsightII 2005 software packages [33], features for its semi-flexible method that the binding-site specified residues and ligands are free to move. The Affinity module, combining Monte Carlo (MC) and simulated annealing (SA) methods, was used to probe the optimal orientations of peptides at the binding site of N1. The solvent effect was accounted by solvating the systems in a sphere of TIP3P [34] water molecules. The non-bonded interactions were described by Cell-Multipole approach.
Energy minimizations were then performed, using the conjugated gradient method with CVFF force-field, until converged to 0.01 kcal mol-1 Å-1. More calculated details are referred elsewhere [31, 35-36].
2.3. Molecular dynamics (MD)
The energy-minimized docked complexes were used as the initial structures during the MD simulations, using GROMACS4.5.5 program [37-38] and Charmm27 force field [39-40].
Each complex was placed in a SPC/E (simple-point-charge) water box of 81×81×81 Å3 [31]. NA+ counter-anions were added to neutralize the charge of the system [41]. 70 ns production MD simulations were performed, after removing bad contacts with molecular mechanics and relaxing water molecules with position-restrained MD. The NPT ensemble was applied at 300 K and 1 Bar [42]. The particle-mesh Ewald (PME) method was used with 8.0 Å direct-space nonbonded cutoff [43]. The covalent bonds involving hydrogen was constrained with the LINCS algorithm [44]. A 1.0 fs time step was used, and coordinates were saved every 10 ps.
3. Results and discussion
The structural and energetic analyses were conducted over the structures obtained from explicit solvent MD simulations. The structural stability of the dynamic trajectories was monitored through time evolutions of potential energies and root-mean-square deviations (RMSD), which is widely used to determine the equilibrium state [45]. Fig. 2 shows that the potential energies of docked complexes are well-behaved during the simulations; the protein backbone RMSD curve reaches the platform since approximately 0.5 ns, with only minor variation afterwards, consistent with the previous results [29-31]. The average protein backbone RMSD for the 16-mer, 10-mer and 6-mer peptide-N1 complexes are 1.8 (0.1), 2.2 (0.1) and 2.2 (0.2) Å, respectively. The data in parentheses are for the standard deviations. Relatively large fluctuations are observed for the ligand positional RMSD until about 60 ns, with the average values of 3.2 (0.6), 2.5 (0.2) and 3.8 (0.8) Å, respectively (Fig. 2). These indicated that the simulations are reasonably converged and
docked complexes have readily been in equilibrium since about 60 ns. Thus, the geometric and energetic analyses are made on the average structures of 60~70 ns MD trajectories. It was found that all of three peptides are in the proximity space within N1, but their binding poses differ somewhat (Fig. 3).
Fig. 2. Time evolution of potential energies and root-mean-square deviations (RMSD) for the peptide-N1 complexes. Backbone-atom (main chain atoms Cα, C, N, O) RMSD plots are represented by solid lines, and those of peptide heavy atoms are represented by dot lines.
3.1. N1 interacting with 16-mer peptide
Albeit virologic and UV spectrum experiments have foreshadowed the effectiveness of 16-mer peptide (EISYIHAEAYRRGELK) against N1 protein [22], however, the interactions still require further exploration. On the time scale of the 70 ns simulations, 16-mer peptide is able to enter and fill the cavity of N1 drug binding site and 150-cavity, which further demonstrated to prevent closure of the 150-loop (Fig. 3). The carboxylate anion of E1 in 16-mer peptide was docked towards guanidino group of residue Arg430, with one H-bond formed (Fig. 4A). Note that the residues of peptides all match 1-letter residue name, while residues of N1 are named with 3-letter one. E1 was further stabilized by the main-chain amides of residues Arg371 and Gly429 via the H-bonding interactions. The main-chain amide of I5 was oriented towards residue Gly348 along with the formation of one H-bond. The imidazole functional group of H6 was stabilized by residue Arg300. The carbonyl group of A7 was docked towards residue Arg292, with one H-bond formed. There is a sandwich-like H-bonding interaction involving E8 and the polar parts of residues Ser279 and Lys350 (Fig. 4A). R11 was observed to engage the H-bond networks with residues Glu227, Glu229 and Glu277. The guanidino group of R12 was docked towards residues Ser179 and Ala180, with one and one H-bonds. It was found that L15 has one H-bond with residue Trp178, and a bidentate H-bond interaction is formed between K16 amino group and the main-chain carbonyl groups of residues Pro197 and Gly200. Detail information is shown in Table 1 and Fig. 4A.
0 10000 20000 30000 40000 50000 60000 70000 -180000
-174000 -168000 -162000 -156000
MD simulation time (ps) Potential energy (kcal mol-1 )
0 10000 20000 30000 40000 50000 60000 70000 0
1 2 3 4 5 6
RMSD (Å)
MD simulation time (ps)
16-mer peptide-N1 10-mer peptide-N1 6-mer peptide-N1 OC-NA PR-NA
Fig. 3. Peptides superposed at the N1 binding site. Protein is shown in property Connolly surface, using the InsightII 2005 scripts. Peptides are represented by stick models.
16-mer peptide 10-mer peptide 6-mer peptide OC-NA PR-NA
Fig. 4. Views of the binding modes of the N1 binding-site residues with the peptides: A) EISYIHAEAYRRGELK, B) YIHAEAYRRG and C) HAEAYR Key residues and peptides are represented by stick, and ball and stick models, respectively. The C, N, O atoms are colored in purple, blue and red for the binding-site residues whereas in green, blue and red for the peptides.
The total interaction energy (Etotal) of 16-mer peptide with N1 was summed to -572.81 kcal mol-1. The electrostatic interactions (Eele) are the main driving force for the binding process, with the value of -385.17 kcal mol-1 (67.24 %), which is consistent with the previous studies of other NA inhibitors [29, 35-36, 46-47]. It was found that binding-site residues Asp198, Glu227 and Glu229 show their primary effects on the electrostatic contributions (Eele), with the values of -50.74, -66.54 and -53.49 kcal mol-1 (Fig. 5A), in agreement with the structural analysis. Besides, the hydrophobic portions of residues Glu227, Glu277, Tyr347 and Lys350 were observed to introduce steadily van der Waals effects (EvdW) to the binding , where the values were -9.74, -7.96, -10.05 and -7.12 kcal mol-1, respectively (Fig. 5A).
A
B C
3.2. N1 interacting with 10-mer/6-mer peptide
As described in above section, the relative stability of 16-mer peptide-N1 complex observed in the 70 ns MD simulations correlates well with the experiments [22], while the 16-mer peptide seems to slightly spin around, with limited compactness to N1 (Figs. 3 and 4A); moreover, it is not facile to synthesize and commercialize, because of the relatively large molecular size. On basis of the characteristics of NA binding site [6, 9, 18, 31, 48], 10-mer peptide YIHAEAYRRG was first used as probe to explore the key portions that make close contact with residues of binding site and 150-cavity. That is, respectively remove three residues from the N- and C-terminuses of 16-mer peptide, and the simplest amino acid Gly is at the C-terminus of the produced 10-mer peptide, which is expected to orient towards the Arg triad (Arg118, Arg292 and Arg371) of the N1 binding site.
Table 1 and Fig. 4B show that the carboxylate anion of G10 in 10-mer peptide was docked towards residues Arg371 and Tyr406, with the formation of two and one H-bonds, respectively. It is encouraging that carboxylate group has been a key feature of all reported NA inhibitors, which forms strong H-bonds with the trio of Arg residues [29, 35-36, 46-47]. The R9 guanidino group was stabilized by the carboxylate group of residue Glu425 through three H-bonds. It was found that the side-chain of Y7 forms one H-bond with residue Glu119, Glu119 also interacts with Y1, which in turn interacts with residues Ser179 and Glu227, associated with numerous H-bonds. The binding pose of 10-mer peptide within N1 protein is rather close to the 150-cavity, in contrast to the case of 16-mer peptide (Fig. 3). The interaction energy (Etotal) of 10-mer peptide and N1 was calculated at -442.21 kcal mol-1, and the electrostatic interactions (Eele) play a larger role, contributing to 84.05 % (Table 1). The calculations predicted that residues Glu119, Glu227, Glu277 and Glu425 have more contributions of electrostatic effects (Eele) with 10-mer peptide, where the values are -122.76, -112.14, -49.49 and -113.69 kcal mol-1 (Fig. 5B). At the same time, residues Tyr 406 and Trp438 were observed to result in strong van der Waals interactions (Evdw), with the value of -7.43 and -7.74 kcal mol-1 (Fig. 5B).
Table 1. The vdW energies (Evdw), electrostatic energies (Eele), total interaction energies (Etotal) and H-bonding information in the peptide-N1 complexes a
Peptide Evdw Eele Etotal
The binding-site residues that form H-bonding with peptide
b
16-mer peptide EISYIHAEAYRRGE LK
-193.04 -385.17 -578.21
Trp178 (1), Ser179 (1), Ala180 (1), Pro197 (1), Gly200 (1), Glu227 (1), Glu229 (1), Glu277 (1), Ser279 (1), Arg292 (1), Arg300 (1), Gly348 (1), Lys350 (1), Arg371 (1) Gly429 (1), Arg430 (1)
10-mer peptide
YIHAEAYRRG -70.53 -371.68 -442.21
Glu119 (2), Ser179 (1), Glu227 (1), Arg371 (2), Tyr406 (1), Glu425 (3)
6-mer peptide
HAEAYR -55.36 -303.44 -358.80
Glu119 (1), Glu227 (1), Glu277 (1), Arg371 (2), Arg430 (1), Tyr438 (1)
a Energy units in kcal mol-1;
b The numbers of H-bonds are given in parentheses.
Fig. 5. The electrostatic (Eele, blue sparse area) and total interaction energies (Etotal) between per binding-site residue of N1 and peptides: A) EISYIHAEAYRRGELK, B)
YIHAEAYRRG and C) HAEAYR
Regarding as 6-mer peptide HAEAYR, its binding pose is characterized by the strong electrostatic interactions involving the carboxylate anion of E3 with the positively charged side chains of binding-site residues Arg371 and Arg430, with the formation of two and one H-bonds (Table 1 and Fig. 4C). The R6 guanidino group was docked towards residue Trp438, with one H-bond formed. It was found that H1 and A2 were stabilized by residues Glu227, Glu277 and residue Glu119, with one H-bond formed with each residue. In addition, Y5 has hydrophobic stacking with the 150-cavity (Figs. 3 and 4C). The binding affinity (Etotal) of 6-mer peptide with N1 is slightly weaken (-358.80 kcal mol-1), in contrast to the above two modes. In 6-mer peptide-N1, the electrostatic interactions (Eele) dominate the binding process with the value (proportion) of -303.44 kcal mol-1 (84.57 %). It was further found that electrostatic effects (Eele) are more helpful for the interactions of residues Glu119, Glu227, Glu277 and Arg371with 6-mer peptide, especially residue Glu227, where the value was -123.74 kcal mol-1 (Fig. 5C). While, residues Tyr 406 and Trp438 show their primary effects on the van der Waals interactions (Evdw), with the values of -6.49 and -5.27 kcal mol-1 (Fig. 5C).
Glu11 9 Trp17
8 Ser17
9 Ala18
0 Ile1
92 Ile1
94 Gly19
6 Pro19
7 Asp19
8 Gly20
0 Ile2
22 Thr22
5 Gln22
6 Glu22
7 Ser22
8 Glu22
9 Glu27
6 Glu27
7 Cys2
78 Ser27
9 Val29
0 Arg29
2 Asp29
3 Arg30
0 Tyr34
7 Gly34
8 Val34
9 Lys35
0 Phe35
2 Arg37
1 Tyr40
6 Gly42
9 Arg43
0 Pro43
1 0
-20 -40 -60 -80
Energy (kcal mol-1 )
A
Glu119 Thr148 Asp151 Trp178 Ser179 Ala180 Leu223 Thr225 Glu227 Glu229 Glu277 Arg371 Tyr406 Glu425 Ile437 Trp438 0
-30 -60 -90 -120
Energy (kcal mol-1)
B C
Glu119 Gln136 Thr148 Trp178 Glu227 Glu277 Arg371 Tyr406 Arg430 Trp438 Thr439 0
-30 -60 -90 -120
Energy (kcal mol-1)
Fig. 6. Time evolution of the H-bonding number of the peptide-N1 complexes.
Our MD simulation results indicated that 16-mer peptide EISYIHAEAYRRGELK is able to fit into both of the N1 binding site and 150-cavity and can lock the 150-loop in an open conformation (Figs. 3 and 4A). Binding-site residues Asp198, Glu227 and Glu229 seem to be important during the binding process (Fig. 5). After removing residues from the N- and C-terminuses, the generated 10-mer and 6-mer peptides still targeted at the 150-cavity, even tightly coordinated to the cavity (Fig. 3). Despite reducing in binding affinity (Etotal), the stability (potential energies) seems to be enhanced, with H-bond population remaining almost constant across the trajectory (Figs. 2 and 6). Especially 6-mer peptide HAEAYR is found to interact stably with the binding site and 150-cavity (Figs. 3 and 6). Overall the three peptides, the portion of Tyr-Arg is more helpful for the bindings, as a result of maintaining favorable contacts with residues around the 150-loop (Fig. 4). This can be comparable to the previous drug design that the long hydrophobic side group of 3-(p-tolyl) allyl-Neu5Ac2en points to the 150-cavity [10].
4. Conclusions
Influenza neuraminidase (NA) is an attractive target for drug design and the 150-cavity of Group-1 NA provides new opportunity to develop dual-site-binding inhibitors. Recently, a novel NA inhibitory peptide (EISYIHAEAYRRGELK) was purified and observed to be active against influenza A/H1N1 virus. In order to explore its interaction mechanisms, docking and molecular dynamics (MD) methods were combined to study the binding modes of N1 protein with peptides EISYIHAEAYRRGELK, YIHAEAYRRG and HAEAYR.
0 10000 20000 30000 40000 50000 60000 70000 0
6 12 18 24
H-bon ding Number
MD simulation time (ps)
16-mer peptide-N1 10-mer peptide-N1 6-mer peptide-N1 OC-NA PR-NA
It was found that 16-mer peptide EISYIHAEAYRRGELK has the binding specificity to N1 and the capability to lock the 150-cavity. The binding-site residues Asp198, Glu227 and Glu229 can form strong interactions with the 16-mer peptide, and more attention should be paid to explore novel agents. 10-mer and 6-mer peptides have also been shown to interact stably with N1 and tightly coordinate to the 150-cavity, especially of the latter. The total interaction energies (Etotal) of the docked complexes of N1 with EISYIHAEAYRRGELK, YIHAEAYRRG, HAEAYR are -578.21, -442.21 and -358.80 kcal mol-1. According to the energetic and geometric analyses, the portion of Tyr-Arg in peptides is more helpful for the bindings and is of great importance in the scaffold modification. We hope that the results will be helpful for designing novel inhibitors to combat the spread of influenza virus.
Acknowledgement
We are grateful for the financial supports from National Natural Science Foundation (No.
31200429) and Scientific Research Foundation of Heilongjiang Provincial Education Department (No. 12531690).
References
[1] G. M. Air, Influenza Other Resp 6, 245 (2012).
[2] S. P. Layne, A. S. Monto, J. K. Taubenberger, Science 323, 1560-1561 (2009).
[3] P. M. Colman, Annu Rev Biochem 78, 95 (2009).
[4] P. J. Collins, L. F. Haire, Y. P. Lin, J. Liu, R. J. Russell, P. A. Walker, J. J. Skehel, S. R. Martin, A. J. Hay, S. J. Gamblin, Nature 453, 1258 (2008).
[5] I. Stephenson, J. Democratis, A. Lackenby, T. McNally, J. Smith, M. Pareek, J. Ellis, A. Bermingham, K. Nicholson, M. Zambon, Clin Infect Dis 48, 389 (2009).
[6] M. von Itzstein, Nat. Rev. Drug. Discov. 6, 967 (2007).
[7] K. Das, J. M. Aramini, L. C. Ma, R. M. Krug, E. Arnold, Nat Struct Mol Biol 17, 530 (2010).
[8] A. Moscona, Annu. Rev. Med. 59, 397 (2008).
[9] R. J. Russell, L. F. Haire, D. J. Stevens, P. J. Collins, Y. P. Lin, G. M. Blackburn, A. J. Hay, S. J. Gamblin, J. J. Skehel, Nature 443, 45 (2006).
[10] S. Rudrawar, J. C. Dyason, M. A. Rameix-Welti, F. J. Rose, P. S. Kerry, R. J. Russell, S. van der Werf, R. J. Thomson, N. Naffakh, M. von Itzstein, Nature communications 1, 113 (2010).
[11] R. E. Amaro, D. D. L. Minh, L. S. Cheng, W. M. Lindstrom, A. J. Olson, J.-H. Lin, W. W. Li, J. A. McCammon, J. Am. Chem. Soc. 129, 7764 (2007).
[12] R. E. Amaro, X. Cheng, I. Ivanov, D. Xu, J. A. McCammon, J. Am. Chem. Soc.
131, 4702 (2009).
[13] R. E. Amaro, R. V. Swift, L. Votapka, W. W. Li, R. C. Walker, R. M. Bush, Nature communications 2, 388 (2011).
[14] Y. Wu, G. Qin, F. Gao, Y. Liu, C. J. Vavricka, J. Qi, H. Jiang, K. Yu, G. F. Gao, Scientific reports 3, 1551 (2013).
[15] L. S. Cheng, R. E. Amaro, D. Xu, W. W. Li, P. W. Arzberger, J. A. McCammon, J. Med.
Chem. 51, 3878 (2008).
[16] A. T. Garcia-Sosa, S. Sild, U. Maran, J Chem Inf Model 48, 2074 (2008).
[17] Y. Li, B. Zhou, R. Wang, J Mol Graph Model 28, 203 (2009).
[18] P. S. Kerry, S. Mohan, R. J. Russell, N. Bance, M. Niikura, B. M. Pinto, Scientific reports 3, 2871 (2013).
[19] T. Matsubara, M. Sumi, H. Kubota, T. Taki, Y. Okahata, T. Sato, J Med Chem 52, 4247 (2009).
[20] M. Rajik, F. Jahanshiri, A. R. Omar, A. Ideris, S. S. Hassan, K. Yusoff, Virol J 6, doi: 10.1186/1743-422X-6-74 (2009).
[21] Z. Yang, G. Yang, Y. Zu, Y. Fu, L. Zhou, Int. J. Mol. Sci. 11, 4932 (2010).
[22] WANG Hai-Tao, WANG Wei, CHEN Ming, YUAN Ning, ZHAO Yuan-Hui, M. Xiang-Zhao, Chem J Chinese U 34, 2540 (2013).
[23] T. Udommaneethanakit, T. Rungrotmongkol, U. Bren, V. Frecer, M. Stanislav, J Chem Inf Model 49, 2323 (2009).
[24] P. Wang, J. Z. Zhang, J Phys Chem B 114, 12958 (2010).
[25] P. Kar, V. Knecht, J Phys Chem B 116, 6137 (2012).
[26] L. Li, Y. Li, L. Zhang, T. Hou, J Chem Inf Model 52, 2715 (2012).
[27] A. Vergara-Jaque, H. Poblete, E. H. Lee, K. Schulten, F. Gonzalez-Nilo, C. Chipot, J Chem Inf Model 52, 2650 (2012).
[28] InisghtII Version 2005. 2005 ed.; Accelrys Inc.: San Diego, USA., 2005.
[29] Z. W. Yang, X. M. Wu, Y. G. Zu, G. Yang, L. J. ZHou, Int. J. Quantum. Chem.
112, 909 (2012).
[30] Z. Yang, Y. Yang, F. Wu, X. Feng, Mol. Simulat. 39, 788-795 (2013).
[31] Z. Yang, G. Yang, L. Zhou, J Comput Aided Mol Des 27, 935-950 (2013).
[32] J. D. Head, M. C. Zerner, Chem. Phys. lett. 122, 264-270 (1985).
[33] Affinity User Guide. Accelrys Inc.: San Diego, USA, 2005.
[34] W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey, M. L. Klein, J. Chem. Phys.
79, 926 (1983).
[35] Z. W. Yang, G. Yang, Y. G. Zu, Y. J. Fu, L. J. Zhou, Phys. Chem. Chem. Phys.
11, 10035 (2009).
[36] Z. Yang, Y. Nie, G. Yang, Y. Zu, Y. Fu, L. Zhou, J. Theor. Biol. 267, 363 (2010).
[37] B. Hess, C. Kutzner, D. van der Spoel, E. Lindahl, J. Chem. Theory Comput.
4, 435 (2008).
[38] S. Pronk, S. Pall, R. Schulz, P. Larsson, P. Bjelkmar, R. Apostolov, M. R. Shirts, J. C. Smith, P. M. Kasson, D. van der Spoel, B. Hess, E. Lindahl, Bioinformatics 29, 845 (2013).
[39] B. R. Brooks, C. L. Brooks, 3rd, A. D. Mackerell, Jr., L. Nilsson, R. J. Petrella, B. Roux, Y. Won, G. Archontis, C. Bartels, S. Boresch, A. Caflisch, L. Caves, Q. Cui, A. R. Dinner, M. Feig, S. Fischer, J. Gao, M. Hodoscek, W. Im, K. Kuczera, T. Lazaridis, J. Ma, V. Ovchinnikov, E. Paci, R. W. Pastor, C. B. Post, J. Z. Pu, M. Schaefer, B. Tidor,
R. M. Venable, H. L. Woodcock, X. Wu, W. Yang, D. M. York, M. Karplus, J Comput Chem 30, 1545 (2009).
[40] P. r. Bjelkmar, P. Larsson, M. A. Cuendet, B. Hess, E. Lindahl, J Chem Theory Comput 6, 459 (2010).
[41] H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, J. Hermans, In: Pullman B, Ed. (Intermolecular forces. Dordrecht: Reidel, 331-342 (1981).
[42] H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, A. DiNola, J. R. Haak, J Chem Phys 81, 3684 (1984).
[43] T. Darden, D. York, L. Pedersen, J Chem Phys 98, 10089 (1993).
[44] B. Hess, H. Bekker, H. J. C. Berendsen, J. G. E. M. Fraaije, J Comput Chem 18, 1463 (1997).
[45] A. Grossfield, D. M. Zuckerman, Annual reports in computational chemistry 5, 23 (2009).
[46] K. M. Masukawa, P. A. Kollman, I. D. Kuntz, J. Med. Chem. 46, 5628 (2003).
[47] N. X. Wang, J. J. Zheng, Protein Sci 18, 707 (2009).
[48] V. Stoll, K. D. Stewart, C. J. Maring, S. Muchmore, V. Giranda, Y.-g. Y. Gu, G. Wang, Y. Chen, M. Sun, C. Zhao, A. L. Kennedy, D. L. Madigan, Y. Xu, A. Saldivar, W. Kati, G. Laver, T. Sowin, H. L. Sham, J. Greer, D. Kempf, Biochemistry 42, 718 (2003).
