Return

The finite difference solution of the Poisson-Boltzmann (PB) equations based upon rigorous physics has been widely used for the calculation of the electrostatic contribution to the solvation energy as the most accurate continuum model [6]. Until recently, use of the PB approach has been largely restricted to calculations involving static representations of molecular structure, but the recent development of simple methods to obtain solvation forces from the PB equation made it possible to use it in molecular dynamics (MD) simulations [7]. However, obtaining accurate numerical solution of PB equation for a large system such as a protein is still too costly to permit useful long time dynamics of biological molecules to be routinely studied. Another simpler and faster implicit continuum model based on PB equation is the so-called generalized Born/solvent-accessible surface area (GB/SA) model introduced by Still and co-workers [8]. This GB/SA model has been developed and applied for molecular dynamics simulations of proteins and nucleic acids [9 - 11]. In more simplified implicit solvent models, all contributions to the solvation free energy are assumed to be proportional to the solvent-accessible surface area [12 - 14] or Gaussian solvent exclusion volume [15].

Recently, we have developed an atomistic Brownian dynamics (BD) with multiple time step algorithm for the simulations that require long time calculations such as protein folding [16, 17]. BD simulation is one of the stochastic dynamics simulations in which the dynamical aspects of the solvent are treated as a dissipative random force in the equations of motion. Furthermore, by eliminating the less interesting high-frequency motions such as those arising from bond vibrations, one can use larger time steps in BD simulations [18 - 20]. So, BD simulation enables much longer simulations to be performed. In the BD algorithm, a protein was described by united-atom model with AMBER91 force field and solvation effect was introduced by distance-dependent dielectric and solvent-accessible surface area (DD/SA) model. This BD simulation of a peptide was 160 times faster than molecular dynamics simulation with explicit solvent and stable.

In this paper, we introduce a new simple implicit solvent model, effective charge (EC) model to improve the force field. In the EC model, the atomic charge is neutralized as a function of the solvent-accessible surface area of the atom. The GB/SA model was also applied for BD algorithm. In order to test the effectiveness of these implicit solvent models, BD simulations of a 28-mer bba folded peptide were carried out. A corresponding molecular dynamics simulation using explicit all-atom water model was also performed for comparison.

(1) DD model

DD model assumes the dielectric constant e = C

(2) SA model

In SA model, the solvation energy is given by:

where

(3) EC model

We introduce a new implicit solvent model, Effective Charge (EC) model. Water is a polar molecule with a dipole moment of 1.85 debye. Therefore, the point charge in water solution is surrounded by oriented water molecules and the point charge is shielded. To represent this shielding effect, atomic charge of atom

Here

This EC model is similar to the solvent-accessibility modified Tanford-Kirkwood method developed by Matthew and Gurd [24]. In their model, the electrostatic energy of pairwise interaction is calculated as the value of the Tanford-Kirkwood electrostatic interaction energy reduced by average solvent-accessibility (normalized by the corresponding accessible surface area of a model tripeptide) over atoms

(4) GB/SA model

In the GB/SA model [8], the total solvation energy (

The first two terms are the non-electrostatic contributions and are approximately linearly related to their solvent-accessible surface area (SA). This is represented by the sum:

where g

The third term,

where is the distance between atoms

In this study, the three implicit solvent models described above were used. However, DD and EC models include only the electrostatic contribution to the free energy of solvation and neglect the non-electrostatic contribution. To consider both contributions, DD and EC models were combined with SA model, that is, DD/SA and DD/SA/EC models.

Figure 1 shows Ca root mean square deviations (RMSD) of the peptide from the NMR structure as a function of the simulation time. A snapshot at 5 ns period of each trajectory is shown in Figure 2. By

Figure 1. Time evolutions of the Ca RMSD during simulations. The RMS deviations from the NMR structure versus time for 5 ns simulations using various solvent models are shown.

Figure 2. Ribbon representations of the snapshots at 5 ns in various solvent model simulations. On the left side, the initial structure of the simulations is also shown as a reference. The figures are generated with MOLMOL [34].

Figure 3 shows the radius of gyration (*R _{g}*) of the peptide as a function of time. As a reference, the average value of

Figure 3. Time evolutions of the radius of gyration (*R _{g}*) during 5 ns simulations using various solvent models. The average value of

Figures 4, 5 show the hydrophilic solvent-accessible surface area (*S*_{phi}) and hydrophobic solvent-accessible surface area (*S*_{pho}) of the peptide as a function of the simulation time, respectively. As a reference, the average values of *S*_{phi} and *S*_{pho} calculated from the 32 structure models determined by NMR were also shown. By the SA and GB/SA models, *S*_{phi} values decreased drastically at early stage. Especially, in the GB/SA model, this decrease in *S*_{phi} was drastic and many salt bridges were observed. On the other hand, *S*_{phi} values generated by the DD/SA and DD/SA/EC models were somewhat larger than the value calculated from NMR structures and by explicit solvent model throughout 5 ns simulation. *S*_{pho} values obtained by the SA, GB/SA models and *in vacuo* simulation were almost the same as the NMR structure value. By the DD/SA and DD/SA/EC models, *S*_{pho} values were a little bit larger than those calculated by the NMR structures but were the closest to the values obtained by explicit solvent MD simulation.

Figure 4. Time evolutions of the total hydrophilic solvent-accessible surface area (*S*_{phi}) during 5 ns simulations using various solvent models. The average values of *S*_{phi} calculated from the 32 structure models determined by NMR are shown as a reference (dashed lines).

Figure 5. Time evolutions of the total hydrophobic solvent-accessible surface area (*S*_{pho}) during 5 ns simulations using various solvent models. The average values of *S*_{pho} calculated from the 32 structure models determined by NMR are shown as a reference (dashed lines).

Finally, we examined the Ca positional root mean square fluctuations (RMSF) of the peptide during simulation using various solvent models (Figure 6). The fluctuations observed in the simulation using the explicit solvent model were large, especially around the N and C termini of the peptide regions and the hairpin turn (residues 7 and 8), while those observed in the simulation using implicit solvent models and *in vacuo* were much smaller throughout the peptide.

Figure 6. The root means square fluctuations (RMSF) of Ca atom positions obtained from simulations using various solvent models. The RMS fluctuation values were calculated as the average values over 2.5 to 5 ns period of each simulation.

Table 1. Computation Time for 1 ns Dynamics of the 28-mer bba folded peptide, pda8d.

Algorithm | Solvent model | Time(min)^{a} | Relative time |
---|---|---|---|

BD^{b} | SA | 12.7 | 1.00 |

BD^{b} | DD/SA | 12.7 | 1.00 |

BD^{b} | DD/SA/EC | 12.7 | 1.00 |

BD^{b} | GB/SA | 17.9 | 1.41 |

MD^{c} | Explicit water | 2,057 | 162.0 |

The main artifacts often observed in

The GB/SA simulation resulted in smaller

Although some artifacts encountered in

In EC model, although the one shielding parameter a

We thank Prof. Akihide Oguchi for helpful suggestions and discussions. This work was supported (I. Y.) by a Grant-in-Aid for Scientific Research on Priority Area (C) Genome Information Science from the Ministry of Education, Culture, Sports, Science and Technology of Japan.

[ 2] A. R. Fersht,

[ 3] A. Warshel and A. Papazyan,

[ 4] B. Roux and T. Simonson,

[ 5] C. J. Cramer and D. G. Truhlar,

[ 6] W. Rocchia, S. Sridharan, A. Nicholls, E. Alexov, A. Chiabrera and B. Honig,

[ 7] A. H. Elcock, M. J. Potter and J. A. McCammon, Application of Poisson-Boltzmann solvation forces to macromolecular simulations,

[ 8] W. C. Still, A. Tempczyk, R. C. Hawley and T. Hendrickson,

[ 9] B. N. Dominy and C. L. Brooks III,

[10] V. Tsui and D. A. Case,

[11] J. Zhu, Y. Shi and H. Liu,

[12] J. Vila, R. L. Williams, M. Vasquez and H. A. Scheraga,

[13] L. Wesson and D. Eisenberg,

[14] F. Fraternali and W. F. van Gunsteren,

[15] T. Lazaridis and M. Karplus,

[16] T. Ando, T. Meguro and I. Yamato,

[17] T. Ando, T. Meguro and I. Yamato,

[18] D. L. Ermak and J. A. McCammon,

[19] C. L. Brooks III, M. Karplus and B. M. Pettitt, Dynamical simulation methods,

[20] W. F. van Gunsteren, Molecular dynamics and stochastic dynamics simulation: a primer,

[21] B. R. Gelin and M. Karplus,

[22] D. Eisenberg and A. D. McLachlan,

[23] W. Hasel, T. F. Hendrickson and W. C. Still,

[24] J. B. Matthew and F. R. N. Gurd,

[25] D. Qui, P. S. Shenkin, F. P. Hollinger and W. C. Still,

[26] Y. Liu and D. L. Beveridge,

[27] G. D. Hawkins, C. J. Cramer and D. G. Truhlar,

[28] B. I. Dahiyat, C. A. Sarisky and S. L. Mayo,

[29] W. D. Cornell, P. Cieplak, C. I. Bayly, I. R. Gould, J. K. M. Merz, D. M. Ferguson, D. C. Spellmeyer, T. Fox, J. W. Caldwell, and P. A. Kollman,

[30] S. J. Weiner, P. A. Kollman, D. A. Case, U. C. Singh, C. Ghio, G. Alagona, S. J. Profeta and P. Weiner,

[31] P. E. Smith and W. F. van Gunsteren, Methods for the evaluation of long range electrostatic force in computer simulations of molecular systems,

[32] R. Zhou,

[33] F. B. Sheinerman and C. L. Brooks III,

[34] R. Koradi, M. Billeter and K. Wuthrich,

Return