Stiff Spring Approximation Revisited: Inertial Effects in Nonequilibrium Trajectories.

Use of harmonic guiding potentials is the most commonly adopted method for implementing steered molecular dynamics (SMD) simulations, performed to obtain potentials of mean force (PMFs) using Jarzynski's equality and other nonequilibrium work (NEW) theorems. The stiff spring approximation (SSA) of Schulten and co-workers enables calculation of the PMF by using the work performed along many SMD trajectories in NEW theorems. We discuss and demonstrate how a high spring constant, k, required for the validity of the SSA can violate another requirement of SSA, the validity of Brownian dynamics in the system under study. These result in skewed work distributions with their width increasing with k. The skew and broadening of work distributions result in biased estimation (through invoking NEW theorems) of the PMF. Remarkably, the skewness and the broadening of work distributions are independent of the average drift velocity and physical asymmetries and can only be attributed to using too-stiff springs. We discuss the proper upper limit for k such that the inertial effects are minimized. In the presence of inertial effects, using the peak value (rather than the statistical mean) of the work distributions vastly reduces the bias in the calculated PMFs and improves the accuracy.

A Permeability Study of O2 and the Trace Amine p-Tyramine through Model Phosphatidylcholine Bilayers.

We study here the permeability of the hydrophobic O2 molecule through a model DPPC bilayer at 323K and 350K, and of the trace amine p-tyramine through PC bilayers at 310K. The tyramine results are compared to previous experimental work at 298K. Nonequilibrium work methods were used in conjunction to simultaneously obtain both the potential of mean force (PMF) and the position dependent transmembrane diffusion coefficient, D(z), from the simulations. These in turn were used to calculate the permeability coefficient, P, through the inhomogeneous solubility-diffusion model. The results for O2 are consistent with previous simulations, and agree with experimentally measured P values for PC bilayers. A temperature dependence in the permeability of O2 through DPPC was obtained, with P decreasing at higher temperatures. Two relevant species of p-tyramine were simulated, from which the PMF and D(z) were calculated. The charged species had a large energetic barrier to crossing the bilayer of ~ 21 kcal/mol, while the uncharged, deprotonated species had a much lower barrier of ~ 7 kcal/mol. The effective in silico permeability for p-tyramine was calculated by applying three approximations, all of which gave nearly identical results (presented here as a function of the pKa). As the permeability value calculated from simulation was highly dependent on the pKa of the amine group, a further pKa study was performed that also varied the fraction of the uncharged and zwitterionic p-tyramine species. Using the experimental P value together with the simulated results, we were able to label the phenolic group as responsible for the pKa1 and the amine for the pKa2, that together represent all of the experimentally measured pKa values for p-tyramine. This agrees with older experimental results, in contrast to more recent work that has suggested there is a strong ambiguity in the pKa values.

Implementation of the forward-reverse method for calculating the potential of mean force using a dynamic restraining protocol.

We present a new sampling and analysis scheme for calculating the potential of mean force (PMF) of systems studied by steered molecular dynamics simulations. This scheme, which we call the bin-passing method, is based on the forward-reverse (FR) method (due to I. Kosztin and co-workers, Kosztin et al. J. Chem. Phys. 2006, 124(6), 064106) and arguments based on the second law of thermodynamics. Applying the bin-passing method results in enhanced sampling, better separation of the reversible and irreversible work distributions, and faster convergence to the underlying PMF of the system under study. Post-simulation analysis is performed using a purpose-built software that we have made publicly available at https://github.com/1particle/bin-passing_analyzer under the terms of the GNU General Public License (version 3). Three examples are provided, for systems of varying sizes and complexities, to demonstrate the efficiency of this method and the quality of the results: for the dissociation PMF of NaCl in water, the bin-passing method obtains PMFs in excellent agreement with that obtained for the same system and using the same force-field through static (equilibrium) methods. The bin-passing method gives a very symmetric PMF for passage of a single water molecule through a DPPC bilayer, and the resultant PMF leads to permeability values in better agreement with experiments than those obtained through previous simulation studies. Finally, we consider the interaction of the antimicrobial peptide HHC-36 with two model membranes and employ the bin-passing method to obtain the PMFs for peptide adsorption to the membranes. The characteristics of these PMFs are consistent with the qualities established for the HHC-36 peptide through in vivo and in vitro experiments, as a non-toxic strong antimicrobial agent.

McMillan-Mayer theory of solutions revisited: simplifications and extensions.

McMillan and Mayer (MM) proved two remarkable theorems in their paper on the equilibrium statistical mechanics of liquid solutions. They first showed that the grand canonical partition function for a solution can be reduced to one with an effectively solute-only form, by integrating out the solvent degrees of freedom. The total effective solute potential in the effective solute grand partition function can be decomposed into components which are potentials of mean force for isolated groups of one, two, three, etc., solute molecules. Second, from the first result, now assuming low solute concentration, MM derived an expansion for the osmotic pressure in powers of the solute concentration, in complete analogy with the virial expansion of gas pressure in powers of the density at low density. The molecular expressions found for the osmotic virial coefficients have exactly the same form as the corresponding gas virial coefficients, with potentials of mean force replacing vacuum potentials. In this paper, we restrict ourselves to binary liquid solutions with solute species A and solvent species B and do three things: (a) By working with a semi-grand canonical ensemble (grand with respect to solvent only) instead of the grand canonical ensemble used by MM, and avoiding graphical methods, we have greatly simplified the derivation of the first MM result, (b) by using a simple nongraphical method developed by van Kampen for gases, we have greatly simplified the derivation of the second MM result, i.e., the osmotic pressure virial expansion; as a by-product, we show the precise relation between MM theory and Widom potential distribution theory, and (c) we have extended MM theory by deriving virial expansions for other solution properties such as the enthalpy of mixing. The latter expansion is proving useful in analyzing ongoing isothermal titration calorimetry experiments with which we are involved. For the enthalpy virial expansion, we have also changed independent variables from semi-grand canonical, i.e., fixed {N(A), µ(B), V, T}, to those relevant to the experiment, i.e., fixed {N(A), N(B), p, T}, where µ denotes chemical potential, N the number of molecules, V the volume, p the pressure, and T the temperature.

Dynamic turn conformation of a short tryptophan-rich cationic antimicrobial peptide and its interaction with phospholipid membranes.

Cationic antimicrobial peptides are promising sources for novel therapeutic agents against multi-drug-resistant bacteria. HHC-36 (KRWWKWWRR) is a simple but effective antimicrobial peptide with similar or superior activity compared with several conventional antibiotics. In this biophysical study, unique conformational properties of this peptide and some of its analogs as well as its interaction with lipid membranes are investigated in detail. Circular dichroism (CD) and molecular dynamics modeling studies of HHC-36 in different environments reveal a dynamic amphipathic structure composed of competing turn conformations with free energies lower than that of the unfolded state, implying a strong influence of tryptophan interactions in formation of the turns. CD spectra and gel electrophoresis also show strong evidence of self-association of this peptide in aqueous milieu and interaction with both neutrally and negatively charged lipid membrane systems. Isothermal titration calorimetry and acrylamide fluorescence quenching experiments emphasize the preference of HHC-36 for negatively charged vesicles. In addition, dye leakage experiments suggest that this peptide functions through a surface-associated mechanism with weak lytic activity against bacterial model membranes.

Membrane permeability of trace amines: evidence for a regulated, activity-dependent, nonexocytotic, synaptic release.

Both pre- and post-synaptic effects of trace amines have been demonstrated. The putative intracellular location of Trace Amine-Associated Receptors necessitate that membrane transport processes be present in order for post-synaptic effects to occur. Here we examine the ability of trace amines to cross synthetic (Fluorosomes) and native (synaptosomes) lipid bilayer membranes. Trace amines readily crossed Fluorosome membranes by simple diffusion, p-tyramine (P = 0.01) and tryptamine (P = 0.0004) showing significantly faster diffusion than dopamine and 5-HT, respectively, with diffusion half-lives of 13.5 ± 4.1 (p-tyramine) and 6.8 ± 0.7 seconds (tryptamine). Similarly, release of [(3) H]p-tyramine and [(3) H]2-phenylethylamine from pre-loaded synaptosomes occurred significantly quicker than did [(3) H]dopamine (P = 0.0001), with half lives of 38.9 (p-tyramine), 7.8 (2-phenylethylamine) and 133.6 seconds (dopamine). This was, however, significantly slower than the diffusion mediated passage across Fluorosome membranes (P = 0.0001), suggesting a role for transporters in mediating trace amine release. Further, a pronounced shoulder region was observed in the synaptosome [(3) H]p-tyramine release curve, suggesting that multiple processes regulate release. No such shoulder region was present for [(3) H]dopamine release. Surprisingly, both [(3) H]p-tyramine (P = 0.001) and [(3) H]2-phenylethylamine (P = 0.0001) release from synaptosomes was significantly decreased under depolarizing conditions. As expected, depolarization significantly increased [(3) H]dopamine release. The data presented indicate that the release of p-tyramine and 2-phenylethylamine from neuronal terminals occurs by a different mechanism than dopamine, and does not involve classical exocytosis. The data are consistent with an initial release of trace amines by simple diffusion, followed by an activity-dependent regulation of synaptic levels via one or more transporter proteins.

Calculating diffusion and permeability coefficients with the oscillating forward-reverse method.

The forward-reverse or FR method is an efficient bidirectional work method for determining the potential of mean force w(z) and also supposedly gives in principle the position-dependent diffusion coefficient D(z). Results from a variation called the OFR (oscillating FR) method suggest inconsistencies in the D(z) values when calculated as prescribed by the FR method. A new steering protocol has thus been developed and applied to the OFR method for the accurate determination of D(z) and also provides greater convergence for w(z) in molecular dynamics simulations. The bulk diffusion coefficient for water was found to be (6.03±0.16)×10(-5) cm2/s at 350 K with system size dependence within the statistical error bars. Using this steering protocol, D(z) and w(z) for water permeating a dipalmitoylphosphatidylcholine (DPPC) bilayer were determined. The potential of mean force is shown to have a barrier of peak height, wmax/(kBT)=8.4, with a width of about 10 Å on either side from the membrane center. The diffusion constant is shown to be highest in the core region of the membrane [peak value ∼(8.0±0.8)×10(-5) cm2/s], lowest in the head-group region [minimum value ∼(2.0±0.3)×10(-5) cm2/s], and to tend toward the bulk value as the water molecule leaves the membrane. The permeability coefficient P for H2O in DPPC was determined using the simulated D(z) and w(z) to give values of (0.129±0.075) cm/s at 323 K and (0.141±0.043) cm/s at 350 K. The results show more spatial detail than results presented in previous work while reducing the computational and user effort.

Scattering density profile model of POPG bilayers as determined by molecular dynamics simulations and small-angle neutron and X-ray scattering experiments.

We combine molecular dynamics (MD) simulations and experiment, both small-angle neutron (SANS) and small-angle X-ray scattering (SAXS), to determine the precise structure of bilayers composed of 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphatidylglycerol (POPG), a lipid commonly encountered in bacterial membranes. Experiment and simulation are used to develop a one-dimensional scattering density profile (SDP) model suitable for the analysis of experimental data. The joint refinement of such data (i.e., SANS and SAXS) results in the area per lipid that is then used in the fixed-area simulations. In the final step, the direct comparison of simulated-to-experimental data gives rise to the detailed structure of POPG bilayers. From these studies we conclude that POPG's molecular area is 66.0 ± 1.3 Å(2), its overall bilayer thickness is 36.7 ± 0.7 Å, and its hydrocarbon region thickness is 27.9 ± 0.6 Å, assuming a simulated value of 1203 Å(3) for the total lipid volume.

Drift-oscillatory steering with the forward-reverse method for calculating the potential of mean force.

We present a method that enables the use of the forward-reverse (FR) method of Kosztin et al. on a broader range of problems in soft matter physics. Our method, which we call the oscillating forward-reverse (OFR) method, adds an oscillatory steering potential to the constant velocity steering potential of the FR method. This enables the calculation of the potential of mean force (PMF) in a single unidirectional oscillatory drift, rather than multiple drifts in both directions as required by the FR method. By following small forward perturbations with small reverse perturbations, the OFR method is able to generate a piecewise reverse path that follows the piecewise forward path much more closely than any practical set of paths used in the FR method. We calculate the PMF for four different systems: a dragged Brownian oscillator, a pair of atoms in a Lennard-Jones liquid, a Na(+)-Cl⁻ ion pair in an aqueous solution, and a deca-alanine molecule being stretched in an implicit solvent. In all cases, the PMF results are in good agreement with those published previously using various other methods, and, to our knowledge, we give for the first time PMFs calculated by nonequilibrium methods for the Lennard-Jones and Na(+)-Cl⁻ systems.

Molecular dynamics-based simulation of trace amine membrane permeability.

Trace amines are endogenous compounds, typified by 2-phenylethylamine (PE) and p-tyramine (TA), found in the vertebrate central nervous system. Although synthesized in pre-synaptic terminals, trace amines do not appear to act as neurotransmitters, but rather modulate responsivity to co-existing neurotransmitters. Trace amines are neither actively accumulated in synaptic vesicles, nor released in an activity-dependent manner. Further, Trace Amine-Associated Receptor 1 (TAAR1), which is selectively activated by PE and TA, is intracellular. As such, PE and TA need to cross cell membranes in order to exert their effects. This has been assumed to occur by simple diffusion, but has not previously been systematically examined. Experimental data were obtained using Fluorosome(®) technology. A permeability coefficient of 25.3 ± 3.8 Å/s (n = 6) was obtained for TA which was not significantly different from that obtained for the monoamine neurotransmitter noradrenaline (20.3 ± 3.8 Å/s, n = 8). PE was unsuitable for use with this system. We have also used molecular dynamics computer simulation techniques to determine the potential of mean force (PMF) associated with trace amine passage across lipid bilayers. A PMF peak barrier of 25 ± 6 kcal/mol (protonated) and 13 ± 1 kcal/mol (deprotonated) was obtained for PE. Protonated TA peak energy barriers were even greater at 31 ± 1 kcal/mol. Application of a homogeneous solubility-diffusion model combining the measured permeability coefficients and simulated PMF allows fitting of the diffusion coefficient for trace amines in the hydrophobic region of the lipid bilayer. The diffusion coefficients in other regions of the membrane were found to make negligible contributions to the permeability coefficient for the calculated PMF. The fit obtained a value for the diffusion coefficient of (163 ± 25) × 10(-10) m(2)/s for TA(+) in the hydrophobic core region. The diffusion coefficient for TA(+) in the aqueous compartment was also calculated directly by simulation yielding a value of (0.62 ± 0.26) × 10(-10) m(2)/s. The adopted simulation methods failed to yield diffusion coefficients in the core region indicating that further work will be required to accurately predict permeability coefficients for trace amines passing through membranes.

Binding free energy and counterion release for adsorption of the antimicrobial peptide lactoferricin B on a POPG membrane.

Molecular dynamics (MD) simulations are used to study the interaction of an anionic palmitoyl-oleoyl-phosphatidylglycerol (POPG) bilayer with the cationic antimicrobial peptide bovine lactoferricin (LFCinB) in a 100 mM NaCl solution at 310 K. The interaction of LFCinB with a POPG bilayer is employed as a model system for studying the details of membrane adsorption selectivity of cationic antimicrobial peptides. Seventy eight 4 ns MD production run trajectories of the equilibrated system, with six restrained orientations of LFCinB at 13 different separations from the POPG membrane, are generated to determine the free energy profile for the peptide as a function of the distance between LFCinB and the membrane surface. To calculate the profile for this relatively large system, a variant of constrained MD and thermodynamic integration is used. A simplified method for relating the free energy profile to the LFCinB-POPG membrane binding constant is employed to predict a free energy of adsorption of -5.4+/-1.3 kcal/mol and a corresponding maximum adsorption binding force of about 58 pN. We analyze the results using Poisson-Boltzmann theory. We find the peptide-membrane attraction to be dominated by the entropy increase due to the release of counterions and polarized water from the region between the charged membrane and peptide, as the two approach each other. We contrast these results with those found earlier for adsorption of LFCinB on the mammalianlike palmitoyl-oleoyl-phosphatidylcholine membrane.

Prediction of binding free energy for adsorption of antimicrobial peptide lactoferricin B on a POPC membrane.

Molecular dynamics (MD) simulations are used to study the interaction of a zwitterionic palmitoyl-oleoyl-phosphatidylcholine (POPC) bilayer with the cationic antimicrobial peptide bovine lactoferricin (LFCinB) in a 100 mM NaCl solution at 310 K. The interaction of LFCinB with POPC is used as a model system for studying the details of membrane-peptide interactions, with the peptide selected because of its antimicrobial nature. Seventy-two 3 ns MD simulations, with six orientations of LFCinB at 12 different distances from a POPC membrane, are carried out to determine the potential of mean force (PMF) or free energy profile for the peptide as a function of the distance between LFCinB and the membrane surface. To calculate the PMF for this relatively large system a new variant of constrained MD and thermodynamic integration is developed. A simplified method for relating the PMF to the LFCinB-membrane binding free energy is described and used to predict a free energy of adsorption (or binding) of -1.05+/-0.39 kcal/mol , and corresponding maximum binding force of about 20 pN, for LFCinB-POPC. The contributions of the ions-LFCinB and the water-LFCinB interactions to the PMF are discussed. The method developed will be a useful starting point for future work simulating peptides interacting with charged membranes and interactions involved in the penetration of membranes, features necessary to understand in order to rationally design peptides as potential alternatives to traditional antibiotics.