- •Foreword
- •Preface
- •Contents
- •Introduction
- •Oren M. Becker
- •Alexander D. MacKerell, Jr.
- •Masakatsu Watanabe*
- •III. SCOPE OF THE BOOK
- •IV. TOWARD A NEW ERA
- •REFERENCES
- •Atomistic Models and Force Fields
- •Alexander D. MacKerell, Jr.
- •II. POTENTIAL ENERGY FUNCTIONS
- •D. Alternatives to the Potential Energy Function
- •III. EMPIRICAL FORCE FIELDS
- •A. From Potential Energy Functions to Force Fields
- •B. Overview of Available Force Fields
- •C. Free Energy Force Fields
- •D. Applicability of Force Fields
- •IV. DEVELOPMENT OF EMPIRICAL FORCE FIELDS
- •B. Optimization Procedures Used in Empirical Force Fields
- •D. Use of Quantum Mechanical Results as Target Data
- •VI. CONCLUSION
- •REFERENCES
- •Dynamics Methods
- •Oren M. Becker
- •Masakatsu Watanabe*
- •II. TYPES OF MOTIONS
- •IV. NEWTONIAN MOLECULAR DYNAMICS
- •A. Newton’s Equation of Motion
- •C. Molecular Dynamics: Computational Algorithms
- •A. Assigning Initial Values
- •B. Selecting the Integration Time Step
- •C. Stability of Integration
- •VI. ANALYSIS OF DYNAMIC TRAJECTORIES
- •B. Averages and Fluctuations
- •C. Correlation Functions
- •D. Potential of Mean Force
- •VII. OTHER MD SIMULATION APPROACHES
- •A. Stochastic Dynamics
- •B. Brownian Dynamics
- •VIII. ADVANCED SIMULATION TECHNIQUES
- •A. Constrained Dynamics
- •C. Other Approaches and Future Direction
- •REFERENCES
- •Conformational Analysis
- •Oren M. Becker
- •II. CONFORMATION SAMPLING
- •A. High Temperature Molecular Dynamics
- •B. Monte Carlo Simulations
- •C. Genetic Algorithms
- •D. Other Search Methods
- •III. CONFORMATION OPTIMIZATION
- •A. Minimization
- •B. Simulated Annealing
- •IV. CONFORMATIONAL ANALYSIS
- •A. Similarity Measures
- •B. Cluster Analysis
- •C. Principal Component Analysis
- •REFERENCES
- •Thomas A. Darden
- •II. CONTINUUM BOUNDARY CONDITIONS
- •III. FINITE BOUNDARY CONDITIONS
- •IV. PERIODIC BOUNDARY CONDITIONS
- •REFERENCES
- •Internal Coordinate Simulation Method
- •Alexey K. Mazur
- •II. INTERNAL AND CARTESIAN COORDINATES
- •III. PRINCIPLES OF MODELING WITH INTERNAL COORDINATES
- •B. Energy Gradients
- •IV. INTERNAL COORDINATE MOLECULAR DYNAMICS
- •A. Main Problems and Historical Perspective
- •B. Dynamics of Molecular Trees
- •C. Simulation of Flexible Rings
- •A. Time Step Limitations
- •B. Standard Geometry Versus Unconstrained Simulations
- •VI. CONCLUDING REMARKS
- •REFERENCES
- •Implicit Solvent Models
- •II. BASIC FORMULATION OF IMPLICIT SOLVENT
- •A. The Potential of Mean Force
- •III. DECOMPOSITION OF THE FREE ENERGY
- •A. Nonpolar Free Energy Contribution
- •B. Electrostatic Free Energy Contribution
- •IV. CLASSICAL CONTINUUM ELECTROSTATICS
- •A. The Poisson Equation for Macroscopic Media
- •B. Electrostatic Forces and Analytic Gradients
- •C. Treatment of Ionic Strength
- •A. Statistical Mechanical Integral Equations
- •VI. SUMMARY
- •REFERENCES
- •Steven Hayward
- •II. NORMAL MODE ANALYSIS IN CARTESIAN COORDINATE SPACE
- •B. Normal Mode Analysis in Dihedral Angle Space
- •C. Approximate Methods
- •IV. NORMAL MODE REFINEMENT
- •C. Validity of the Concept of a Normal Mode Important Subspace
- •A. The Solvent Effect
- •B. Anharmonicity and Normal Mode Analysis
- •VI. CONCLUSIONS
- •ACKNOWLEDGMENT
- •REFERENCES
- •Free Energy Calculations
- •Thomas Simonson
- •II. GENERAL BACKGROUND
- •A. Thermodynamic Cycles for Solvation and Binding
- •B. Thermodynamic Perturbation Theory
- •D. Other Thermodynamic Functions
- •E. Free Energy Component Analysis
- •III. STANDARD BINDING FREE ENERGIES
- •IV. CONFORMATIONAL FREE ENERGIES
- •A. Conformational Restraints or Umbrella Sampling
- •B. Weighted Histogram Analysis Method
- •C. Conformational Constraints
- •A. Dielectric Reaction Field Approaches
- •B. Lattice Summation Methods
- •VI. IMPROVING SAMPLING
- •A. Multisubstate Approaches
- •B. Umbrella Sampling
- •C. Moving Along
- •VII. PERSPECTIVES
- •REFERENCES
- •John E. Straub
- •B. Phenomenological Rate Equations
- •II. TRANSITION STATE THEORY
- •A. Building the TST Rate Constant
- •B. Some Details
- •C. Computing the TST Rate Constant
- •III. CORRECTIONS TO TRANSITION STATE THEORY
- •A. Computing Using the Reactive Flux Method
- •B. How Dynamic Recrossings Lower the Rate Constant
- •IV. FINDING GOOD REACTION COORDINATES
- •A. Variational Methods for Computing Reaction Paths
- •B. Choice of a Differential Cost Function
- •C. Diffusional Paths
- •VI. HOW TO CONSTRUCT A REACTION PATH
- •A. The Use of Constraints and Restraints
- •B. Variationally Optimizing the Cost Function
- •VII. FOCAL METHODS FOR REFINING TRANSITION STATES
- •VIII. HEURISTIC METHODS
- •IX. SUMMARY
- •ACKNOWLEDGMENT
- •REFERENCES
- •Paul D. Lyne
- •Owen A. Walsh
- •II. BACKGROUND
- •III. APPLICATIONS
- •A. Triosephosphate Isomerase
- •B. Bovine Protein Tyrosine Phosphate
- •C. Citrate Synthase
- •IV. CONCLUSIONS
- •ACKNOWLEDGMENT
- •REFERENCES
- •Jeremy C. Smith
- •III. SCATTERING BY CRYSTALS
- •IV. NEUTRON SCATTERING
- •A. Coherent Inelastic Neutron Scattering
- •B. Incoherent Neutron Scattering
- •REFERENCES
- •Michael Nilges
- •II. EXPERIMENTAL DATA
- •A. Deriving Conformational Restraints from NMR Data
- •B. Distance Restraints
- •C. The Hybrid Energy Approach
- •III. MINIMIZATION PROCEDURES
- •A. Metric Matrix Distance Geometry
- •B. Molecular Dynamics Simulated Annealing
- •C. Folding Random Structures by Simulated Annealing
- •IV. AUTOMATED INTERPRETATION OF NOE SPECTRA
- •B. Automated Assignment of Ambiguities in the NOE Data
- •C. Iterative Explicit NOE Assignment
- •D. Symmetrical Oligomers
- •VI. INFLUENCE OF INTERNAL DYNAMICS ON THE
- •EXPERIMENTAL DATA
- •VII. STRUCTURE QUALITY AND ENERGY PARAMETERS
- •VIII. RECENT APPLICATIONS
- •REFERENCES
- •II. STEPS IN COMPARATIVE MODELING
- •C. Model Building
- •D. Loop Modeling
- •E. Side Chain Modeling
- •III. AB INITIO PROTEIN STRUCTURE MODELING METHODS
- •IV. ERRORS IN COMPARATIVE MODELS
- •VI. APPLICATIONS OF COMPARATIVE MODELING
- •VII. COMPARATIVE MODELING IN STRUCTURAL GENOMICS
- •VIII. CONCLUSION
- •ACKNOWLEDGMENTS
- •REFERENCES
- •Roland L. Dunbrack, Jr.
- •II. BAYESIAN STATISTICS
- •A. Bayesian Probability Theory
- •B. Bayesian Parameter Estimation
- •C. Frequentist Probability Theory
- •D. Bayesian Methods Are Superior to Frequentist Methods
- •F. Simulation via Markov Chain Monte Carlo Methods
- •III. APPLICATIONS IN MOLECULAR BIOLOGY
- •B. Bayesian Sequence Alignment
- •IV. APPLICATIONS IN STRUCTURAL BIOLOGY
- •A. Secondary Structure and Surface Accessibility
- •ACKNOWLEDGMENTS
- •REFERENCES
- •Computer Aided Drug Design
- •Alexander Tropsha and Weifan Zheng
- •IV. SUMMARY AND CONCLUSIONS
- •REFERENCES
- •Oren M. Becker
- •II. SIMPLE MODELS
- •III. LATTICE MODELS
- •B. Mapping Atomistic Energy Landscapes
- •C. Mapping Atomistic Free Energy Landscapes
- •VI. SUMMARY
- •REFERENCES
- •Toshiko Ichiye
- •II. ELECTRON TRANSFER PROPERTIES
- •B. Potential Energy Parameters
- •IV. REDOX POTENTIALS
- •A. Calculation of the Energy Change of the Redox Site
- •B. Calculation of the Energy Changes of the Protein
- •B. Calculation of Differences in the Energy Change of the Protein
- •VI. ELECTRON TRANSFER RATES
- •A. Theory
- •B. Application
- •REFERENCES
- •Fumio Hirata and Hirofumi Sato
- •Shigeki Kato
- •A. Continuum Model
- •B. Simulations
- •C. Reference Interaction Site Model
- •A. Molecular Polarization in Neat Water*
- •B. Autoionization of Water*
- •C. Solvatochromism*
- •F. Tautomerization in Formamide*
- •IV. SUMMARY AND PROSPECTS
- •ACKNOWLEDGMENTS
- •REFERENCES
- •Nucleic Acid Simulations
- •Alexander D. MacKerell, Jr.
- •Lennart Nilsson
- •D. DNA Phase Transitions
- •III. METHODOLOGICAL CONSIDERATIONS
- •A. Atomistic Models
- •B. Alternative Models
- •IV. PRACTICAL CONSIDERATIONS
- •A. Starting Structures
- •C. Production MD Simulation
- •D. Convergence of MD Simulations
- •WEB SITES OF INTEREST
- •REFERENCES
- •Membrane Simulations
- •Douglas J. Tobias
- •II. MOLECULAR DYNAMICS SIMULATIONS OF MEMBRANES
- •B. Force Fields
- •C. Ensembles
- •D. Time Scales
- •III. LIPID BILAYER STRUCTURE
- •A. Overall Bilayer Structure
- •C. Solvation of the Lipid Polar Groups
- •IV. MOLECULAR DYNAMICS IN MEMBRANES
- •A. Overview of Dynamic Processes in Membranes
- •B. Qualitative Picture on the 100 ps Time Scale
- •C. Incoherent Neutron Scattering Measurements of Lipid Dynamics
- •F. Hydrocarbon Chain Dynamics
- •ACKNOWLEDGMENTS
- •REFERENCES
- •Appendix: Useful Internet Resources
- •B. Molecular Modeling and Simulation Packages
- •Index
44 |
Becker and Watanabe |
C. Molecular Dynamics: Computational Algorithms
Solving Newton’s equation of motion requires a numerical procedure for integrating the differential equation. A standard method for solving ordinary differential equations, such as Newton’s equation of motion, is the finite-difference approach. In this approach, the molecular coordinates and velocities at a time t ∆t are obtained (to a sufficient degree of accuracy) from the molecular coordinates and velocities at an earlier time t. The equations are solved on a step-by-step basis. The choice of time interval ∆t depends on the properties of the molecular system simulated, and ∆t must be significantly smaller than the characteristic time of the motion studied (Section V.B).
A good starting point for understanding finite-difference methods is the Taylor expansion about time t of the position at time t ∆t,
r(t ∆t) r(t) r˙(t)∆t |
1 |
r¨(t)∆t2 |
(11) |
||
|
|
||||
2 |
|
|
|||
Alternatively, this can be written as |
|
|
|||
r(t ∆t) r(t) v(t)∆t |
1 |
a(t)∆t2 |
(12) |
||
|
|||||
2 |
|
|
where v(t) is the velocity vector and a(t) is the acceleration. Because the integration proceeds in a stepwise fashion, and recalling Eq. (6), it is convenient to rewrite the above expansion in a discrete form. Using rn to indicate the position at step n (at time t) and rn 1 to indicate the position at the next step, n 1 (at time t ∆t), Eq. (12) can be written as
rn 1 rn vn ∆t |
1 |
|
Fn |
∆t2 O(∆t3) |
(13) |
2 |
m |
where O(∆tn) is the terms of order ∆tn or smaller. With this information the velocity vn 1 at time n 1 can be crudely estimated, for example, as
vn 1 (rn 1 rn)/2 |
(14) |
Together, Eqs. (13) and (14) form an integration algorithm. Given the position rn, the velocity vn, and the force Fn at step n, these equations allow one to calculate (actually, estimate) the position rn 1 and velocity vn 1 at step n 1. The formulation is highly trivial and results in a low quality integration algorithm (large errors). Other, more accurate, algorithms have been developed using the same kind of reasoning. In the following subsections we survey some of the more commonly used finite-difference integration algorithms, highlighting their advantages and disadvantages.
1. Verlet Integrator
The most common integration algorithm used in the study of biomolecules is due to Verlet [11]. The Verlet integrator is based on two Taylor expansions, a forward expansion (t ∆t) and a backward expansion (t ∆t),
rn 1 rn vn ∆t |
1 |
|
Fn |
∆t2 O(∆t3) |
(15a) |
2 |
m |
Dynamics Methods |
|
|
|
|
|
|
45 |
|
rn 1 |
rn vn ∆t |
1 |
|
Fn |
∆t2 O(∆t3) |
(15b) |
||
2 |
m |
|||||||
The sum of the two expansions yields an algorithm for propagating the position, |
|
|||||||
rn 1 |
2rn rn 1 |
Fn |
∆t2 O(∆t4) |
(16) |
||||
|
||||||||
|
|
|
m |
|
Translated into a stream of commands, this algorithm is executed in two steps:
1.Use the current position rn to calculate the current force Fn.
2.Use the current and previous positions rn and rn 1 together with the current force Fn (calculated in step 1) to calculate the position in the next step, rn 1, according to Eq. (16).
These two steps are repeated for every time step for each atom in the molecule. Subtracting Eq. (15b) from Eq (15a) yields a complementary algorithm for propagating the velocities,
vn |
rn 1 rn 1 |
O(∆t2) |
(17) |
|
|||
|
2∆t |
|
Figure 1a gives a graphical representation of the steps involved in a Verlet propagation. The algorithm embodied in Eqs. (16) and (17) provides a stable numerical method for solving Newton’s equation of motion for systems ranging from simple fluids to biopolymers. Like any algorithm, the Verlet algorithm has advantages as well as disadvantages.
Figure 1 A stepwise view of the Verlet integration algorithm and its variants. (a) The basic Verlet method. (b) Leap-frog integration. (c) Velocity Verlet integration. At each algorithm dark and light gray cells indicate the initial and calculating variables, respectively. The numbers in the cells represent the orders in the calculation procedures. The arrows point from the data that are used in the calculation of the variable that is being calculated at each step.
46 |
Becker and Watanabe |
Advantages of the Verlet algorithm are
1.The position integration is quite accurate [errors on the order of O(∆t4)] and is independent of the velocity propagation, a fact that simplifies the position integration and reduces memory requirements.
2.The algorithm requires only a single force evaluation per integration cycle (computationally, force evaluations are the most ‘‘expensive’’ part of the simulation).
3.This formulation, which is based on forward and backward expansions, guarantees time reversibility (a property of the equation of motion).
Disadvantages of the Verlet algorithm are
1. The velocity propagation is subject to relatively large errors, on the order O(∆t2). Recall that an accurate estimate of the velocity is required for the kinetic energy evaluations. An added inconvenience is that vn can be computed only if rn 1 is already known.
2.Further numerical errors are introduced when an O(∆t2) term is added to an O(∆t0) term.
3.The Verlet algorithm is not ‘‘self-starting.’’ A lower order Taylor expansion [e.g., Eq. (13)] is often used to initiate the propagation.
4.It must be modified to incorporate velocity-dependent forces or temperature scaling.
2.Leap-Frog Integrator
Modifications to the basic Verlet scheme have been proposed to tackle the above deficiencies, particularly to improve the velocity evaluation. One of these modifications is the leap-frog algorithm, so called for its half-step scheme: Velocities are evaluated at the midpoint of the position evaluation and vice versa [12,13]. The algorithm can be written as
rn 1 rn vn 1/2 ∆t |
(18a) |
||||
vn 1/2 |
vn 1/2 |
|
Fn |
∆t |
(18b) |
|
|||||
|
|
|
m |
|
where vn 1/2 stands for the velocity at the mid-step time [t (1/2)∆t]. Elimination of the velocities from these equations shows that the method is algebraically equivalent to the Verlet algorithm. Cast in the form of execution instructions, the leap-frog algorithm involves three steps:
1.Use the current position rn to calculate the current force Fn.
2.Use the current force Fn and previous mid-step velocity vn 1/2 to calculate the next mid-step velocity vn 1/2.
3.Use the current position rn and the next mid-step velocity vn 1/2 (from step 2) to calculate the position in the next step, rn 1.
Figure 1b gives a graphical representation of the steps involved in the leap-frog propagation. The current velocity vn, which is necessary for calculating the kinetic energy, can be calculated as
vn (vn 1/2 vn 1/2)/2 |
(19) |
Dynamics Methods |
47 |
Advantages of the leap-frog algorithm are
1.It improves velocity evaluation.
2.The direct evaluation of velocities gives a useful handle for controlling the simulation temperature (via velocity scaling).
3.It reduces the numerical errors, since here O(∆t1) terms are added to O(∆t0) terms.
Disadvantages of the leap-frog algorithm are
1.It still does not handle the velocities in a completely satisfactory manner, because the velocities at time t are only approximated by Eq. (19).
2.This algorithm is computationally a little more expensive than the Verlet algorithm.
3.Velocity Verlet Integrator
An even better handling of the velocities is obtained by another variant of the basic Verlet integrator, known as the ‘‘velocity Verlet’’ algorithm. This is a Verlet-type algorithm that stores positions, velocities, and accelerations all at the same time t and minimizes roundoff errors [14]. The velocity Verlet algorithm is written
rn 1 |
rn vn ∆t |
1 |
|
Fn |
∆t2 |
(20a) |
||||||
2 |
m |
|||||||||||
vn 1 |
vn |
1 |
|
Fn |
|
Fn 1 |
∆t |
(20b) |
||||
2 |
m |
m |
Again, elimination of the velocities from these equations recovers the Verlet algorithm. In practice, the velocity Verlet algorithm consists of the following steps:
1.Calculate the position rn 1 at time t ∆t from Eq. (20a).
2.Calculate the velocity at mid-step vn 1/2 using the equation
vn 1/2 vn |
1 |
|
Fn |
∆t |
(21) |
2 |
m |
3.Calculate the force Fn 1 at time t ∆t.
4.Finally, complete the velocity move to vn by using
vn 1 vn 1/2 |
1 |
|
Fn 1 |
∆t |
(22) |
|
|||||
2 |
m |
At this point, the kinetic energy at time t ∆t is available. Figure 1c gives a graphical representation of the steps involved in the velocity Verlet propagation.
Advantages of the velocity Verlet algorithm are
1.It is numerically very stable.
2.It is convenient and simple [the code of this method is a straightforward transcription of Eqs. (20)–(22)].
3.It provides an accurate evaluation of velocities and hence of the kinetic energy.