SOFTWARE NEWS AND UPDATES
WWW.C-CHEM.ORG
Introduction of Periodic Boundary Conditions into UNRES Force Field Adam K. Sieradzan* In this article, implementation of periodic boundary conditions (PBC) into physics-based coarse-grained UNited RESidue (UNRES) force field is presented, which replaces droplet-like restraints previously used. Droplet-like restraints are necessary to keep multichain systems together and prevent them from dissolving to infinitely low concentration. As an alternative for droplet-like restrains cuboid PBCs with imaging of the molecules were introduced. Owing to this modification, artificial forces which arose from restraints keeping a droplet together were eliminated what leads to more realistic trajectories. Due
to computational reasons cutoff and smoothing functions were introduced on the long range interactions. The UNRES force field with PBC was tested by performing microcanonical simulations. Moreover, to asses the behavior of the thermostat in PBCs Langevin and Berendsen thermostats were studied. The influence of PBCs on association pattern was compared with droplet-like restraints on the bba hetero tetramer 1 proC 2015 Wiley Periodicals, Inc. tein system. V
Introduction
molecular.[34,35] Moreover, there are also studies to treat system as a nonperiodic for instance continuum-molecular system. An interesting study about nonperiodic systems has been presented in Refs [36,37]. In this work, a cuboid PBC with cutoff on long range interactions has been implemented to run coarse-grained molecular simulations with the UNRES force field.
Protein–protein interactions play an extremely important role in living organisms and understanding of these interactions is a necessary condition to understand how the organisms are functioning.[1] Amyloid growth[2,3] and their intermediates behavior,[2] enzyme-inhibitor interactions,[4] and modulatory functions of proteins[5] are only a few examples showing the importance of protein–protein interactions. Theoretical studies of such interactions have provided valuable information and explanation of experimental results.[6] Moreover, there has been tremendous advance in all-atom simulation recently,[7–11] allowing to simulate small proteins even in ms time-scale.[11] However, for large systems coarse-graining is still required. The UNited RESidue (UNRES)[12–28] force field is a carefully derived physic-based force-field which has proven to be a powerful tool in protein structure prediction.[29] The previous version of UNRES force field was extended to simulate oligomeric proteins.[3,30–32] To keep the monomers together, the restrains C-terminus of chain i and N-terminus of chain i 1 1 were imposed what prevents monomers to go too far away from each other. This methodology keeps the angular momentum constant, however, the concentration of the proteins in the solution can only be approximated. In some cases, the restrains are pronounced with “collidingcourse” association pattern especially in trajectories simulated in high temperature. It must, however, be noted that this pattern did not always occur, and influence association temperature and folding pattern only in a certain degree. Use of periodic boundary conditions (PBC) is necessary to reproduce the experimental observables that describe motion such as, for example, the diffusion constant in cross-linked protein membranes.[33] Implementation of PBC in molecular simulation algorithms is of constant interest,[34,35] there has been improvement in treatment of electrostatic interactions in 940
Journal of Computational Chemistry 2015, 36, 940–946
DOI: 10.1002/jcc.23864
Methodology UNRES model of the polypeptide chain In the UNRES,[12–28] force field polypeptide chain geometry is represented by sequence of Ca atoms connected by virtual bonds. There are two sites of interaction per residue. The peptide group (p), which is located half-way between two consecutive Ca atoms, and the side-chain (SC) group which is attached to Ca atom (Fig. 1). The geometry of backbone is ! defined by vectors ( dC ) connecting two consecutive Ca atoms or by internal coordinates, namely the valence angles h between three consecutive Ca atoms, the torsional angle c between four consecutive Ca atoms, and the d virtual-bond lengths distance between two consecutive Ca atoms (dCi). The position of a virtual side chain is described by either the vec! tor dX from Ca to the side chain or by internal coordinates (the virtual-bond length d, and the spherical angles a and b) positioning the SC with respect to backbone.
A. K. Sieradzan Department of Physics and Astronomy, Uppsala University, A˚ngstr€ omlaboratoriet, L€ agerhyddsv€ agen, 1,751 20 Uppsala, Sweden; Faculty of Chemistry, University of Gda nsk, Wita Stwosza 63, 80–308 Gda nsk, Poland E-mail: [email protected] Contract grant sponsor: Government of the Republic of Poland bugdet for years 2013–2014; Contract grant number: 0558/IP3/2013/72 C 2015 Wiley Periodicals, Inc. V
WWW.CHEMISTRYVIEWS.COM
SOFTWARE NEWS AND UPDATES
WWW.C-CHEM.ORG
fn ðTÞ5
Figure 1. The UNRES model of polypeptide chains. The interaction sites are peptide-group centers (p), and SC centers attached to the corresponding a-carbons with different Ca. . .SC bond lengths, dSC. The peptide groups are represented as dark gray circles and the side chains are represented as light gray ellipsoids of different size. The a-carbon atoms are represented by small open circles. The geometry of the chain can be described either by the virtual-bond vectors dCi (from Cai to Cai11 ), i51; 2; . . . ; n21, and dXi (from Cai to SCi), i52; . . . ; n21, represented by thick lines, where n is the number of residues, or in terms of virtual-bond lengths, backbone virtualbond angles hi ; i51; 2; . . . ; n22, backbone virtual-bond-dihedral angles ci ; i51; 2; . . . ; n23, and the angles ai and bi ; i52; 3; . . . ; n21 that describe the location of a side chain with respect to the coordinate frame defined by Cai21 ; Cai , and Cai ; Cai11 .
The effective energy function is a sum of factors which are represented by analytical approximation of a restricted free energy or a potential of mean force of a given conformation ensemble restricted to coarse-grained conformation defined by Ca atoms and SC groups and is expressed by [eq. (1)]: U5wSC
X
USCi SCj 1wSCp
i < XX U5 scaleðrlm Þ3Ulm > l i521 j521 k521 m : Ulm
if rlm > cutoff 1 if cutoff 2 < rlm < cutoff 1 if rlm < cutoff 2
(6) where scaleðrlm Þ511k2 ð2k23Þ
(7)
k5ðrlm 2cutoff 2 Þ=ðcutoff 1 2cutoff 2 Þ where the sum over i, j, k is sum over neighboring cells (when i 5 j 5 k 5 0 describes the simulated system - primary box). The scaling function is used for smoothing the energy function in the values between cutoff1 and cutoff2 to prevent sharp changes in a gradient (in contrary to case when single cutoff would be used). The smoothing function methodology is similar to the one used in CHARMM force filed.[42] 942
Journal of Computational Chemistry 2015, 36, 940–946
An introduction of the PBC or introduction of cutoff on long range interactions may lead to occurrence of artificial forces and energy drift. Therefore, the UNRES force field was tested if it keeps total energy constant. As a test for keeping constant energy, a microcanonical simulation was performed with two variants of the velocity integration algorithm: the variable time step (VTS)[43] or adaptive-multiple time step (A-MTS).[44] The simulation was performed for 8,000,000 steps with a time step length of 0.489 fs resulting in 4 ns UNRES time (corresponding to 4 ls of real time). As a starting structure octamer of minimized segment MVGGVVIA from the amyloid-beta peptide (Ab, residues 35–42, PDB code:2Y3K) was used. The cutoff1 and cutoff2 on long range interactions were set to 15 A˚ and 14.7 A˚, respectively. Canonical simulations with periodic boundary conditions Introduction of PBC disturbs the ratio of collision energy with respect to the kinetic energy.[45] This ratio is equal to N/ (N – 1), therefore, the smaller system is the larger the ratio becomes what may lead to disturbance of thermostat behavior. To verify if there was influence of the PBC, we measured temperature distribution and average temperature of the system. To verify the thermostat behavior canonical molecular dynamics simulation with the Berendsen thermostat[43,46] and Langevin dynamics[43] were performed. As a model peptide the same octamer system (Ab, residues 35–42, PDB code:2Y3K) was used as in microcanonical simulation. Both simulations were performed at 200, 300, 400, and 500 K with approximately 5,000,000 step with time step 4.89 fs. For Berendsen thermostat,[46] the coupling parameter s 5 48.9 fs was used. For Langevin dynamics, water friction was scaled by the factor of 0.01 to speed up the calculations, as in our earlier work.[47] The cutoff1 and cutoff2 on long range interactions were set to ˚ and 14.7 A ˚ , respectively. The influence of the PBC 15 A on thermostat behavior was tested. The results were compared with theoretical studies of PBC influence on kinetic energy.[45] WWW.CHEMISTRYVIEWS.COM
WWW.C-CHEM.ORG
Figure 3. The influence of the cutoff on long range interaction energy. The red line represents the energy component error of interactions between SC and peptide group for the charged system (PDB code: 3SGO); the green line represents energy component error of interactions between two peptide groups for the charged system (PDB code: 3SGO); the blue line represents the energy component error of interactions between SC and peptide group for the noncharged system (PDB code: 2Y3k); the black line represents energy component error of interactions between two peptide groups for the noncharged system (PDB code: 2Y3K).
Folding and association of bba hetero tetramer To determine the influence of the PBC on the association pathway and the heat capacity, the bba heterotetramer (BBAhet1; PDB code: 1XOF)[48] simulation was performed, analyzed and compared with droplet-like restraints.[49] The multiplexed-
Figure 4. The kinetic (black), potential (blue) and total energy (red) plot for microcanonical simulation with use of UNRES force field with PBC for VTS velocity integrator algorithm. Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]
SOFTWARE NEWS AND UPDATES
Figure 5. The kinetic (black), potential (blue) and total energy (red) plot for microcanonical simulation with use of UNRES force field with PBC for AMTS velocity integrator algorithm. Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]
replica exchange molecular dynamics[50] was used at the following 36 temperatures: 210, 220, 230, 240, 250, 255, 260, 265, 270, 275, 280, 285, 290, 295, 300, 305, 310, 315, 320, 325, 330, 335, 340, 345, 350, 360, 380, 390, 400, 410, 420, 430, 440, 460, 480, and 500K with 2 replicas for each temperature. The other technical details were the same as in our previous paper,[49] however, because of introduction of the PBC, the cutoff on the long range interactions was introduced. The cutoff1 and cutoff2 on long ˚ and 14.7 A ˚ , respectively. range interactions were set to 15 A
Figure 6. The total energy plot for VTS algorithm (red) and A-MTS algorithm (blue) in microcanonical simulation. Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]
Journal of Computational Chemistry 2015, 36, 940–946
943
SOFTWARE NEWS AND UPDATES
WWW.C-CHEM.ORG
Table 1. Comparision between the temperature set in thermostat and averages obtained in simulations. Temperature (K)
200 300 400 500
Langevin
Berndsen
Thermostat average (K)
Themostat average (K)
201.74 302.30 403.30 504.24
200.11 299.93 397.12 496.80
Results and Discussion Infulence of introduction of cutoff on long range interactions on energy The influence of introduction of cutoff on long range interactions was evaluated and results are shown in Figure 3. As it can be seen introduction of cutoff on long range interaction has a significant influence especially when cutoff less than 10 A˚ is used. Moreover, as shown in Figure 3 the SC peptide group and peptide group peptide group interactions are most prone to error associated with cutoff. The error of the other energy components was significantly smaller then these two. It was found that for both systems the cutoff of 15 A˚ is sufficitent for correct description of long range interactions with error less than 1%. Microcanonical simulations with periodic boundary conditions The results for microcanonical simulations in UNRES force field with PBC are shown in Figures 4–6. As it can be seen both VTS and A-MTS perform well as the fluctuations of total energy (Fig. 6) are significantly smaller than fluctuations of both, kinetic and potential energy (Figs. 4 and 5). It should be noted that both the A-MTS algorithm and VTS algorithm perform well in 4ns range in keeping constant energy as the total energy fluctuations are less than 0.5 kcal/mol. This indicates that no artificial forces are created and no energy drift is observed.
Figure 7. The comparison between theoretical thermal distribution (green) with Berendsen thermostat (red) and Langevin thermostat (blue). Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]
Canonical simulations with periodic boundary conditions To verify the influence of the PBC on the thermostat behavior the average temperatures were calculated (Table 1). As it can be seen in the Table 1 the Berendsen thermostat is not subject to artificial increase in the average temperature. This is due to coupling the kinetic energy with temperature. Conversely, Langevin thermostat seems to be influenced by the PBC. The temperature error increase linearly with temperature (R2 5 0.99). This result is consistent with theoretically predicted error.[45] Nevertheless, even for such a small system the effect of PBC on thermostat is negligible and the temperature is kept within 1% error with the thermostat-set temperature showing the thermostat behaviour is correct. As the implementation of PBC disturbs the thermostat behavior the distribution of the temperatures was studied. The obtained temperature distributions for Langevin and Berendsen thermostat were compared with
Figure 8. The comparison between heat capacity curves obtained in droplet-like environment A) and PBCs B). For the heat capacity curve in PBCs for bba heterotetramer different snapshots are ordered in rainbow style: red, orange, yellow, green, and blue. Figure 1A is reproduced with permission from F. C AIP Publishing LLC, 2014. Color figure can be viewed in the online issue, Yasar, A. K. Sieradzan, U. H. E. Hansmann, J. Chem. Phys. 2014, 140, 105103. V which is available at wileyonlinelibrary.com.]
944
Journal of Computational Chemistry 2015, 36, 940–946
WWW.CHEMISTRYVIEWS.COM
SOFTWARE NEWS AND UPDATES
WWW.C-CHEM.ORG
Figure 9. The comparision of the fraction of bba heterotetramer monomers (red), dimers (green), trimers (blue), and tetramers (violet) as a function of temperature for droplet-like environment A) and PBCs B). Figure 1A is reproduced with permission from F. Yasar, A. K. Sieradzan, U. H. E. Hansmann, J. Chem. C AIP Publishing LLC, 2014. Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.] Phys. 2014, 140, 105103. V
theoretical distribution {eq. (55) in Ref. [51]} and are shown in Figure 7. The introduction of PBC do not influence significantly the thermostat behavior; the Berendsen thermostat (red) gives the temperature distribution too narrow in comparison with the theoretical (green) as previously shown.[52] Whereas the Langevin (blue) thermostat gives similar temperature distribution as theoretically predicted. Folding of bba heterotetramer protein To evaluate the influence of previously used droplet-like restraints the heat capacity curves obtained in droplet-like environment (Fig. 8A) were compared with the ones obtained in PBC (Fig. 8B). As it can be seen, both methods yield virtually the same heat capacity curves. However, in the PBC the equilibrium (measured by heat capacity curve convergence) was obtained faster than in droplet-like restraint environment. It should be noted that no algorithm increasing replica exchange flow such as maximum replica exchange flow algorithm[53] was used. Therefore it shows that in PBC the equilibrium is obtained faster when no artificial forces are introduced showing additional advantage of PBC over droplet-like restraints. Conversely, when the association pattern is concerned, the differences are noticeable (Fig. 9). The association temperature is significantly lower Tassoc 325K compared with previously obtained Tassoc 375 K. Therefore, the high temperature shoulder (T 400–450K) seems to be rather dimerization associated, not tetramer association as previously[44] suggested. Moreover, even in the low temperature the trimer fraction is still significant and only 80% of the molecules are in form of tetramers as this is concentration dependent process and in PBC the concentration can be precisely controlled.
Conclusions In this work, PBCs were implemented in the UNRES force field. The influence of the cutoff on long range interactions was studied and shown that its influence is minimal when cutoff of ˚ is used. The UNRES force field with PBC keeps constant 15 A
energy for microcanonical simulations even for very long simulations proving that no artificial forces arise during the simulations. Moreover, the thermostat behaviour was verified. The Berendsen thermostat produces too narrow temperature distribution whereas the Langevin thermostat gives theoretically consistent thermal fluctuations. Conversely, the Langevin is affected by introduction of PBC as the average temperature is associated with error, however it was shown that this effect in negligible. The PBC give similar folding pathway as droplet-like restraints, however, the convergence is obtained much faster and association temperature is much closer to realistic vales. The newly implemented PBC are currently used for studying the arginine residue association with arginine-binding protein as well as studying of aggregation of amyloid proteins. In our future work, the Ewald’s summation[54] and isotropic periodic sum[55] will be tested and the speed efficiency will be compared. Moreover, the work on the other than only translational symmetry is currently under way.
Acknowledgments The author would like to thank Dr Bernard Brooks from NIH for scientific consultation in this topic. Computational resources were provided by (a) the supercomputer resources at the Informatics Center of the Metropolitan Academic Network (IC MAN) in Gdansk, (b) the 624-processor Beowulf cluster at the Baker Laboratory of Chemistry, Cornell University, and (c) our 184-processor Beowulf cluster at the Faculty of Chemistry, University of Gdansk. Keywords: multi-chain systems periodicity coarse-grain force field electrostatic interaction cutoff proteins
How to cite this article: A. K. Sieradzan J. Comput. Chem. 2015, 36, 940–946. DOI: 10.1002/jcc.23864
[1] H. Jeong, S. P. Mason, A.-L. Barabasi, and Z. N. Oltvai, Nature 2001, 411, 41.
Journal of Computational Chemistry 2015, 36, 940–946
945
SOFTWARE NEWS AND UPDATES
WWW.C-CHEM.ORG
[2] T. R. Jahn, M. J. Parker, S. W. Homans, S. E. Radford, Nat. Struct. Mol. Biol. 2006, 13, 195. [3] A. Rojas, A. Liwo, D. Browne, H. A. Scheraga, J. Mol. Biol. 2010, 404, 537. [4] G. Song, M. Zhou, W. Chen, T. Chen, B. Walker, C. Shaw, Biochem. Biophys. Res. Commun. 2008, 372, 191. [5] M. Scheffner, U. Nuber, J. M. Huibregtse, Nature 1995, 373, 81. [6] G. Zhao, J. R. Perilla, E. L. Yufenyuy, X. Meng, B. Chen, J. Ning, J. Ahn, A. M. Gronenborn, K. Schulten, C. Aiken, Nature 2013, 497, 643. [7] V. S. Pande, I. Baker, J. Chapman, S. Elmer, S. Kaliqand, S. M. Larson, Y. M. Rhee, M. R. Shirts, C. D. Snow, E. J. Sorinand, B. Zagrovic, Biopolymers 2003, 68, 91. [8] J. C. Phillips, R. Braun, W. Wang, J. Gumbart, E. Tajkhorshid, E. Villa, C. Chipot, R. D. Skeel, L. Kal e, K. Schulten, J. Comput. Chem. 2005, 26, 1781. [9] B. Hess, C. Kutzner, D. van der Spoel, E. Lindahl, J. Chem. Theory Comput. 2008, 4, 435. [10] D. E. Shaw, P. Maragakis, K. Lindorff-Larsen, S. Piana, R. O. Drorand, M. P. Eastwood, J. A. Bank, J. M. Jumper, J. K. Salmon, Y. Shan, W. Wriggers, Science 2010, 330, 341. [11] K. Lindorff-Larsen, S. Piana, R. O. Dror, D. E. Shaw, Science 2011, 334, 517. [12] A. Liwo, S. Ołdziej, M. R. Pincus, R. J. Wawak, S. Rackovsky, H. A. Scheraga, J. Comput. Chem. 1997, 18, 849. [13] A. Liwo, M. R. Pincus, R. J. Wawak, S. Rackovsky, S. Ołdziej, H. A. Scheraga, J. Comput. Chem. 1997, 18, 874. [14] A. Liwo, R. Kazmierkiewicz, C. Czaplewski, M. Groth, S. Ołdziej, R. J. Wawak, S. Rackovsky, M. R. Pincus, H. A. Scheraga, J. Comput. Chem. 1998, 19, 259. [15] A. Liwo, C. Czaplewski, J. Pillardy, H. A. Scheraga, J. Chem. Phys. 2001, 115, 2323. [16] A. Liwo, P. Arłukowicz, C. Czaplewski, S. Ołdziej, J. Pillardy, H. Scheraga, Proc. Natl. Acad. Sci. USA 2002, 99, 1937. [17] S. Ołdziej, U. Kozłowska, A. Liwo, H. A. Scheraga, J. Phys. Chem. A 2003, 107, 8035. [18] A. Liwo, S. Ołdziej, C. Czaplewski, U. Kozłowska, H. A. Scheraga, J. Phys. Chem. B 2004, 108, 9421. [19] S. Ołdziej, J. Ła˛giewka, A. Liwo, C. Czaplewski, M. Chinchio, M. Nanias, H. A. Scheraga, J. Phys. Chem. B 2004, 108, 16950. [20] U. Kozłowska, A. Liwo, H. A. Scheraga, J. Phys. Condens Matter 2007, 19, 285203. [21] A. Liwo, M. Khalili, C. Czaplewski, S. Kalinowski, S. Ołdziej, K. Wachucik, H. Scheraga, J. Phys. Chem. B 2007, 111, 260. [22] A. Liwo, C. Czaplewski, S. Ołdziej, A. V. Rojas, R. Kazmierkiewicz, M. Makowski, R. K. Murarka, H. A. Scheraga, In Coarse-Graining of Condensed Phase and Biomolecular Systems, Chapter 8; G. Voth, Ed.; CRC Press, 2008; pp. 1391–1411. [23] Y. He, Y. Xiao, A. Liwo, H. A. Scheraga, J. Comput. Chem. 2009, 30, 2127. [24] U. Kozłowska, A. Liwo, H. A. Scheraga, J. Comput. Chem. 2010, 31, 1143. [25] U. Kozłowska, G. G. Maisuradze, A. Liwo, H. A. Scheraga, J. Comput. Chem. 2010, 31, 1154. [26] A. K. Sieradzan, H. A. Scheraga, A. Liwo, J. Chem. Theor. Comput. 2012, 8, 1334. Available at: http://pubs.acs.org/doi/pdf/10.1021/ct2008439, http://pubs.acs.org/doi/abs/10.1021/ct2008439. [27] A. K. Sieradzan, U. H. E. Hansmann, H. A. Scheraga, A. Liwo, J. Chem. Theory Comput. 2012, 8, 4746. Available at: http://pubs.acs.org/doi/pdf/ 10.1021/ct3005563, http://pubs.acs.org/doi/abs/10.1021/ct3005563. [28] P. Krupa, A. K. Sieradzan, S. Rackovsky, M. Baranowski, S. Ołdziej, H. A. Scheraga, A. Liwo, C. Czaplewski, J. Chem. Theory Comput. 2013, 9, 4620. Available at: http://pubs.acs.org/doi/pdf/10.1021/ct4004977, http://pubs.acs.org/doi/abs/10.1021/ct4004977.
946
Journal of Computational Chemistry 2015, 36, 940–946
[29] Y. He, M. A. Mozolewska, P. Krupa, A. K. Sieradzan, T. K. Wirecki, A. Liwo, K. Kachlishvili, S. Rackovsky, D. Jagieła, R. Slusarz, C. R. Czaplewski, S. Ołdziej, H. A. Scheraga, Proc. Natl. Acad. Sci. USA 2013, 110, 14936. Available at: http://www.pnas.org/content/110/37/14936. full.pdf1html, http://www.pnas.org/content/110/37/14936.abstract. [30] A. Rojas, A. Liwo, H. Scheraga, J. Phys. Chem. 2011, B 115, 12978. [31] A. L. Y. He, H. Weinstein, H. Scheraga, J. Mol. Biol. 2011, 405, 298. [32] A. K. Sieradzan, A. Liwo, U. H. E. Hansmann, J. Chem. Theory Comput. 2012, 8, 3416. Available at: http://pubs.acs.org/doi/pdf/10.1021/ ct300528r, http://pubs.acs.org/doi/abs/10.1021/ct300528r. [33] G. Nair, J. F. Gargiuli, N. R. Shiju, Z. Rong, E. Shapiro, D. Drikakis, P. Vadgama, Chem. Bio. Chem. 2006, 7, 1683. [34] P. A. Schultz, Phys. Rev. B 1999, 60, 1551. [35] C. A. Rozzi, D. Varsano, A. Marini, E. K. U. Gross, A. Rubio, Phys. Rev. B 2006, 73, 205119. [36] M. Kalweit, D. Drikakis, Proc. Inst. Mech. Eng. Part C: J. Mech. Eng. Sci. 2008, 222, 797. [37] K. M. Stocker, Ph.D. Thesis, University of Notre Dame, 2014. [38] R. Kubo, J. Phys. Soc. Jpn. 1962, 17, 1100. [39] H. Shen, A. Liwo, H. A. Scheraga, J. Phys. Chem. B 2009, 113, 8738. [40] A. Kolinski, J. Skolnick, J. Chem. Phys. 1992, 97, 9412. [41] E. A. Dolan, R. M. Venable, R. W. Pastor, B. R. Brooks, Biophys. J. 2002, 82, 2317. [42] P. J. Steinbach, B. R. Brooks, J. Comput. Chem. 1994, 15, 667. ISSN 1096-987X, Available at: http://dx.doi.org/10.1002/jcc.540150702. [43] M. Khalili, A. Liwo, F. Rakowski, P. Grochowski, H. Scheraga, J. Phys. Chem. B 2005, 109, 13785. [44] F. Rakowski, P. Grochowski, B. Lesyng, A. Liwo, H. A. Scheraga, J. Chem. Phys. 2006, 125, 204107. [45] R. B. Shirts, S. R. Burt, A. M. Johnson, J. Chem. Phys. 2006, 125, 164102. Available at: http://scitation.aip.org/content/aip/journal/jcp/125/16/10. 1063/1.2359432. [46] H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, A. DiNola, J. R. Haak, J. Chem. Phys. 1984, 81, 3684. [47] M. Khalili, A. Liwo, A. Jagielska, H. Scheraga, J. Phys. Chem. B 2005, 109, 13798. [48] M. H. Ali, C. M. Taylor, G. Grigoryan, K. N. Allen, B. Imperiali, A. E. Keating, Structure 2005, 13, 225. [49] F. Yasar, A. K. Sieradzan, U. H. E. Hansmann, J. Chem. Phys. 2014, 140, 105103. Available at: http://scitation.aip.org/content/aip/journal/jcp/ 140/10/10.1063/1.4868140. [50] M. Nanias, C. Czaplewski, H. A. Scheraga, J. Chem. Theor. Comput. 2006, 2, 513. [51] D. S. Kleinerman, C. Czaplewski, A. Liwo, H. A. Scheraga, J. Chem. Phys. 2008, 128, 245103. [52] M. D’Alessandro, A. Tenenbaum, A. Amadei, J. Phys. Chem. B 2002, 106, 5050. Available at: http://pubs.acs.org/doi/pdf/10.1021/jp013689i, http://pubs.acs.org/doi/abs/10.1021/jp013689i. [53] W. Nadler, J. H. Meinke, U. H. E. Hansmann, Phys. Rev. E 2008, 78, 061905. Available at: http://link.aps.org/doi/10.1103/PhysRevE.78.061905. [54] T. I. Cheatham, J. Miller, T. Fox, T. Darden, P. Kollman, J. Am. Chem. Soc. 1995, 117, 4193. [55] X. Wu, B. R. Brooks, J. Chem. Phys. 2008, 129, 154115.
Received: 10 November 2014 Revised: 15 January 2015 Accepted: 17 January 2015 Published online on 8 March 2015
WWW.CHEMISTRYVIEWS.COM
Copyright of Journal of Computational Chemistry is the property of John Wiley & Sons, Inc. and its content may not be copied or emailed to multiple sites or posted to a listserv without the copyright holder's express written permission. However, users may print, download, or email articles for individual use.
WWW.C-CHEM.ORG
Introduction of Periodic Boundary Conditions into UNRES Force Field Adam K. Sieradzan* In this article, implementation of periodic boundary conditions (PBC) into physics-based coarse-grained UNited RESidue (UNRES) force field is presented, which replaces droplet-like restraints previously used. Droplet-like restraints are necessary to keep multichain systems together and prevent them from dissolving to infinitely low concentration. As an alternative for droplet-like restrains cuboid PBCs with imaging of the molecules were introduced. Owing to this modification, artificial forces which arose from restraints keeping a droplet together were eliminated what leads to more realistic trajectories. Due
to computational reasons cutoff and smoothing functions were introduced on the long range interactions. The UNRES force field with PBC was tested by performing microcanonical simulations. Moreover, to asses the behavior of the thermostat in PBCs Langevin and Berendsen thermostats were studied. The influence of PBCs on association pattern was compared with droplet-like restraints on the bba hetero tetramer 1 proC 2015 Wiley Periodicals, Inc. tein system. V
Introduction
molecular.[34,35] Moreover, there are also studies to treat system as a nonperiodic for instance continuum-molecular system. An interesting study about nonperiodic systems has been presented in Refs [36,37]. In this work, a cuboid PBC with cutoff on long range interactions has been implemented to run coarse-grained molecular simulations with the UNRES force field.
Protein–protein interactions play an extremely important role in living organisms and understanding of these interactions is a necessary condition to understand how the organisms are functioning.[1] Amyloid growth[2,3] and their intermediates behavior,[2] enzyme-inhibitor interactions,[4] and modulatory functions of proteins[5] are only a few examples showing the importance of protein–protein interactions. Theoretical studies of such interactions have provided valuable information and explanation of experimental results.[6] Moreover, there has been tremendous advance in all-atom simulation recently,[7–11] allowing to simulate small proteins even in ms time-scale.[11] However, for large systems coarse-graining is still required. The UNited RESidue (UNRES)[12–28] force field is a carefully derived physic-based force-field which has proven to be a powerful tool in protein structure prediction.[29] The previous version of UNRES force field was extended to simulate oligomeric proteins.[3,30–32] To keep the monomers together, the restrains C-terminus of chain i and N-terminus of chain i 1 1 were imposed what prevents monomers to go too far away from each other. This methodology keeps the angular momentum constant, however, the concentration of the proteins in the solution can only be approximated. In some cases, the restrains are pronounced with “collidingcourse” association pattern especially in trajectories simulated in high temperature. It must, however, be noted that this pattern did not always occur, and influence association temperature and folding pattern only in a certain degree. Use of periodic boundary conditions (PBC) is necessary to reproduce the experimental observables that describe motion such as, for example, the diffusion constant in cross-linked protein membranes.[33] Implementation of PBC in molecular simulation algorithms is of constant interest,[34,35] there has been improvement in treatment of electrostatic interactions in 940
Journal of Computational Chemistry 2015, 36, 940–946
DOI: 10.1002/jcc.23864
Methodology UNRES model of the polypeptide chain In the UNRES,[12–28] force field polypeptide chain geometry is represented by sequence of Ca atoms connected by virtual bonds. There are two sites of interaction per residue. The peptide group (p), which is located half-way between two consecutive Ca atoms, and the side-chain (SC) group which is attached to Ca atom (Fig. 1). The geometry of backbone is ! defined by vectors ( dC ) connecting two consecutive Ca atoms or by internal coordinates, namely the valence angles h between three consecutive Ca atoms, the torsional angle c between four consecutive Ca atoms, and the d virtual-bond lengths distance between two consecutive Ca atoms (dCi). The position of a virtual side chain is described by either the vec! tor dX from Ca to the side chain or by internal coordinates (the virtual-bond length d, and the spherical angles a and b) positioning the SC with respect to backbone.
A. K. Sieradzan Department of Physics and Astronomy, Uppsala University, A˚ngstr€ omlaboratoriet, L€ agerhyddsv€ agen, 1,751 20 Uppsala, Sweden; Faculty of Chemistry, University of Gda nsk, Wita Stwosza 63, 80–308 Gda nsk, Poland E-mail: [email protected] Contract grant sponsor: Government of the Republic of Poland bugdet for years 2013–2014; Contract grant number: 0558/IP3/2013/72 C 2015 Wiley Periodicals, Inc. V
WWW.CHEMISTRYVIEWS.COM
SOFTWARE NEWS AND UPDATES
WWW.C-CHEM.ORG
fn ðTÞ5
Figure 1. The UNRES model of polypeptide chains. The interaction sites are peptide-group centers (p), and SC centers attached to the corresponding a-carbons with different Ca. . .SC bond lengths, dSC. The peptide groups are represented as dark gray circles and the side chains are represented as light gray ellipsoids of different size. The a-carbon atoms are represented by small open circles. The geometry of the chain can be described either by the virtual-bond vectors dCi (from Cai to Cai11 ), i51; 2; . . . ; n21, and dXi (from Cai to SCi), i52; . . . ; n21, represented by thick lines, where n is the number of residues, or in terms of virtual-bond lengths, backbone virtualbond angles hi ; i51; 2; . . . ; n22, backbone virtual-bond-dihedral angles ci ; i51; 2; . . . ; n23, and the angles ai and bi ; i52; 3; . . . ; n21 that describe the location of a side chain with respect to the coordinate frame defined by Cai21 ; Cai , and Cai ; Cai11 .
The effective energy function is a sum of factors which are represented by analytical approximation of a restricted free energy or a potential of mean force of a given conformation ensemble restricted to coarse-grained conformation defined by Ca atoms and SC groups and is expressed by [eq. (1)]: U5wSC
X
USCi SCj 1wSCp
i < XX U5 scaleðrlm Þ3Ulm > l i521 j521 k521 m : Ulm
if rlm > cutoff 1 if cutoff 2 < rlm < cutoff 1 if rlm < cutoff 2
(6) where scaleðrlm Þ511k2 ð2k23Þ
(7)
k5ðrlm 2cutoff 2 Þ=ðcutoff 1 2cutoff 2 Þ where the sum over i, j, k is sum over neighboring cells (when i 5 j 5 k 5 0 describes the simulated system - primary box). The scaling function is used for smoothing the energy function in the values between cutoff1 and cutoff2 to prevent sharp changes in a gradient (in contrary to case when single cutoff would be used). The smoothing function methodology is similar to the one used in CHARMM force filed.[42] 942
Journal of Computational Chemistry 2015, 36, 940–946
An introduction of the PBC or introduction of cutoff on long range interactions may lead to occurrence of artificial forces and energy drift. Therefore, the UNRES force field was tested if it keeps total energy constant. As a test for keeping constant energy, a microcanonical simulation was performed with two variants of the velocity integration algorithm: the variable time step (VTS)[43] or adaptive-multiple time step (A-MTS).[44] The simulation was performed for 8,000,000 steps with a time step length of 0.489 fs resulting in 4 ns UNRES time (corresponding to 4 ls of real time). As a starting structure octamer of minimized segment MVGGVVIA from the amyloid-beta peptide (Ab, residues 35–42, PDB code:2Y3K) was used. The cutoff1 and cutoff2 on long range interactions were set to 15 A˚ and 14.7 A˚, respectively. Canonical simulations with periodic boundary conditions Introduction of PBC disturbs the ratio of collision energy with respect to the kinetic energy.[45] This ratio is equal to N/ (N – 1), therefore, the smaller system is the larger the ratio becomes what may lead to disturbance of thermostat behavior. To verify if there was influence of the PBC, we measured temperature distribution and average temperature of the system. To verify the thermostat behavior canonical molecular dynamics simulation with the Berendsen thermostat[43,46] and Langevin dynamics[43] were performed. As a model peptide the same octamer system (Ab, residues 35–42, PDB code:2Y3K) was used as in microcanonical simulation. Both simulations were performed at 200, 300, 400, and 500 K with approximately 5,000,000 step with time step 4.89 fs. For Berendsen thermostat,[46] the coupling parameter s 5 48.9 fs was used. For Langevin dynamics, water friction was scaled by the factor of 0.01 to speed up the calculations, as in our earlier work.[47] The cutoff1 and cutoff2 on long range interactions were set to ˚ and 14.7 A ˚ , respectively. The influence of the PBC 15 A on thermostat behavior was tested. The results were compared with theoretical studies of PBC influence on kinetic energy.[45] WWW.CHEMISTRYVIEWS.COM
WWW.C-CHEM.ORG
Figure 3. The influence of the cutoff on long range interaction energy. The red line represents the energy component error of interactions between SC and peptide group for the charged system (PDB code: 3SGO); the green line represents energy component error of interactions between two peptide groups for the charged system (PDB code: 3SGO); the blue line represents the energy component error of interactions between SC and peptide group for the noncharged system (PDB code: 2Y3k); the black line represents energy component error of interactions between two peptide groups for the noncharged system (PDB code: 2Y3K).
Folding and association of bba hetero tetramer To determine the influence of the PBC on the association pathway and the heat capacity, the bba heterotetramer (BBAhet1; PDB code: 1XOF)[48] simulation was performed, analyzed and compared with droplet-like restraints.[49] The multiplexed-
Figure 4. The kinetic (black), potential (blue) and total energy (red) plot for microcanonical simulation with use of UNRES force field with PBC for VTS velocity integrator algorithm. Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]
SOFTWARE NEWS AND UPDATES
Figure 5. The kinetic (black), potential (blue) and total energy (red) plot for microcanonical simulation with use of UNRES force field with PBC for AMTS velocity integrator algorithm. Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]
replica exchange molecular dynamics[50] was used at the following 36 temperatures: 210, 220, 230, 240, 250, 255, 260, 265, 270, 275, 280, 285, 290, 295, 300, 305, 310, 315, 320, 325, 330, 335, 340, 345, 350, 360, 380, 390, 400, 410, 420, 430, 440, 460, 480, and 500K with 2 replicas for each temperature. The other technical details were the same as in our previous paper,[49] however, because of introduction of the PBC, the cutoff on the long range interactions was introduced. The cutoff1 and cutoff2 on long ˚ and 14.7 A ˚ , respectively. range interactions were set to 15 A
Figure 6. The total energy plot for VTS algorithm (red) and A-MTS algorithm (blue) in microcanonical simulation. Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]
Journal of Computational Chemistry 2015, 36, 940–946
943
SOFTWARE NEWS AND UPDATES
WWW.C-CHEM.ORG
Table 1. Comparision between the temperature set in thermostat and averages obtained in simulations. Temperature (K)
200 300 400 500
Langevin
Berndsen
Thermostat average (K)
Themostat average (K)
201.74 302.30 403.30 504.24
200.11 299.93 397.12 496.80
Results and Discussion Infulence of introduction of cutoff on long range interactions on energy The influence of introduction of cutoff on long range interactions was evaluated and results are shown in Figure 3. As it can be seen introduction of cutoff on long range interaction has a significant influence especially when cutoff less than 10 A˚ is used. Moreover, as shown in Figure 3 the SC peptide group and peptide group peptide group interactions are most prone to error associated with cutoff. The error of the other energy components was significantly smaller then these two. It was found that for both systems the cutoff of 15 A˚ is sufficitent for correct description of long range interactions with error less than 1%. Microcanonical simulations with periodic boundary conditions The results for microcanonical simulations in UNRES force field with PBC are shown in Figures 4–6. As it can be seen both VTS and A-MTS perform well as the fluctuations of total energy (Fig. 6) are significantly smaller than fluctuations of both, kinetic and potential energy (Figs. 4 and 5). It should be noted that both the A-MTS algorithm and VTS algorithm perform well in 4ns range in keeping constant energy as the total energy fluctuations are less than 0.5 kcal/mol. This indicates that no artificial forces are created and no energy drift is observed.
Figure 7. The comparison between theoretical thermal distribution (green) with Berendsen thermostat (red) and Langevin thermostat (blue). Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]
Canonical simulations with periodic boundary conditions To verify the influence of the PBC on the thermostat behavior the average temperatures were calculated (Table 1). As it can be seen in the Table 1 the Berendsen thermostat is not subject to artificial increase in the average temperature. This is due to coupling the kinetic energy with temperature. Conversely, Langevin thermostat seems to be influenced by the PBC. The temperature error increase linearly with temperature (R2 5 0.99). This result is consistent with theoretically predicted error.[45] Nevertheless, even for such a small system the effect of PBC on thermostat is negligible and the temperature is kept within 1% error with the thermostat-set temperature showing the thermostat behaviour is correct. As the implementation of PBC disturbs the thermostat behavior the distribution of the temperatures was studied. The obtained temperature distributions for Langevin and Berendsen thermostat were compared with
Figure 8. The comparison between heat capacity curves obtained in droplet-like environment A) and PBCs B). For the heat capacity curve in PBCs for bba heterotetramer different snapshots are ordered in rainbow style: red, orange, yellow, green, and blue. Figure 1A is reproduced with permission from F. C AIP Publishing LLC, 2014. Color figure can be viewed in the online issue, Yasar, A. K. Sieradzan, U. H. E. Hansmann, J. Chem. Phys. 2014, 140, 105103. V which is available at wileyonlinelibrary.com.]
944
Journal of Computational Chemistry 2015, 36, 940–946
WWW.CHEMISTRYVIEWS.COM
SOFTWARE NEWS AND UPDATES
WWW.C-CHEM.ORG
Figure 9. The comparision of the fraction of bba heterotetramer monomers (red), dimers (green), trimers (blue), and tetramers (violet) as a function of temperature for droplet-like environment A) and PBCs B). Figure 1A is reproduced with permission from F. Yasar, A. K. Sieradzan, U. H. E. Hansmann, J. Chem. C AIP Publishing LLC, 2014. Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.] Phys. 2014, 140, 105103. V
theoretical distribution {eq. (55) in Ref. [51]} and are shown in Figure 7. The introduction of PBC do not influence significantly the thermostat behavior; the Berendsen thermostat (red) gives the temperature distribution too narrow in comparison with the theoretical (green) as previously shown.[52] Whereas the Langevin (blue) thermostat gives similar temperature distribution as theoretically predicted. Folding of bba heterotetramer protein To evaluate the influence of previously used droplet-like restraints the heat capacity curves obtained in droplet-like environment (Fig. 8A) were compared with the ones obtained in PBC (Fig. 8B). As it can be seen, both methods yield virtually the same heat capacity curves. However, in the PBC the equilibrium (measured by heat capacity curve convergence) was obtained faster than in droplet-like restraint environment. It should be noted that no algorithm increasing replica exchange flow such as maximum replica exchange flow algorithm[53] was used. Therefore it shows that in PBC the equilibrium is obtained faster when no artificial forces are introduced showing additional advantage of PBC over droplet-like restraints. Conversely, when the association pattern is concerned, the differences are noticeable (Fig. 9). The association temperature is significantly lower Tassoc 325K compared with previously obtained Tassoc 375 K. Therefore, the high temperature shoulder (T 400–450K) seems to be rather dimerization associated, not tetramer association as previously[44] suggested. Moreover, even in the low temperature the trimer fraction is still significant and only 80% of the molecules are in form of tetramers as this is concentration dependent process and in PBC the concentration can be precisely controlled.
Conclusions In this work, PBCs were implemented in the UNRES force field. The influence of the cutoff on long range interactions was studied and shown that its influence is minimal when cutoff of ˚ is used. The UNRES force field with PBC keeps constant 15 A
energy for microcanonical simulations even for very long simulations proving that no artificial forces arise during the simulations. Moreover, the thermostat behaviour was verified. The Berendsen thermostat produces too narrow temperature distribution whereas the Langevin thermostat gives theoretically consistent thermal fluctuations. Conversely, the Langevin is affected by introduction of PBC as the average temperature is associated with error, however it was shown that this effect in negligible. The PBC give similar folding pathway as droplet-like restraints, however, the convergence is obtained much faster and association temperature is much closer to realistic vales. The newly implemented PBC are currently used for studying the arginine residue association with arginine-binding protein as well as studying of aggregation of amyloid proteins. In our future work, the Ewald’s summation[54] and isotropic periodic sum[55] will be tested and the speed efficiency will be compared. Moreover, the work on the other than only translational symmetry is currently under way.
Acknowledgments The author would like to thank Dr Bernard Brooks from NIH for scientific consultation in this topic. Computational resources were provided by (a) the supercomputer resources at the Informatics Center of the Metropolitan Academic Network (IC MAN) in Gdansk, (b) the 624-processor Beowulf cluster at the Baker Laboratory of Chemistry, Cornell University, and (c) our 184-processor Beowulf cluster at the Faculty of Chemistry, University of Gdansk. Keywords: multi-chain systems periodicity coarse-grain force field electrostatic interaction cutoff proteins
How to cite this article: A. K. Sieradzan J. Comput. Chem. 2015, 36, 940–946. DOI: 10.1002/jcc.23864
[1] H. Jeong, S. P. Mason, A.-L. Barabasi, and Z. N. Oltvai, Nature 2001, 411, 41.
Journal of Computational Chemistry 2015, 36, 940–946
945
SOFTWARE NEWS AND UPDATES
WWW.C-CHEM.ORG
[2] T. R. Jahn, M. J. Parker, S. W. Homans, S. E. Radford, Nat. Struct. Mol. Biol. 2006, 13, 195. [3] A. Rojas, A. Liwo, D. Browne, H. A. Scheraga, J. Mol. Biol. 2010, 404, 537. [4] G. Song, M. Zhou, W. Chen, T. Chen, B. Walker, C. Shaw, Biochem. Biophys. Res. Commun. 2008, 372, 191. [5] M. Scheffner, U. Nuber, J. M. Huibregtse, Nature 1995, 373, 81. [6] G. Zhao, J. R. Perilla, E. L. Yufenyuy, X. Meng, B. Chen, J. Ning, J. Ahn, A. M. Gronenborn, K. Schulten, C. Aiken, Nature 2013, 497, 643. [7] V. S. Pande, I. Baker, J. Chapman, S. Elmer, S. Kaliqand, S. M. Larson, Y. M. Rhee, M. R. Shirts, C. D. Snow, E. J. Sorinand, B. Zagrovic, Biopolymers 2003, 68, 91. [8] J. C. Phillips, R. Braun, W. Wang, J. Gumbart, E. Tajkhorshid, E. Villa, C. Chipot, R. D. Skeel, L. Kal e, K. Schulten, J. Comput. Chem. 2005, 26, 1781. [9] B. Hess, C. Kutzner, D. van der Spoel, E. Lindahl, J. Chem. Theory Comput. 2008, 4, 435. [10] D. E. Shaw, P. Maragakis, K. Lindorff-Larsen, S. Piana, R. O. Drorand, M. P. Eastwood, J. A. Bank, J. M. Jumper, J. K. Salmon, Y. Shan, W. Wriggers, Science 2010, 330, 341. [11] K. Lindorff-Larsen, S. Piana, R. O. Dror, D. E. Shaw, Science 2011, 334, 517. [12] A. Liwo, S. Ołdziej, M. R. Pincus, R. J. Wawak, S. Rackovsky, H. A. Scheraga, J. Comput. Chem. 1997, 18, 849. [13] A. Liwo, M. R. Pincus, R. J. Wawak, S. Rackovsky, S. Ołdziej, H. A. Scheraga, J. Comput. Chem. 1997, 18, 874. [14] A. Liwo, R. Kazmierkiewicz, C. Czaplewski, M. Groth, S. Ołdziej, R. J. Wawak, S. Rackovsky, M. R. Pincus, H. A. Scheraga, J. Comput. Chem. 1998, 19, 259. [15] A. Liwo, C. Czaplewski, J. Pillardy, H. A. Scheraga, J. Chem. Phys. 2001, 115, 2323. [16] A. Liwo, P. Arłukowicz, C. Czaplewski, S. Ołdziej, J. Pillardy, H. Scheraga, Proc. Natl. Acad. Sci. USA 2002, 99, 1937. [17] S. Ołdziej, U. Kozłowska, A. Liwo, H. A. Scheraga, J. Phys. Chem. A 2003, 107, 8035. [18] A. Liwo, S. Ołdziej, C. Czaplewski, U. Kozłowska, H. A. Scheraga, J. Phys. Chem. B 2004, 108, 9421. [19] S. Ołdziej, J. Ła˛giewka, A. Liwo, C. Czaplewski, M. Chinchio, M. Nanias, H. A. Scheraga, J. Phys. Chem. B 2004, 108, 16950. [20] U. Kozłowska, A. Liwo, H. A. Scheraga, J. Phys. Condens Matter 2007, 19, 285203. [21] A. Liwo, M. Khalili, C. Czaplewski, S. Kalinowski, S. Ołdziej, K. Wachucik, H. Scheraga, J. Phys. Chem. B 2007, 111, 260. [22] A. Liwo, C. Czaplewski, S. Ołdziej, A. V. Rojas, R. Kazmierkiewicz, M. Makowski, R. K. Murarka, H. A. Scheraga, In Coarse-Graining of Condensed Phase and Biomolecular Systems, Chapter 8; G. Voth, Ed.; CRC Press, 2008; pp. 1391–1411. [23] Y. He, Y. Xiao, A. Liwo, H. A. Scheraga, J. Comput. Chem. 2009, 30, 2127. [24] U. Kozłowska, A. Liwo, H. A. Scheraga, J. Comput. Chem. 2010, 31, 1143. [25] U. Kozłowska, G. G. Maisuradze, A. Liwo, H. A. Scheraga, J. Comput. Chem. 2010, 31, 1154. [26] A. K. Sieradzan, H. A. Scheraga, A. Liwo, J. Chem. Theor. Comput. 2012, 8, 1334. Available at: http://pubs.acs.org/doi/pdf/10.1021/ct2008439, http://pubs.acs.org/doi/abs/10.1021/ct2008439. [27] A. K. Sieradzan, U. H. E. Hansmann, H. A. Scheraga, A. Liwo, J. Chem. Theory Comput. 2012, 8, 4746. Available at: http://pubs.acs.org/doi/pdf/ 10.1021/ct3005563, http://pubs.acs.org/doi/abs/10.1021/ct3005563. [28] P. Krupa, A. K. Sieradzan, S. Rackovsky, M. Baranowski, S. Ołdziej, H. A. Scheraga, A. Liwo, C. Czaplewski, J. Chem. Theory Comput. 2013, 9, 4620. Available at: http://pubs.acs.org/doi/pdf/10.1021/ct4004977, http://pubs.acs.org/doi/abs/10.1021/ct4004977.
946
Journal of Computational Chemistry 2015, 36, 940–946
[29] Y. He, M. A. Mozolewska, P. Krupa, A. K. Sieradzan, T. K. Wirecki, A. Liwo, K. Kachlishvili, S. Rackovsky, D. Jagieła, R. Slusarz, C. R. Czaplewski, S. Ołdziej, H. A. Scheraga, Proc. Natl. Acad. Sci. USA 2013, 110, 14936. Available at: http://www.pnas.org/content/110/37/14936. full.pdf1html, http://www.pnas.org/content/110/37/14936.abstract. [30] A. Rojas, A. Liwo, H. Scheraga, J. Phys. Chem. 2011, B 115, 12978. [31] A. L. Y. He, H. Weinstein, H. Scheraga, J. Mol. Biol. 2011, 405, 298. [32] A. K. Sieradzan, A. Liwo, U. H. E. Hansmann, J. Chem. Theory Comput. 2012, 8, 3416. Available at: http://pubs.acs.org/doi/pdf/10.1021/ ct300528r, http://pubs.acs.org/doi/abs/10.1021/ct300528r. [33] G. Nair, J. F. Gargiuli, N. R. Shiju, Z. Rong, E. Shapiro, D. Drikakis, P. Vadgama, Chem. Bio. Chem. 2006, 7, 1683. [34] P. A. Schultz, Phys. Rev. B 1999, 60, 1551. [35] C. A. Rozzi, D. Varsano, A. Marini, E. K. U. Gross, A. Rubio, Phys. Rev. B 2006, 73, 205119. [36] M. Kalweit, D. Drikakis, Proc. Inst. Mech. Eng. Part C: J. Mech. Eng. Sci. 2008, 222, 797. [37] K. M. Stocker, Ph.D. Thesis, University of Notre Dame, 2014. [38] R. Kubo, J. Phys. Soc. Jpn. 1962, 17, 1100. [39] H. Shen, A. Liwo, H. A. Scheraga, J. Phys. Chem. B 2009, 113, 8738. [40] A. Kolinski, J. Skolnick, J. Chem. Phys. 1992, 97, 9412. [41] E. A. Dolan, R. M. Venable, R. W. Pastor, B. R. Brooks, Biophys. J. 2002, 82, 2317. [42] P. J. Steinbach, B. R. Brooks, J. Comput. Chem. 1994, 15, 667. ISSN 1096-987X, Available at: http://dx.doi.org/10.1002/jcc.540150702. [43] M. Khalili, A. Liwo, F. Rakowski, P. Grochowski, H. Scheraga, J. Phys. Chem. B 2005, 109, 13785. [44] F. Rakowski, P. Grochowski, B. Lesyng, A. Liwo, H. A. Scheraga, J. Chem. Phys. 2006, 125, 204107. [45] R. B. Shirts, S. R. Burt, A. M. Johnson, J. Chem. Phys. 2006, 125, 164102. Available at: http://scitation.aip.org/content/aip/journal/jcp/125/16/10. 1063/1.2359432. [46] H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, A. DiNola, J. R. Haak, J. Chem. Phys. 1984, 81, 3684. [47] M. Khalili, A. Liwo, A. Jagielska, H. Scheraga, J. Phys. Chem. B 2005, 109, 13798. [48] M. H. Ali, C. M. Taylor, G. Grigoryan, K. N. Allen, B. Imperiali, A. E. Keating, Structure 2005, 13, 225. [49] F. Yasar, A. K. Sieradzan, U. H. E. Hansmann, J. Chem. Phys. 2014, 140, 105103. Available at: http://scitation.aip.org/content/aip/journal/jcp/ 140/10/10.1063/1.4868140. [50] M. Nanias, C. Czaplewski, H. A. Scheraga, J. Chem. Theor. Comput. 2006, 2, 513. [51] D. S. Kleinerman, C. Czaplewski, A. Liwo, H. A. Scheraga, J. Chem. Phys. 2008, 128, 245103. [52] M. D’Alessandro, A. Tenenbaum, A. Amadei, J. Phys. Chem. B 2002, 106, 5050. Available at: http://pubs.acs.org/doi/pdf/10.1021/jp013689i, http://pubs.acs.org/doi/abs/10.1021/jp013689i. [53] W. Nadler, J. H. Meinke, U. H. E. Hansmann, Phys. Rev. E 2008, 78, 061905. Available at: http://link.aps.org/doi/10.1103/PhysRevE.78.061905. [54] T. I. Cheatham, J. Miller, T. Fox, T. Darden, P. Kollman, J. Am. Chem. Soc. 1995, 117, 4193. [55] X. Wu, B. R. Brooks, J. Chem. Phys. 2008, 129, 154115.
Received: 10 November 2014 Revised: 15 January 2015 Accepted: 17 January 2015 Published online on 8 March 2015
WWW.CHEMISTRYVIEWS.COM
Copyright of Journal of Computational Chemistry is the property of John Wiley & Sons, Inc. and its content may not be copied or emailed to multiple sites or posted to a listserv without the copyright holder's express written permission. However, users may print, download, or email articles for individual use.
