Accelerating Convergence in Molecular Dynamics Simulations of Solutes in Lipid Membranes by Conducting a Random Walk along the Bilayer Normal

Abstract

All molecular dynamics simulations are susceptible to sampling errors, which degrade the accuracy and precision of observed values. The statistical convergence of simulations containing atomistic lipid bilayers is limited by the slow relaxation of the lipid phase, which can exceed hundreds of nanoseconds. These long conformational autocorrelation times are exacerbated in the presence of charged solutes, which can induce significant distortions of the bilayer structure. Such long relaxation times represent hidden barriers that induce systematic sampling errors in simulations of solute insertion. To identify optimal methods for enhancing sampling efficiency, we quantitatively evaluate convergence rates using generalized ensemble sampling algorithms in calculations of the potential of mean force for the insertion of the ionic side chain analog of arginine in a lipid bilayer. Umbrella sampling (US) is used to restrain solute insertion depth along the bilayer normal, the order parameter commonly used in simulations of molecular solutes in lipid bilayers. When US simulations are modified to conduct random walks along the bilayer normal using a Hamiltonian exchange algorithm, systematic sampling errors are eliminated more rapidly and the rate of statistical convergence of the standard free energy of binding of the solute to the lipid bilayer is increased 3-fold. We compute the ratio of the replica flux transmitted across a defined region of the order parameter to the replica flux that entered that region in Hamiltonian exchange simulations. We show that this quantity, the transmission factor, identifies sampling barriers in degrees of freedom orthogonal to the order parameter. The transmission factor is used to estimate the depth-dependent conformational autocorrelation times of the simulation system, some of which exceed the simulation time, and thereby identify solute insertion depths that are prone to systematic sampling errors and estimate the lower bound of the amount of sampling that is required to resolve these sampling errors. Finally, we extend our simulations and verify that the conformational autocorrelation times estimated by the transmission factor accurately predict correlation times that exceed the simulation time scale-something that, to our knowledge, has never before been achieved.

Permeation of the three aromatic dipeptides through lipid bilayers: Experimental and computational study

“…As Neale et al discovered, there are rare sampling barriers at the lipid interface with time scales on the order of 10 µs-far too long for our simulations to detect. 93 Whether our amino acids fall into the small or large category remains unknown and would require much longer simulation times. Finally, the increased electric field strength within the membrane may be assisting our amino acids to adopt favorable orientations on a much faster time scale due to the increased forces on the molecular dipole.…”

Sphingolipids contribute to acetic acid resistance in Zygosaccharomyces bailii

“…The distance between the acetic acid molecule and the center of the bilayer was restrained using a harmonic spring with a 1000 kJ mol −1  nm −2 force constant, and the instantaneous distance was recorded every 10 ps. The simulations with acetic acid molecules at different restrained distances were coupled, and exchange between neighbors was attempted every 4 ps to improve sampling (Neale et al, 2013). The umbrella sampling simulations were 30 ns in length at each point (simulation parameters as described above), and data from the initial 5 ns were throughout discarded as equilibration.…”

Thermodynamic and kinetic characterization of transmembrane helix association

“…Previously, slow membrane dynamics has been implicated to contribute to insufficient sampling in biased simulations along a single reaction coordinate. 42,43 However, this is of less concern in our unbiased simulations since the PMF is built by projection of the data onto the 1D reaction coordinate. In particular, the sampling of multiple association-dissociation events suggests a complete representation of the underlying coarse-grain free energy landscape, and the main features of the PMF of association can be considered converged.…”