Role of lipid charge in organization of water/lipid bilayer interface: insights via computer simulations

Anionic unsaturated lipid bilayers represent suitable model systems that mimic real cell membranes: they are fluid and possess a negative surface charge. Understanding of detailed molecular organization of water-lipid interfaces in such systems may provide an important insight into the mechanisms of proteins' binding to membranes. Molecular dynamics (MD) of full-atom hydrated lipid bilayers is one of the most powerful tools to address this problem in silico. Unfortunately, wide application of computational methods for such systems is limited by serious technical problems. They are mainly related to correct treatment of long-range electrostatic effects. In this study a physically reliable model of an anionic unsaturated bilayer of 1,2-dioleoyl-sn-glycero-3-phosphoserine (DOPS) was elaborated and subjected to long-term MD simulations. Electrostatic interactions were treated with two different algorithms: spherical cutoff function and particle-mesh Ewald summation (PME). To understand the role of lipid charge in the system behavior, similar calculations were also carried out for zwitterionic bilayer composed of 1,2-dioleoyl-sn-glycero-3-phosphocholine (DOPC). It was shown that, for the charged DOPS bilayer, the PME protocol performs much better than the cutoff scheme. In the last case a number of artifacts in the structural organization of the bilayer were observed. All of them were attributed to inadequate treatment of electrostatic interactions of lipid headgroups with counterions. Electrostatic properties, along with structural and dynamic parameters, of both lipid bilayers were investigated. Comparative analysis of the MD data reveals that the water-lipid interface of the DOPC bilayer is looser than that for DOPS. This makes possible deeper penetration of water molecules inside the zwitterionic (DOPC) bilayer, where they strongly interact with carbonyls of lipids. This can lead to thickening of the membrane interface in zwitterionic as compared to negatively charged bilayers.     

Molecular dynamics simulation of nonsteroidal antiinflammatory drugs,naproxen and relafen,in a lipid bilayer membrane

Naproxen and relafen, as nonsteroidal antiinflammatory drugs, were simulated in neutral and charged forms and their effects on a lipid bilayer membrane were investigated by molecular dynamics simulation using Groningen machine for chemical simulations software (GROMACS). Simulation of 10 systems was performed, which included different dosages of the drug molecules, naproxen and Relafen, in charged and neutral forms, and a mixture of naproxen and Relafen in neutral forms. The effects of the mixture and the individual drugs' dosages on membrane properties, such as electrostatic potential, order parameter, diffusion coefficients, and hydrogen bond formation, were analyzed. Hydration of the drugs in the membrane system was investigated using radial distribution function analysis. Using fully hydrated dimyristoylphosphatidylcholine (DMPC) as a reference system, 128 lipid molecules and water molecules were simulated exclusively, and the same simulation technique was performed on 10 other systems, including drug mixtures and a DMPC membrane. Angular distributions of lipid chains of the membrane were calculated, and the effects of the drug insertion and chain orientation in the membrane were evaluated. © 2013 Wiley Periodicals, Inc.     

Combined Monte Carlo and molecular dynamics simulation of hydrated 18:0 sphingomyelin-cholesterol lipid bilayers

We have carried out atomic level molecular dynamics and Monte Carlo simulations of hydrated 18:0 sphingomyelin (SM)-cholesterol (CHOL) bilayers at temperatures of 20 and 50 degrees C. The simulated systems each contained 266 SM, 122 CHOL, and 11861 water molecules. Each simulation was run for 10 ns under semi-isotropic pressure boundary conditions. The particle-mesh Ewald method was used for long-range electrostatic interactions. Properties of the systems were calculated over the final 3 ns. We compare the properties of 20 and 50 degrees C bilayer systems with each other, with experimental data, and with experimental and simulated properties of pure SM bilayers and dipalmitoyl phospatidyl choline (DPPC)-CHOL bilayers. The simulations reveal an overall similarity of both systems, despite the 30 degrees C temperature difference which brackets the pure SM main phase transition. The area per molecule, lipid chain order parameter profiles, atom distributions, and electron density profiles are all very similar for the two simulated systems. Consistent with simulations from our lab and others, we find strong intramolecular hydrogen bonding in SM molecules between the phosphate ester oxygen and the hydroxyl hydrogen atoms. We also find that cholesterol hydroxyl groups tend to form hydrogen bonds primarily with SM carbonyl, methyl, and amide moieties and to a lesser extent methyl and hydroxyl oxygens.     

Enhanced Lipid Diffusion and Mixing in Accelerated Molecular Dynamics

Accelerated molecular dynamics (aMD) is an enhanced sampling technique that expedites conformational space sampling by reducing the barriers separating various low-energy states of a system. Here, we present the first application of the aMD method on lipid membranes. Altogether, ～1.5 μs simulations were performed on three systems: a pure POPC bilayer, a pure DMPC bilayer, and a mixed POPC:DMPC bilayer. Overall, the aMD simulations are found to produce significant speedup in trans-gauche isomerization and lipid lateral diffusion versus those in conventional MD (cMD) simulations. Further comparison of a 70-ns aMD run and a 300-ns cMD run of the mixed POPC:DMPC bilayer shows that the two simulations yield similar lipid mixing behaviors, with aMD generating a 2-3-fold speedup compared to cMD. Our results demonstrate that the aMD method is an efficient approach for the study of bilayer structural and dynamic properties. On the basis of simulations of the three bilayer systems, we also discuss the impact of aMD parameters on various lipid properties, which can be used as a guideline for future aMD simulations of membrane systems.     

Force field dependence of phospholipid headgroup and acyl chain properties: Comparative molecular dynamics simulations of DMPC bilayers

The reliability of molecular simulations largely depends on the quality of the empirical force field parameters. Force fields used in lipid simulations continue to be improved to enhance the agreement with experiments for a number of different properties. In this work, we have carried out molecular dynamics simulations of neat DMPC bilayers using united‐atom Berger force field and three versions of all‐atom CHARMM force fields. Three different systems consisting of 48, 72, and 96 lipids were studied. Both particle mesh Ewald (PME) and spherical cut‐off schemes were used to evaluate the long‐range electrostatic interactions. In total, 21 simulations were carried out and analyzed to find out the dependence of lipid properties on force fields, system size, and schemes to calculate long‐range interactions. The acyl chain order parameters calculated from Berger and the recent versions of CHARMM simulations have shown generally good agreement with the experimental results. However, both sets of force fields deviate significantly from the experimentally observed P‐C dipolar coupling values for the carbon atoms that link the choline and glycerol groups with the phosphate groups. Significant differences are also observed in several headgroup parameters between CHARMM and Berger simulations. Our results demonstrate that when changes were introduced to improve CHARMM force field using PME scheme, all the headgroup parameters have not been reoptimized. The headgroup properties are likely to play a significant role in lipid–lipid, protein–lipid, and ligand–lipid interactions and hence headgroup parameters in phospholipids require refinement for both Berger and CHARMM force fields. © 2009 Wiley Periodicals, Inc.J Comput Chem, 2010     

Study of frictional properties of a phospholipid bilayer in a liquid environment with lateral force microscopy as a function of NaCl concentration

Friction properties of 1,2-dimyristoyl-sn-glycero-3-phosphocholine (DMPC)-supported planar bilayers deposited on mica were tested in a liquid environment by lateral force microscopy. The presence of these bilayers was detected by imaging and force measurements with atomic force microscopy. To test how the presence of NaCl affects the frictional properties of the phospholipid bilayers, four DMPC bilayers were prepared on mica in saline media ranging from 0 to 0.1 M NaCl. Changes in the lateral vs vertical force curves were recorded as a function of NaCl concentration and related to structural changes induced in the DMPC bilayer by electrolyte ions. Three friction regimes were observed as the vertical force exerted by the tip on the bilayer increased. To relate the friction response to the structure of the DMPC bilayer, topographic images were recorded at the same time as friction data. Ions in solution screened charges present in DMPC polar heads, leading to more compact bilayers. As a consequence, the vertical force at which the bilayer broke during friction experiments increased with NaCl concentration. In addition, the topographic images showed that low-NaCl-concentration bilayers recover more easily due to the low cohesion between phospholipid molecules.     

Effect of different treatments of long-range interactions and sampling conditions in molecular dynamic simulations of rhodopsin embedded in a dipalmitoyl phosphatidylcholine bilayer

The present study analyzes the effect of the simulation conditions on the results of molecular dynamics simulations of G-protein coupled receptors (GPCRs) performed with an explicit lipid bilayer. Accordingly, the present work reports the analysis of different simulations of bovine rhodopsin embedded in a dipalmitoyl phosphatidylcholine (DPPC) lipid bilayer using two different sampling conditions and two different approaches for the treatment of long-range electrostatic interactions. Specifically, sampling was carried out either by using the statistical ensembles NVT or NPT (constant number of atoms, a pressure of 1 atm in all directions and fixed temperature), and the electrostatic interactions were treated either by using a twin-cutoff, or the particle mesh Ewald summation method (PME). The results of the present study suggest that the use of the NPT ensemble in combination with the PME method provide more realistic simulations. The use of NPT during the equilibration avoids the need of an a priori estimation of the box dimensions, giving the correct area per lipid. However, once the system is equilibrated, the simulations are irrespective of the sampling conditions used. The use of an electrostatic cutoff induces artifacts on both lipid thickness and the ion distribution, but has no direct effect on the protein and water molecules.     

Nonintercalating nanosubstrates create asymmetry between bilayer leaflets

The physical properties of lipid bilayers can be remodeled by a variety of environmental factors. Here we investigate using molecular dynamics simulations the specific effects of nanoscopic substrates or external contact points on lipid membranes. We expose palmitoyl-oleoyl phosphatidylcholine bilayers unilaterally and separately to various model nanosized substrates differing in surface hydroxyl densities. We find that a surface hydroxyl density as low as 10% is sufficient to keep the bilayer juxtaposed to the substrate. The bilayer interacts with the substrate indirectly through multiple layers of water molecules; however, despite such buffered interaction, the bilayers exhibit certain properties different from unsupported bilayers. The substrates modify transverse lipid fluctuations, charge density profiles, and lipid diffusion rates, although differently in the two leaflets, which creates an asymmetry between bilayer leaflets. Other properties that include lipid cross-sectional areas, component volumes, and order parameters are minimally affected. The extent of asymmetry that we observe between bilayer leaflets is well beyond what has been reported for bilayers adsorbed on infinite solid supports. This is perhaps because the bilayers are much closer to our nanosized finite supports than to infinite solid supports, resulting in a stronger support-bilayer electrostatic coupling. The exposure of membranes to nanoscopic contact points, therefore, cannot be considered as a simple linear interpolation between unsupported membranes and membranes supported on infinite supports. In the biological context, this suggests that the exposure of membranes to nonintercalating proteins, such as those belonging to the cytoskeleton, should not always be considered as passive nonconsequential interactions.     

Interface water dynamics and porating electric fields for phospholipid bilayers

Lipid bilayers, normally a barrier to charged species and large molecules, are permeabilized by electric fields, a phenomenon exploited by cell biologists and geneticists for porating and transfecting cells and tissues. Recent molecular simulation studies have advanced our understanding of electroporation, but the relative contributions of atomically local details (interface water and headgroup dipole and counterion configurations) and medium- and long-range electrostatic gradients and changes in membrane structural shifts to the initiating conditions and mechanisms of pore formation remain unclear. Molecular dynamics simulations of electroporation in several lipid systems presented here reveal the effects of lipid hydrocarbon tail length and composition on the magnitude of the field required for poration and on the location of the initial sites of field-driven water intrusion into the bilayer. Minimum porating external fields of 260 mV nm(-1), 280 mV nm(-1), 320 mV nm(-1), and 380 mV nm(-1) were found for 1,2-dilauroyl-sn-glycero-3-phosphatidylcholine (DLPC), 1,2-dipalmitoyl-sn-glycero-3-phosphatidylcholine (DPPC), 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphatidylcholine (POPC), and 1,2-dioleoyl- sn-glycero-3-phosphatidylcholine (DOPC) bilayers, respectively, and correlated most strongly with the bilayer thickness. These phospholipid systems share several common features including a wide, dynamic distribution of the headgroup dipole angle with the bilayer normal ranging from 0 to 155 degrees that is only slightly shifted in a porating electric field, and similar electric field-induced shifts in water dipole orientation, although the mean water dipole moment profile at the aqueous-membrane interface is more sensitive to the electric field for DOPC than for the other phospholipids. The location of pore initiation, at the anode- or cathode-facing leaflet, varies with the composition of the bilayer and correlates with a change in the polarity of the localized electric field at the interface.     

Structure and dynamics of water at the interface with phospholipid bilayers

We have performed two molecular-dynamics simulations to study the structural and dynamical properties of water at the interface with phospholipid bilayers. In one of the simulations the bilayer contained neutral phospholipid molecules, dioleoylphosphatidylcholine (DOPC); in the second simulation the bilayer contained charged lipid molecules, dioleoylphosphatidylserine (DOPS). From the density profile of water we observe that water next to the DOPS bilayer is more perturbed as compared to water near the DOPC bilayer. Using an energetic criterion for the determination of hydrogen bonding we find that water molecules create strong hydrogen bonds with the headgroups of the phospholipid molecules. Due to the presence of these bonds and also due to the confinement of water, the translational and orientational dynamics of water at the interface are slowed down. The degree of slowing down of the dynamics depends upon the location of water molecules near a lipid headgroup.     

Dynamical motions of lipids and a finite size effect in simulations of bilayers

Molecular dynamics (MD) simulations of dipalmitoylphosphatidylcholine bilayers composed of 72 and 288 lipids are used to examine system size dependence on dynamical properties associated with the particle mesh Ewald (PME) treatment of electrostatic interactions. The lateral diffusion constant Dl is 2.92 x 10(-7) and 0.95 x 10(-7) cm2/s for 72 and 288 lipids, respectively. This dramatic finite size effect originates from the correlation length of lipid diffusion, which extends to next-nearest neighbors in the 288 lipid system. Consequently, diffusional events in smaller systems can propagate across the boundaries of the periodic box. The internal dynamics of lipids calculated from the PME simulations are independent of the system size. Specifically, reorientational correlation functions for the slowly relaxing phosphorus-glycerol hydrogen, phosphorus-nitrogen vectors, and more rapidly relaxing CH vectors in the aliphatic chains are equivalent for the 72 and 288 lipid simulations. A third MD simulation of a bilayer with 72 lipids using spherical force-shift electrostatic cutoffs resulted in interdigitated chains, thereby rendering this cutoff method inappropriate.     

Poly(amidoamine) dendrimers on lipid bilayers I: Free energy and conformation of binding

Third-generation (G3) poly(amidoamine) (PAMAM) dendrimers are simulated approaching 1,2-dimyristoyl- sn-glycero-3-phosphocholine (DMPC) bilayers with fully atomistic molecular dynamics, which enables the calculation of a free energy profile along the approach coordinate. Three different dendrimer terminations are examined: protonated primary amine, uncharged acetamide, and deprotonated carboxylic acid. As the dendrimer and lipids become closer, their attractive force increases (up to 240 pN) and the dendrimer becomes deformed as it interacts with the lipids. The total energy release upon binding of a G3-NH3+, G3-Ac, or G3-COO- dendrimer to a DMPC bilayer is, respectively, 36, 26, or 47 kcal/mol or, equivalently, 5.2, 3.2, or 4.7x10(-3) kcal/g. These results are analyzed in terms of the dendrimers' size, shape, and atomic distributions as well as proximity of individual lipid molecules and particular lipid atoms to the dendrimer. For example, an area of 9.6, 8.2, or 7.9 nm2 is covered on the bilayer for the G3-NH3+, G3-Ac, or G3-COO- dendrimers, respectively, while interacting strongly with 18-13 individual lipid molecules.     

A multiscale coarse-graining method for biomolecular systems

A new approach is presented for obtaining coarse-grained (CG) force fields from fully atomistic molecular dynamics (MD) trajectories. The method is demonstrated by applying it to derive a CG model for the dimyristoylphosphatidylcholine (DMPC) lipid bilayer. The coarse-graining of the interparticle force field is accomplished by an application of a force-matching procedure to the force data obtained from an explicit atomistic MD simulation of the biomolecular system of interest. Hence, the method is termed a "multiscale" CG (MS-CG) approach in which explicit atomistic-level forces are propagated upward in scale to the coarse-grained level. The CG sites in the lipid bilayer application were associated with the centers-of-mass of atomic groups because of the simplicity in the evaluation of the forces acting on them from the atomistic data. The resulting CG lipid bilayer model is shown to accurately reproduce the structural properties of the phospholipid bilayer.     

Perturbation of phospholipid bilayer properties by ethanol at a high concentration

Atomistic molecular dynamics (MD) simulations have been carried out at 30 degrees C on a fully hydrated liquid crystalline lamellar phase of dimyrystoylphosphatidylcholine (DMPC) lipid bilayer with embedded ethanol molecules at 1:1 composition, as well as on the pure bilayer phase. The ethanol molecules are found to exhibit a preference to occupy regions near the upper part of the lipid acyl chains and the phosphocholine headgroups. The calculations revealed that the phosphocholine headgroup dipoles (P- --> N+) of the lipids prefer to orient more toward the aqueous layer in the presence of ethanol. It is noticed that the ethanol molecules modify the dynamic properties of both lipids as well as the water molecules in the hydration layer of the lipid headgroups. Both the in-plane "rattling" and out-of-plane "protrusion" motions of the lipids have been found to increase in the presence of ethanol. Most importantly, it is observed that the water molecules within the hydration layer of the lipid headgroups exhibit faster translational and rotational motions in the presence of ethanol. This arises due to faster dynamics of hydrogen bonds between lipid headgroups and water in the presence of ethanol.     

Modulating membrane properties: the effect of trehalose and cholesterol on a phospholipid bilayer

The protective properties of trehalose on cholesterol-containing lipid dipalmitoylphosphatidylcholine (DPPC) bilayers are studied through molecular simulations. The ability of the disaccharide to interact with the phospholipid headgroups and stabilize the membrane persists even at high cholesterol concentrations and restricts some of the changes to the structure that would otherwise be imposed by cholesterol molecules. Predictions of bilayer properties such as area per lipid, tail ordering, and chain conformation support the notion that the disaccharide decreases the main melting transition in these multicomponent model membranes, which correspond more closely to common biological systems than pure bilayers. Molecular simulations indicate that the membrane dynamics are slowed considerably by the presence of trehalose, indicating that high sugar concentrations would serve to avert possible phase separations that could arise in mixed phospholipid systems. Various time correlation functions suggest that the character of the modifications in lipid dynamics induced by trehalose and cholesterol is different in the hydrophilic and hydrophobic regions of the membrane.     

Computational analysis of local membrane properties

In the field of biomolecular simulations, dynamics of phospholipid membranes is of special interest. A number of proteins, including channels, transporters, receptors and short peptides are embedded in lipid bilayers and tightly interact with phospholipids. While the experimental measurements report on the spatial and/or temporal average membrane properties, simulation results are not restricted to the average properties. In the current study, we present a collection of methods for an efficient local membrane property calculation, comprising bilayer thickness, area per lipid, deuterium order parameters, Gaussian and mean curvature. The local membrane property calculation allows for a direct mapping of the membrane features, which subsequently can be used for further analysis and visualization of the processes of interest. The main features of the described methods are highlighted in a number of membrane systems, namely: a pure dimyristoyl-phosphatidyl-choline (DMPC) bilayer, a fusion peptide interacting with a membrane, voltage-dependent anion channel protein embedded in a DMPC bilayer, cholesterol enriched bilayer and a coarse grained simulation of a curved palmitoyl-oleoyl-phosphatidyl-choline lipid membrane. The local membrane property analysis proves to provide an intuitive and detailed view on the observables that are otherwise interpreted as averaged bilayer properties.     

Effect of ions on a dipalmitoyl phosphatidylcholine bilayer. a molecular dynamics simulation study

The effect of physiological concentrations of different chlorides on the structure of a dipalmitoyl phosphatidylcholine (DPPC) bilayer has been investigated through atomistic molecular dynamics simulations. These calculations provide support to the concept that Li+, Na+, Ca2+, Mg2+, Sr2+, Ba2+, and Ac3+, but not K+, bind to the lipid-head oxygens. Ion binding exhibits an influence on lipid order, area per lipid, orientation of the lipid head dipole, the charge distribution in the system, and therefore the electrostatic potential across the head-group region of the bilayer. These structural effects are sensitive to the specific characteristics of each cation, i.e., radius, charge, and coordination properties. These results provide evidence aimed at shedding some light into the apparent contradictions among different studies reported recently regarding the ordering effect of ions on zwitterionic phosphatidylcholine lipid bilayers.     

Molecular dynamics study on the stabilization of dehydrated lipid bilayers with glucose and trehalose

In this study, we use molecular dynamics simulations to investigate and compare the interactions of DPPC bilayers with and without saccharides (glucose or trehalose) under dehydrated conditions. Results from the simulations indicate that unilamellar bilayers lose their structural integrity under dehydrated conditions in the absence of saccharides; however, in the presence of either glucose or trehalose, the bilayers maintain their stability. Hydrogen bond analysis shows that the saccharide molecules displace a significant amount of water surrounding the lipid headgroups. At the same time, the additional hydrogen bonds formed between water and saccharide molecules help to maintain a hydration layer on the lipid bilayer interface. On the basis of the hydrogen bond distributions, trehalose forms more hydrogen bonds with the lipids than glucose, and it is less likely to interact with neighboring saccharide molecules. These results suggest that the interaction between the saccharide and lipid molecules through hydrogen bonds is an essential component of the mechanism for the stabilization of lipid bilayers.     

Molecular dynamics simulations of cardiolipin bilayers

Cardiolipin is a key lipid component in the inner mitochondrial membrane, where the lipid is involved in energy production, cristae structure, and mechanisms in the apoptotic pathway. In this article we used molecular dynamics computer simulations to investigate cardiolipin and its effect on the structure of lipid bilayers. Three cardiolipin/POPC bilayers with different lipid compositions were simulated: 100, 9.2, and 0% cardiolipin. We found strong association of sodium counterions to the carbonyl groups of both lipid types, leaving in the case of 9.2% cardiolipin virtually no ions in the aqueous compartment. Although binding occurred primarily at the carbonyl position, there was a preference to bind to the carbonyl groups of cardiolipin. Ion binding and the small headgroup of cardiolipin gave a strong ordering of the hydrocarbon chains. We found significant effects in the water dipole orientation and water dipole potential which can compensate for the electrostatic repulsion that otherwise should force charged lipids apart. Several parameters relevant for the molecular structure of cardiolipin were calculated and compared with results from analyses of coarse-grained simulations and available X-ray structural data.     

Comparative molecular dynamics study of ether- and ester-linked phospholipid bilayers

The lipid membranes found in archaea have high bilayer stability and low permeability. The molecular structure of their constituent lipids is characterized by ether-linked, branched hydrophobic chains, whereas the conventional lipids obtained from eukaryotic or eubacterial sources have ester linked straight chains. In order to elucidate the influence of the ether linkage, instead of an ester one, on the physical properties of the lipid bilayers, we have carried out comparative 10 ns molecular dynamics simulations of diphytanyl phosphatidylcholine (ether-DPhPC) and diphytanoyl phosphatidylcholine (ester-DPhPC) bilayers in water, respectively. We analyze bilayer structures, hydration of the lipids, membrane dipole potentials, and free energy profiles of water and oxygen across the bilayers. We observe that the membrane dipole potential for the ether-DPhPC bilayer, which arises mainly from the ether linkage, is about half of that of the ester-DPhPC. The calculated free energy barrier for a water molecule in the ether-DPhPC bilayer system is slightly higher than that in the ester-DPhPC counterpart, which is in accord with experimental data.