Molecular dynamics simulations of lipid bilayers: major artifacts due to truncating electrostatic interactions

SMARTINI3 parametrization of multi-scale membrane models via unsupervised learning methods

In this study, we utilize genetic algorithms to develop a realistic implicit solvent ultra-coarse-grained (ultra-CG) membrane model comprising only three interaction sites. The key philosophy of the ultra-CG membrane model SMARTINI3 is its compatibility with realistic membrane proteins, for example, modeled within the Martini coarse-grained (CG) model, as well as with the widely used GROMACS software for molecular simulations. Our objective is to parameterize this ultra-CG model to accurately reproduce the experimentally observed structural and thermodynamic properties of Phosphatidylcholine (PC) membranes in real units, including properties such as area per lipid, area compressibility, bending modulus, line tension, phase transition temperature, density profile, and radial distribution function. In our example, we specifically focus on the properties of a POPC membrane, although the developed membrane model could be perceived as a generic model of lipid membranes. To optimize the performance of the model (the fitness), we conduct a series of evolutionary runs with diverse random initial population sizes (ranging from 96 to 384). We demonstrate that the ultra-CG membrane model we developed exhibits authentic lipid membrane behaviors, including self-assembly into bilayers, vesicle formation, membrane fusion, and gel phase formation. Moreover, we demonstrate compatibility with the Martini coarse-grained model by successfully reproducing the behavior of a transmembrane domain embedded within a lipid bilayer. This facilitates the simulation of realistic membrane proteins within an ultra-CG bilayer membrane, enhancing the accuracy and applicability of our model in biophysical studies.

Molecular insights into the interactions of theaflavin and epicatechin with different lipid bilayer membranes by molecular dynamics simulation

At present, consumers increasingly favored the natural food preservatives with fewer side-effects on health. The green tea catechins and black tea theaflavins attracted considerable interest, and their antibacterial effects were extensively reported in the literature. Epicatechin (EC), a green tea catechin without a gallate moiety, showed no bactericidal activity, whereas the theaflavin (TF), also lacking a gallate moiety, exhibited potent bactericidal activity, and the antibacterial effects of green tea catechins and black tea theaflavins were closely correlated with their abilities to disrupt the bacterial cell membrane. In our present study, the mechanisms of membrane interaction modes and behaviors of TF and EC were explored by molecular dynamics simulations. It was demonstrated that TF exhibited markedly stronger affinity for the POPG bilayer compared to EC. Additionally, the hydrophobic interactions of tropolone/catechol rings with the acyl chain part could significantly contribute to the penetration of TF into the POPG bilayer. It was also found that the resorcinol/pyran rings were the key functional groups in TF for forming hydrogen bonds with the POPG bilayer. We believed that the findings from our current study could offer useful insights to better understand the stronger antibacterial effects of TF compared to EC.

Molecular dynamics simulations reliably identify vibrational modes in far-IR spectra of phospholipids

The properties of self-assembled phospholipid membranes are of essential importance in biochemistry and physical chemistry, providing a platform for many cellular life functions. Far-infrared (far-IR) vibrational spectroscopy, on the other hand, is a highly information-rich method to characterize intermolecular interactions and collective behaviour of lipids that can help explain, e.g., chain packing, thermodynamic phase behaviour, and sequestration. However, reliable interpretation of the far-IR spectra is still lacking. Here we present a molecular dynamics (MD) based approach to simulate vibrational modes of individual lipids and in an ensemble. The results are a good match to synchrotron far-IR measurements and enable identification of the molecular motions corresponding to each vibrational mode, thus allowing the correct interpretation of membrane spectra with high accuracy and resolving the longstanding ambiguities in the literature in this regard. Our results demonstrate the feasibility of using MD simulations for interpreting far-IR spectra broadly, opening new avenues for practical use of this powerful method.

Analytic Gradients for the Electrostatic Embedding QM/MM Model in Periodic Boundary Conditions Using Particle-Mesh Ewald Sums and Electrostatic Potential Fitted Charge Operators

Long-range electrostatic effects are fundamental for describing chemical reactivity in the condensed phase. Here, we present the methodology of an efficient quantum mechanical/molecular mechanical (QM/MM) model in periodic boundary conditions (PBC) compatible with QM/MM boundaries at chemical bonds. The method combines electrostatic potential fitted charge operators and electrostatic potentials derived from the smooth particle-mesh Ewald (PME) sum approach. The total energy and its analytic first derivatives with respect to QM, MM, and lattice vectors allow QM/MM molecular dynamics (MD) in the most common thermodynamic ensembles. We demonstrate the robustness of the method by performing a QM/MM MD equilibration of methanol in water. We simulate the cis/trans isomerization free-energy profiles in water of proline amino acid and a proline-containing oligopeptide, showing a correct description of the reaction barrier. Our PBC-compatible QM/MM model can efficiently be used to study the chemical reactivity in the condensed phase and enzymatic catalysis.

Unsolved problem of long-range interactions: dipolar spin-ice study

Long-range interactions derive various strange phenomena. As illustrated by cutoff simulations of water, increasing cutoff length does not improve the simulation result necessarily; on the contrary, it makes the result worse. In the extreme situation, the structure of water transforms into a layer structure. In this study, to explore the underlying mechanism of this phenomenon, we performed Monte Carlo simulations on dipolar spins arranged on a pyrochlore spin-ice lattice. Like the water case, the present dipolar spin system also showed cutoff-induced dipole ordering and layer formation. The width of the layers depended on the cutoff length; and longer cutoff length led to a broader layer. These features are certainly consistent with the previous water case. This indicates that layer formation is the general behavior of dipolar systems whose interactions are truncated within a finite distance. The result is important for future exploration of the relationship between long-range interactions and resulting structures. In addition, it emphasizes the necessity of rigorous treatment of long-range interactions because increasing the cutoff length prevents convergence and provides an entirely different result from the rigorous Ewald calculation.

Molecular insights into the structure destabilization effects of ECG and EC on the Aβ protofilament: An all-atom molecular dynamics simulation study

The formation of Aβ into amyloid fibrils was closely connected to AD, therefore, the Aβ aggregates were the primary therapeutic targets against AD. Previous studies demonstrated that epicatechin-3-gallate (ECG), which possessed a gallate moiety, exhibited a greater ability to disrupt the preformed Aβ amyloid fibrils than epicatechin (EC), indicating that the gallate moiety was crucial. In the present study, the molecular mechanisms were investigated. Our results demonstrated that ECG had more potent disruptive impacts on the β-sheet structure and K28-A42 salt bridges than EC. We found that ECG significantly interfered the interactions between Peptide-4 and Peptide-5. However, EC could not. The disruption of K28-A42 salt bridges by ECG was mainly due to the interactions between ECG and the hydrophobic residues located at C-terminus. Interestingly, EC disrupted the K28-A42 salt bridges by the interactions with C-terminal hydrophobic residues and the cation-π interactions with K28. Moreover, our results indicated that hydrophobic interactions, H-bonds, π-π interactions and cation-π interactions between ECG and the bend of L-shaped region caused the disaggregation of interactions between Peptide-4 and Peptide-5. Significantly, gallate moiety in ECG had contributed tremendously to the disaggregation. We believed that our findings could be useful for designing prospective drug candidates targeting AD.

Probing the Hydrophobic Region of a Lipid Bilayer at Specific Depths Using Vibrational Spectroscopy

A novel spectroscopic approach for studying the flexibility and mobility in the hydrophobic interior of lipid bilayers at specific depths is proposed. A set of test compounds featuring an azido moiety and a cyano or carboxylic acid moiety, connected by an alkyl chain of different lengths, was synthesized. FTIR data and molecular dynamics calculations indicated that the test compounds in a bilayer are oriented so that the cyano or carboxylic acid moiety is located in the lipid head-group region, while the azido group stays inside the bilayer at the depth determined by its alkyl chain length. We found that the asymmetric stretching mode of the azido group (νN3) can serve as a reporter of the membrane interior dynamics. FTIR and two-dimensional infrared (2DIR) studies were performed at different temperatures, ranging from 22 to 45 °C, covering the Lβ-Lα phase transition temperature of dipalmitoylphosphatidylcholine (∼41 °C). The width of the νN3 peak was found to be very sensitive to the phase transition and to the temperature in general. We introduced an order parameter, SN3, which characterizes restrictions to motion inside the bilayer. 2DIR spectra of νN3 showed different extents of inhomogeneity at different depths in the bilayer, with the smallest inhomogeneity in the middle of the leaflet. The spectral diffusion dynamics of the N3 peak was found to be dependent on the depth of the N3 group location in the bilayer. The obtained results enhance our understanding of the bilayer dynamics and can be extended to investigate membranes with more complex compositions.

SPEADI: Accelerated Analysis of IDP-Ion Interactions from MD-Trajectories

The disordered nature of Intrinsically Disordered Proteins (IDPs) makes their structural ensembles particularly susceptible to changes in chemical environmental conditions, often leading to an alteration of their normal functions. A Radial Distribution Function (RDF) is considered a standard method for characterizing the chemical environment surrounding particles during atomistic simulations, commonly averaged over an entire or part of a trajectory. Given their high structural variability, such averaged information might not be reliable for IDPs. We introduce the Time-Resolved Radial Distribution Function (TRRDF), implemented in our open-source Python package SPEADI, which is able to characterize dynamic environments around IDPs. We use SPEADI to characterize the dynamic distribution of ions around the IDPs Alpha-Synuclein (AS) and Humanin (HN) from Molecular Dynamics (MD) simulations, and some of their selected mutants, showing that local ion-residue interactions play an important role in the structures and behaviors of IDPs.

Conformational dynamics of A30G α-synuclein that causes familial Parkinson disease

The first gene shown to be responsible for autosomal-dominant Parkinson's disease (PD) is the SNCA gene, which encodes for alpha synuclein (α-Syn). Recently, a novel heterozygous A30G mutation of the SNCA gene associated with familial PD has been reported. However, little research has been done on how the A30G mutation affects the structure of α-Syn. So, using atomistic molecular dynamics (MD) simulation, we demonstrate here the key structural characteristics of A30G α-Syn in the free monomer form and membrane associated state. From the MD trajectory analysis, the structure of A30G α-Syn was noticed to exhibit rapid conformational change, increase in backbone flexibility near the site of mutation and decrease in α-helical propensity. The typical torsion angles in residues (Val26 and Glu28) near the mutation site were observed to deviate significantly in A30G α-Syn. In the case of membrane bound A30G α-Syn, the regions that were submerged in the lipid bilayer (N-helix (3-37) and turn region (38-44)) found to contain higher helical content than the elevated region above the lipid surface. The bending angle in the helix-N and helix-C regions were noticed to be relatively higher in the free form of A30G α-Syn (38.50) than in the membrane bound form (370). The A30G mutation in α-Syn was predicted to have an impact on the stability and function of the protein based on ΔΔG values obtained from the online servers. Our results demonstrate that the A30G mutation in α-Syn altered the protein's α-helical structure and slightly altered the membrane binding.Communicated by Ramaswamy H. Sarma.

APL@voro-interactive visualization and analysis of cell membrane simulations

SUMMARY

Molecular dynamics (MD) simulations of cell membranes allow for a better understanding of complex processes such as changing membrane dynamics, lipid rafts and the incorporation/passing of macromolecules into/through membranes. To explore and understand cell membrane compositions, dynamics and processes, visual analytics can help to interpret MD simulation data. APL@Voro is a software for the interactive visualization and analysis of cell membrane simulations. Here, we present the new APL@Voro, which has been continuously developed since its initial release in 2013. We discuss newly implemented algorithms, methodologies and features, such as the interactive comparison of related simulations and methods to assign lipids to either the upper or lower leaflet.

AVAILABILITY AND IMPLEMENTATION

The current open-source version of APL@Voro can be downloaded from http://aplvoro.com.

Effect of Local Stress on Accurate Modeling of Bacterial Outer Membranes Using All-Atom Molecular Dynamics

Biological membranes are fundamental components of living organisms that play an undeniable role in their survival. Molecular dynamics (MD) serves as an essential computational tool for studying biomembranes on molecular and atomistic scales. The status quo of MD simulations of biomembranes studies a nanometer-sized membrane patch periodically extended under periodic boundary conditions (PBCs). In nature, membranes are usually composed of different lipids in their two layers (referred to as leaflets). This compositional asymmetry imposes a fixed ratio of lipid numbers between the two leaflets in a periodically constrained membrane, which needs to be set appropriately. The widely adopted methods of defining a leaflet lipid ratio suffer from the lack of control over the mechanical tension of each leaflet, which could significantly influence research findings. In this study, we investigate the role of membrane-building protocol and the resulting initial stress state on the interaction between small molecules and asymmetric membranes. We model the outer membrane of Pseudomonas aeruginosa bacteria using two different building protocols and probe their interactions with the Pseudomonas quinolone signal (PQS). Our results show that differential stress could shift the position of free energy minimum for the PQS molecule between the two leaflets of the asymmetric membrane. This work provides critical insights into the relationship between the initial per-leaflet tension and the spontaneous intercalation of PQS.

Lamellarity-Driven Differences in Surface Structural Features of DPPS Lipids: Spectroscopic, Calorimetric and Computational Study

Although single-lipid bilayers are usually considered models of eukaryotic plasma membranes, their research drops drastically when it comes to exclusively anionic lipid membranes. Being a major anionic phospholipid in the inner leaflet of eukaryote membranes, phosphatidylserine-constituted lipid membranes were occasionally explored in the form of multilamellar liposomes (MLV), but their inherent instability caused a serious lack of efforts undertaken on large unilamellar liposomes (LUVs) as more realistic model membrane systems. In order to compensate the existing shortcomings, we performed a comprehensive calorimetric, spectroscopic and MD simulation study of time-varying structural features of LUV made from 1,2-dipalmitoyl-sn-glycero-3-phospho-L-serine (DPPS), whereas the corresponding MLV were examined as a reference. A substantial uncertainty of UV/Vis data of LUV from which only T_m was unambiguously determined (53.9 ± 0.8 °C), along with rather high uncertainty on the high-temperature range of DPPS melting profile obtained from DSC (≈50-59 °C), presumably reflect distinguished surface structural features in LUV. The FTIR signatures of glycerol moiety and those originated from carboxyl group serve as a strong support that in LUV, unlike in MLV, highly curved surfaces occur continuously, whereas the details on the attenuation of surface features in MLV were unraveled by molecular dynamics.

Bergamottin: location, aggregation and interaction with the plasma membrane

Bioactive furanocoumarins, a group of natural secondary metabolites common in higher plants, are recognized for their benefits to human health and have been shown to have numerous biological properties. However, the knowledge of its biomolecular mechanism is not known. One of the main furanocoumarins is bergamottin (BGM), which is characterized by a planar three-ringed structure and a hydrocarbon chain, which give BGM its high lipid/water partition coefficient. Because of that, and although the biological mechanism of BGM is not known, BGM bioactive properties could be ascribed to its potential to interact with the biological membrane, modulating its structure, changing its dynamics and at the same time that it might interact with lipids. For our goal, we have applied molecular dynamics to determine the position of BGM in a complex membrane and discern the possibility of certain interactions with membrane lipids. Our findings establish that BGM tends to locate in the middle of the hydrocarbon layer of the membrane, inserts in between the hydrocarbon chains of the phospholipids in an oblique position with respect to the membrane plane, increasing the fluidity of the membrane. Significantly, BGM tends to be surrounded by POPC molecules but exclude the molecule of CHOL. Outstandingly, BGM molecules associate spontaneously creating aggregates, which does not preclude them from interacting with and inserting into the membrane. The bioactive properties of BGM could be ascribed to its membranotropic effects and support the improvement of these molecules as therapeutic molecules, giving place to new opportunities for potential medical improvements.Communicated by Ramaswamy H. Sarma.

Fuzzy spherical truncation-based multi-linear protein descriptors: From their definition to application in structural-related predictions

This study introduces a set of fuzzy spherically truncated three-dimensional (3D) multi-linear descriptors for proteins. These indices codify geometric structural information from kth spherically truncated spatial-(dis)similarity two-tuple and three-tuple tensors. The coefficients of these truncated tensors are calculated by applying a smoothing value to the 3D structural encoding based on the relationships between two and three amino acids of a protein embedded into a sphere. At considering, the geometrical center of the protein matches with center of the sphere, the distance between each amino acid involved in any specific interaction and the geometrical center of the protein can be computed. Then, the fuzzy membership degree of each amino acid from an spherical region of interest is computed by fuzzy membership functions (FMFs). The truncation value is finally a combination of the membership degrees from interacting amino acids, by applying the arithmetic mean as fusion rule. Several fuzzy membership functions with diverse biases on the calculation of amino acids memberships (e.g., Z-shaped (close to the center), PI-shaped (middle region), and A-Gaussian (far from the center)) were considered as well as traditional truncation functions (e.g., Switching). Such truncation functions were comparatively evaluated by exploring: 1) the frequency of membership degrees, 2) the variability and orthogonality analyses among them based on the Shannon Entropy’s and Principal Component’s methods, respectively, and 3) the prediction performance of alignment-free prediction of protein folding rates and structural classes. These analyses unraveled the singularity of the proposed fuzzy spherically truncated MDs with respect to the classical (non-truncated) ones and respect to the MDs truncated with traditional functions. They also showed an improved prediction power by attaining an external correlation coefficient of 95.82% in the folding rate modelling and an accuracy of 100% in distinguishing structural protein classes. These outcomes are better than the ones attained by existing approaches, justifying the theoretical contribution of this report. Thus, the fuzzy spherically truncated-based protein descriptors from MuLiMs-MCoMPAs (http://tomocomd.com/mulims-mcompas) are promising alignment-free predictors for modeling protein functions and properties.

Reversible Unwrapping Algorithm for Constant-Pressure Molecular Dynamics Simulations

Molecular simulation technologies have afforded researchers a unique look into the nanoscale interactions driving physical processes. However, a limitation for molecular dynamics (MD) simulations is that they must be performed on finite-sized systems in order to map onto computational resources. To minimize artifacts arising from finite-sized simulation systems, it is common practice for MD simulations to be performed with periodic boundary conditions (PBCs). However, in order to calculate specific physical properties, such as mean square displacements to calculate diffusion coefficients, continuous particle trajectories where the atomic movements are continuous and do not jump between cell faces are required. In these cases, modifying atomic coordinates through unwrapping schemes is an essential post-processing tool to remove these jumps. Here, two established trajectory unwrapping schemes are applied to 1 μs wrapped trajectories for a small water box and lysozyme in water. The existing schemes can result in spurious diffusion coefficients, long bonds within unwrapped molecules, and inconsistent atomic coordinates when coordinates are rewrapped after unwrapping. We determine that prior unwrapping schemes do not account for changing periodic box dimensions and introduce an additional correction term to the existing displacement unwrapping scheme to correct for these artifacts. We also demonstrate that the resulting algorithm is a hybrid between the existing heuristic and displacement unwrapping schemes. After treatment using this new unwrapping scheme, molecular geometries are correct even after long simulations. In anticipation for longer MD trajectories, we develop implementations for this new scheme in multiple PBC handling tools.

Review of Electrostatic Force Calculation Methods and Their Acceleration in Molecular Dynamics Packages Using Graphics Processors

Molecular dynamics (MD) simulations probe the conformational repertoire of macromolecular systems using Newtonian dynamic equations. The time scales of MD simulations allow the exploration of biologically relevant phenomena and can elucidate spatial and temporal properties of the building blocks of life, such as deoxyribonucleic acid (DNA) and protein, across microsecond (μs) time scales using femtosecond (fs) time steps. A principal bottleneck toward extending MD calculations to larger time scales is the long-range electrostatic force measuring component of the naive nonbonded force computation algorithm, which scales with a complexity of (N, number of atoms). In this review, we present various methods to determine electrostatic interactions in often-used open-source MD packages as well as the implementation details that facilitate acceleration of the electrostatic interaction calculation.

Fast Prediction of Binding Affinities of SARS-CoV-2 Spike Protein and Its Mutants with Antibodies through Intermolecular Interaction Modeling-Based Machine Learning

Since the introduction of the novel SARS-CoV-2 virus (COVID-19) in late 2019, various new variants have appeared with mutations that confer resistance to the vaccines and monoclonal antibodies that were developed in response to the wild-type virus. As we continue through the pandemic, an accurate and efficient methodology is needed to help predict the effects certain mutations will have on both our currently produced therapeutics and those that are in development. Using published cryo-electron microscopy and X-ray crystallography structures of the spike receptor binding domain region with currently known antibodies, in the present study, we created and cross-validated an intermolecular interaction modeling-based multi-layer perceptron machine learning approach that can accurately predict the mutation-caused shifts in the binding affinity between the spike protein (wild-type or mutant) and various antibodies. This validated artificial intelligence (AI) model was used to predict the binding affinity (Kd) of reported SARS-CoV-2 antibodies with various variants of concern, including the most recently identified "Deltamicron" (or "Deltacron") variant. This AI model may be employed in the future to predict the Kd of developed novel antibody therapeutics to overcome the challenging antibody resistance issue and develop structural bases for the effects of both current and new mutants of the spike protein. In addition, the similar AI strategy and approach based on modeling of the intermolecular interactions may be useful in development of machine learning models predicting binding affinities for other protein-protein binding systems, including other antibodies binding with their antigens.

The role of structural order in heterogeneous ice nucleation

The freezing of water into ice is a key process that is still not fully understood. It generally requires an impurity of some description to initiate the heterogeneous nucleation of the ice crystals. The molecular structure, as well as the extent of structural order within the impurity in question, both play an essential role in determining its effectiveness. However, disentangling these two contributions is a challenge for both experiments and simulations. In this work, we have systematically investigated the ice-nucleating ability of the very same compound, cholesterol, from the crystalline (and thus ordered) form to disordered self-assembled monolayers. Leveraging a combination of experiments and simulations, we identify a “sweet spot” in terms of the surface coverage of the monolayers, whereby cholesterol maximises its ability to nucleate ice (which remains inferior to that of crystalline cholesterol) by enhancing the structural order of the interfacial water molecules. These findings have practical implications for the rational design of synthetic ice-nucleating agents.

The freezing of water into ice is still not fully understood. Here, we investigate the role of structural disorder within the biologically relevant impurities that facilitate this fundamental phase transition.

Envelope E protein of dengue virus and phospholipid binding to the late endosomal membrane

Flaviviruses include many significant human pathogens, comprising dengue, West Nile, Yellow fever, Japanese encephalitis, Zika and tick-borne encephalitis viruses and many others, affecting millions of people in the world. These viruses have produced important epidemics in the past, they continue to do it and they will undoubtedly continue to do so in the future. Flaviviruses enter into the cells via receptor-mediated endocytosis by fusing its membrane with the endosomal membrane in a pH-dependent manner with the help of the envelope E protein, a prototypical class II membrane fusion protein. The envelope E protein has a conserved fusion peptide at its distal end, which is responsible in the first instance of inserting the protein into the host membrane. Since the participation of other segments of the E protein in the fusion process should not be ruled out, we have used atomistic molecular dynamics to study the binding of the distal end of domain II of the envelope E protein from Dengue virus (DENV) with a complex membrane similar to the late-endosome one. Our work shows that not only the fusion peptide participates directly in the fusion, but also two other sequences of the protein, next to the fusion peptide it in the three-dimensional structure, are jointly wrapped in the fusion process. Overall, these three sequences represent a new target that would make it possible to obtain effective antivirals against DENV in particular and Flaviviruses in general.

Generalized Methodology for the Quick Prediction of Variant SARS-CoV-2 Spike Protein Binding Affinities with Human Angiotensin-Converting Enzyme II

Variants of the SARS-CoV-2 virus continue to remain a threat 2 years from the beginning of the pandemic. As more variants arise, and the B.1.1.529 (Omicron) variant threatens to create another wave of infections, a method is needed to predict the binding affinity of the spike protein quickly and accurately with human angiotensin-converting enzyme II (ACE2). We present an accurate and convenient energy minimization/molecular mechanics Poisson-Boltzmann surface area methodology previously used with engineered ACE2 therapeutics to predict the binding affinity of the Omicron variant. Without any additional data from the variants discovered after the publication of our first model, the methodology can accurately predict the binding of the spike/ACE2 variant complexes. From this methodology, we predicted that the Omicron variant spike has a Kd of ∼22.69 nM (which is very close to the experimental Kd of 20.63 nM published during the review process of the current report) and that spike protein of the new "Stealth" Omicron variant (BA.2) will display a Kd of ∼12.9 nM with the wild-type ACE2 protein. This methodology can be used with as-yet discovered variants, allowing for quick determinations regarding the variant's infectivity versus either the wild-type virus or its variants.

Molecular Modeling is an Enabling Approach to Complement and Enhance Channelopathy Research

Hundreds of human membrane proteins form channels that transport necessary ions and compounds, including drugs and metabolites, yet details of their normal function or how function is altered by genetic variants to cause diseases are often unknown. Without this knowledge, researchers are less equipped to develop approaches to diagnose and treat channelopathies. High-resolution computational approaches such as molecular modeling enable researchers to investigate channelopathy protein function, facilitate detailed hypothesis generation, and produce data that is difficult to gather experimentally. Molecular modeling can be tailored to each physiologic context that a protein may act within, some of which may currently be difficult or impossible to assay experimentally. Because many genomic variants are observed in channelopathy proteins from high-throughput sequencing studies, methods with mechanistic value are needed to interpret their effects. The eminent field of structural bioinformatics integrates techniques from multiple disciplines including molecular modeling, computational chemistry, biophysics, and biochemistry, to develop mechanistic hypotheses and enhance the information available for understanding function. Molecular modeling and simulation access 3D and time-dependent information, not currently predictable from sequence. Thus, molecular modeling is valuable for increasing the resolution with which the natural function of protein channels can be investigated, and for interpreting how genomic variants alter them to produce physiologic changes that manifest as channelopathies. © 2022 American Physiological Society. Compr Physiol 12:3141-3166, 2022.

A deep potential model with long-range electrostatic interactions

Machine learning models for the potential energy of multi-atomic systems, such as the deep potential (DP) model, make molecular simulations with the accuracy of quantum mechanical density functional theory possible at a cost only moderately higher than that of empirical force fields. However, the majority of these models lack explicit long-range interactions and fail to describe properties that derive from the Coulombic tail of the forces. To overcome this limitation, we extend the DP model by approximating the long-range electrostatic interaction between ions (nuclei + core electrons) and valence electrons with that of distributions of spherical Gaussian charges located at ionic and electronic sites. The latter are rigorously defined in terms of the centers of the maximally localized Wannier distributions, whose dependence on the local atomic environment is modeled accurately by a deep neural network. In the DP long-range (DPLR) model, the electrostatic energy of the Gaussian charge system is added to short-range interactions that are represented as in the standard DP model. The resulting potential energy surface is smooth and possesses analytical forces and virial. Missing effects in the standard DP scheme are recovered, improving on accuracy and predictive power. By including long-range electrostatics, DPLR correctly extrapolates to large systems the potential energy surface learned from quantum mechanical calculations on smaller systems. We illustrate the approach with three examples: the potential energy profile of the water dimer, the free energy of interaction of a water molecule with a liquid water slab, and the phonon dispersion curves of the NaCl crystal.

Fast Prediction of Binding Affinities of the SARS-CoV-2 Spike Protein Mutant N501Y (UK Variant) with ACE2 and Miniprotein Drug Candidates

A recently identified variant of SARS-CoV-2 virus, known as the United Kingdom (UK) variant (lineage B.1.1.7), has an N501Y mutation on its spike protein. SARS-CoV-2 spike protein binds with angiotensin-converting enzyme 2 (ACE2), a key protein for the viral entry into the host cells. Here, we report an efficient computational approach, including the simple energy minimizations and binding free energy calculations, starting from an experimental structure of the binding complex along with experimental calibration of the calculated binding free energies, to rapidly and reliably predict the binding affinities of the N501Y mutant with human ACE2 (hACE2) and recently reported miniprotein and hACE2 decoy (CTC-445.2) drug candidates. It has been demonstrated that the N501Y mutation markedly increases the ACE2-spike protein binding affinity (Kd) from 22 to 0.44 nM, which could partially explain why the UK variant is more infectious. The miniproteins are predicted to have ∼10,000- to 100,000-fold diminished binding affinities with the N501Y mutant, creating a need for design of novel therapeutic candidates to overcome the N501Y mutation-induced drug resistance. The N501Y mutation is also predicted to decrease the binding affinity of a hACE2 decoy (CTC-445.2) binding with the spike protein by ∼200-fold. This convenient computational approach along with experimental calibration may be similarly used in the future to predict the binding affinities of potential new variants of the spike protein.

Influence of the gallate moiety on the interactions between green tea polyphenols and lipid membranes elucidated by molecular dynamics simulations

Previous studies suggested that naturally occurring EGCG primarily acted on the bacterial cell membrane then damaged the membrane and the gallate moiety in EGCG was very important to its anti-bacterial activity. However, the detailed mechanisms were still poorly understood. In this paper, EGCG and EGC were selected to study the great contribution of gallate moiety on the anti-bacterial activities of polyphenols. The results indicated that EGCG could penetrate deeper into the POPG lipid bilayer and possess more potent structure-perturbing potency on the POPG lipid bilayer than EGC. We also found that EGCG had the ability to form hydrogen bonds with the deeper inside oxygen atoms in the POPG lipid bilayer and the gallate moiety was the key functional group for EGCG forming hydrogen bonds with the POPG lipid bilayer. Moreover, results from the binding free energy analysis demonstrated that the gallate moiety made great contribution to the high affinity between EGCG and the POPG lipid bilayer. We believed that these findings could yield useful insights into the influence mechanisms of gallate moiety on the anti-bacterial activities of polyphenols.

Not only in silico drug discovery: Molecular modeling towards in silico drug delivery formulations

The use of methods at molecular scale for the discovery of new potential active ligands, as well as previously unknown binding sites for target proteins, is now an established reality. Literature offers many successful stories of active compounds developed starting from insights obtained in silico and approved by Food and Drug Administration (FDA). One of the most famous examples is raltegravir, a HIV integrase inhibitor, which was developed after the discovery of a previously unknown transient binding area thanks to molecular dynamics simulations. Molecular simulations have the potential to also improve the design and engineering of drug delivery devices, which are still largely based on fundamental conservation equations. Although they can highlight the dominant release mechanism and quantitatively link the release rate to design parameters (size, drug loading, et cetera), their spatial resolution does not allow to fully capture how phenomena at molecular scale influence system behavior. In this scenario, the "computational microscope" offered by simulations at atomic scale can shed light on the impact of molecular interactions on crucial parameters such as release rate and the response of the drug delivery device to external stimuli, providing insights that are difficult or impossible to obtain experimentally. Moreover, the new paradigm brought by nanomedicine further underlined the importance of such computational microscope to study the interactions between nanoparticles and biological components with an unprecedented level of detail. Such knowledge is a fundamental pillar to perform device engineering and to achieve efficient and safe formulations. After a brief theoretical background, this review aims at discussing the potential of molecular simulations for the rational design of drug delivery systems.

Dielectric continuum model examination of real-space electrostatic treatments

Electrostatic interaction is long ranged; thus, the accurate calculation is not an easy task in molecular dynamics or Monte Carlo simulations. Though the rigorous Ewald method based on the reciprocal space has been established, real-space treatments have recently become an attractive alternative because of the efficient calculation. However, the construction is not yet completed and is now a challenging subject. In an earlier theoretical study, Neumann and Steinhauser employed the Onsager dielectric continuum model to explain how simple real-space cutoff produces artificial dipolar orientation. In the present study, we employ this continuum model to explore the fundamental properties of the recently developed real-space treatments of three shifting schemes. The result of the distance-dependent Kirkwood function G_K(R) showed that the simple bare cutoff produces a well-known hole-shaped artifact, whereas the shift treatments do not. Two-dimensional mapping of electric field well explained how these shift treatments remove the hole-shaped artifact. Still, the shift treatments are not sufficient because they do not produce a flat G_K(R) profile unlike ideal no-cutoff treatment. To test the continuum model results, we also performed Monte Carlo simulations of dipolar particles. The results found that the continuum model could predict the qualitative tendency as to whether each electrostatic treatment produces the hole-shaped artifact of G_K(R) or not. We expect that the present study using the continuum model offers a stringent criterion to judge whether the primitive electrostatic behavior is correctly described or not, which will be useful for future construction of electrostatic treatments.

Molecular insights into the inhibitory mechanisms of gallate moiety on the Aβ1-40 amyloid aggregation: A molecular dynamics simulation study

Alzheimer's disease is the most common form of neurodegenerative disease and the formation of Aβ amyloid aggregates has been widely demonstrated to be the principal cause of Alzheimer's disease. Our previous study and other studies suggested that the gallate moiety played an obligatory role in the inhibition process of naturally occurring polyphenols on Aβ amyloid fibrils formation. However, the detailed mechanisms were still unknown. Thus, in the present study, the gallic acid (GA) was specially selected and the molecular recognition mechanisms between GA molecules and Aβ_1-40 monomer were examined and analyzed by molecular dynamics simulation. The in silico experiments revealed that GA significantly prevented the conformational changes of Aβ_1-40 monomer with no β-sheet structure during the whole 100 ns. By analyzing the binding sites of GA molecules to Aβ_1-40 monomer, we found that both hydrophilic and hydrophobic amino acid residues were participated in the binding of GA molecules to Aβ_1-40 monomer. Moreover, results from the binding free energy analysis further demonstrated that the strength of polar interactions was significantly stronger than that of nonpolar interactions. We believed that our results could help to elucidate the underlying mechanisms of gallate moiety on the anti-amyloidogenic effects of polyphenols at the atomic level.

Identification of the phospholipid binding regions of the envelope E protein of flaviviruses by molecular dynamics

The Flavivirus genus comprise several important human pathogens, including dengue, West Nile, Yellow fever, Japanese encephalitis, Zika, and tick-borne encephalitis viruses. These enveloped viruses affect more than 2 billion people in the world, mainly in less developed countries. Although some vaccines exist for some flaviviruses, these vaccines are not universally available due to many factors and since their infections are a world-wide public health issue, the development of antiviral molecules is fundamental. Flavivirus membranes, through the help of the envelope E glycoprotein, fuse with endosomal compartments in a pH-dependent way to release their genome into the cytoplasm and require specific lipids, such as bis(monoacylglycero)phosphate (BMP), for efficient fusion. The fundamental role the envelope E protein has on viral entry and membrane fusion suggest that it is an essential antiviral target. In this work, we have used atomistic molecular dynamics simulations to study the binding of the head-group of BMP to the tip of the envelope E proteins of ZIKV, DENV, TBEV and JEV viruses whose three-dimensional structures are known. Our results indicate that, apart from the fusion loop, there are different amino acid residues in different regions of the envelope E proteins of flaviviruses capable of binding the head-group of BMP. These regions should work together to accomplish the binding and fusion of the envelope and endosomal membranes and represent a new target to develop and design potent and effective antiviral agents capable of blocking flavivirus-endosome membrane fusion. [Formula: see text].

Smoothed Spherical Truncation based on Fuzzy Membership Functions: Application to the Molecular Encoding

A novel spherical truncation method, based on fuzzy membership functions, is introduced to truncate interatomic (or interaminoacid) relations according to smoothing values computed from fuzzy membership degrees. In this method, the molecules are circumscribed into a sphere, so that the geometric centers of the molecules are the centers of the spheres. The fuzzy membership degree of each atom (or aminoacid) is computed from its distance with respect to the geometric center of the molecule, by using a fuzzy membership function. So, the smoothing value to be applied in the truncation of a relation (or interaction) is computed by averaging the fuzzy membership degrees of the atoms (or aminoacids) involved in the relation. This truncation method is rather different from the existing ones, at considering the geometric center for the whole molecule and not only for atom-groups, as well as for using fuzzy membership functions to compute the smoothing values. A variability study on a set comprised of 20,469 compounds (15,050 drug-like compounds, 2994 drugs approved, 880 natural products from African sources, and 1545 plant-derived natural compounds exhibiting anti-cancerous activity) demonstrated that the truncation method proposed allows to determine molecular encodings with better ability for discriminating among structurally different molecules than the encodings obtained without applying truncation or applying non-fuzzy truncation functions. Moreover, a principal component analysis revealed that orthogonal chemical information of the molecules is encoded by using the method proposed. Lastly, a modeling study proved that the truncation method improves the modeling ability of existing geometric molecular descriptors, at allowing to develop more robust models than the ones built only using non-truncated descriptors. In this sense, a comparison and statistical assessment were performed on eight chemical datasets. As a result, the models based on the truncated molecular encodings yielded statistically better results than 12 procedures considered from the literature. It can thus be stated that the proposed truncation method is a relevant strategy for obtaining better molecular encodings, which will be ultimately useful in enhancing the modeling ability of existing encodings both on small-to-medium size molecules and biomacromolecules. © 2019 Wiley Periodicals, Inc.

Hybrid Atomistic and Coarse-Grained Model for Surfactants in Apolar Solvents

Here, we develop and verify the performance of a hybrid molecular modeling approach that combines coarse-grained apolar solvent and atomistic solute or polar solvent description, for example, for description of reverse micellar systems. The coarse-grained solvent model is directly applicable to organic solvents encompassing alkane, alkene, and fatty acid ester functional groups and connects directly to both standard united-atom GROMOS 53A6 and all-atom CHARMM27 force fields, as well as the atomistic detail water models compatible with these force fields. The different levels of description are coupled via explicit, unscaled electrostatics, and scaled mixing rules for dispersive interactions. The hybrid model is in near-quantitative agreement with fully atomistic simulations when combined with the CHARMM27 model but underestimates modestly surfactant aggregation when using GROMOS 53A6 united-atom description. The use of truncated electrostatics affords up to a 9-fold increase in computational speed without significant loss of accuracy. However, long-range electrostatic calculations and load imbalance at high core counts can significantly degrade the performance. We demonstrate the usability of the hybrid model by assessing the reverse micelle formation of a homologous series of nonionic glycerolipids via large-scale self-assembly simulations. The presented model is demonstrated here for accurate description of surfactant systems in apolar solvents, with and without also polar solvent (water) in the system. The formulation can be expected to describe well also other solute species or interfaces with an apolar solvent in an apolar environment.

Amphiphilic copolymers change the nature of the ordered-to-disordered phase transition of lipid membranes from discontinuous to continuous

The phase behaviour as a function of temperature is explored for pure phospholipid (DPPC) and hybrid lipid-polymer (DPPC/Pluronic L64) bilayers with the aid of atomistic MD simulations. The range of the fixed-temperature simulations includes temperatures below and above the known melting temperature (Tm) for DPPC membranes. For the pure lipid bilayer, the main phase transition is discontinuous, as verified by the abrupt changes observed in the membrane structure, elasticity and the lipid diffusivity near the critical temperature Tm, which lies in the region 298.15-303.15 K. A pre-transition step is detected at 298.15 K which has been identified as the ripple phase (Pβ'), where ordered and disordered lipids coexist, causing thickness fluctuations. In the ordered gel phase, the positional ordering as assessed by the lipid radial distribution functions is long-range and some degree of hexagonal packing is measured. The hybrid bilayers on the other hand, transform from a more ordered to a disordered phase in a continuous manner, without finite jumps in their properties. No signs of the ripple phase are identified and the ordered phase exhibits very limited hexagonal packing and some positional ordering that decays fast. The effect of the inserted polymers in the two phases is reversed; at low temperatures, they render the membrane thinner, less cohesive and less ordered compared to the pure one, with the lipids assuming faster diffusion rates, whereas at high temperatures, the polymer interaction with the lipids acts reducing their diffusivity, but also increasing the lipid tail ordering and the membrane stiffness. The ability of the amphiphilic L64 copolymers to change the nature of the main phase transition of lipid membranes and their properties both in the ordered and the disordered phase is of vital importance for the prediction of both the efficacy of hybrid lipid/polymer nanoparticles as drug delivery vehicles as well as their potential adverse implications during interactions with healthy cell membranes.

Computational Modeling of Realistic Cell Membranes

Cell membranes contain a large variety of lipid types and are crowded with proteins, endowing them with the plasticity needed to fulfill their key roles in cell functioning. The compositional complexity of cellular membranes gives rise to a heterogeneous lateral organization, which is still poorly understood. Computational models, in particular molecular dynamics simulations and related techniques, have provided important insight into the organizational principles of cell membranes over the past decades. Now, we are witnessing a transition from simulations of simpler membrane models to multicomponent systems, culminating in realistic models of an increasing variety of cell types and organelles. Here, we review the state of the art in the field of realistic membrane simulations and discuss the current limitations and challenges ahead.

Shedding light on the structural properties of lipid bilayers using molecular dynamics simulation: a review study

In recent years, a massive increase has been observed in the number of published articles describing accurate and reliable molecular dynamics simulations of lipid bilayers. This is due to several reasons, including the development of fast and efficient methods for treating long-range electrostatic interactions, significant progress in computer hardware, progress in atomistic simulation algorithms and the development of well-validated empirical molecular mechanical force fields. Although molecular dynamics is an effective approach for investigating different aspects of lipid bilayers, to the best of our knowledge, there is no review in the literature that explains the different analyses that can be carried out with membrane simulation. This review gives an overview about the some of the most important possible analyses, technical challenges, and existing protocols that can be performed on the biological membrane by molecular dynamics simulation. The reviewed analyses include the degree of membrane disruption, average area per lipid, probability distributions for the area per lipid molecule, membrane thickness, membrane area compressibility, lateral diffusion, rotational diffusion, order parameters, head group tilt, electron density profile, mass density profile, electrostatic potential profile, ordering of vicinity waters, number of hydrogen bonds, and radial distribution function.

Multiscale Modeling and Simulation Approaches to Lipid-Protein Interactions

Lipid membranes play a crucial role in living systems by compartmentalizing biological processes and forming a barrier between these processes and the environment. Naturally, a large apparatus of biomolecules is responsible for construction, maintenance, transport, and degradation of these lipid barriers. Additional classes of biomolecules are tasked with transport of specific substances or transduction of signals from the environment across lipid membranes. In this article, we intend to describe a set of techniques that enable one to build accurate models of lipid systems and their associated proteins, and to simulate their dynamics over a variety of time and length scales. We discuss the methods and challenges that allow us to derive structural, mechanistic, and thermodynamic information from these models. We also show how these models have recently been applied in research to study some of the most complex lipid-protein systems to date, including bacterial and viral envelopes, neuronal membranes, and mammalian signaling systems.

Coarse-grained simulation studies on the adsorption of polyelectrolyte complexes upon lipid membranes

Conformations of a polyelectrolyte complex irreversibly bound to a zwitterionic lipid bilayer.

Epigallocatechin-3-gallate location and interaction with late endosomal and plasma membrane model membranes by molecular dynamics

Epigallocatechin-3-gallate (EGCG) is the most abundant polyphenol in green tea and it has been reported to have many beneficial properties against many different types of illnesses and infections. However, the exact mechanism/s underlying its biological effects are unknown. It has been previously shown that EGCG is capable of binding to and disrupting the membrane, so that some of its effects on biological systems could be ascribed to its capacity to incorporate into the biological membrane and modulate its structure. In this work, we have used atomistic molecular dynamics (MD) to discern the location and orientation of EGCG in model membranes and the possible existence of specific interactions with membrane lipids. For that goal, we have used in our simulation two complex model membranes, one resembling the plasma membrane (PM) and the other one the late endosome (LE) membrane. Our results support that EGCG tends to associate with the membrane and exists inside it in a relatively stable and steady location with a low propensity to be associated with other EGCG molecules. Interestingly, EGCG forms hydrogen bonds with POPC and POPE in the PM system but POPC and BMP and no POPE in the LE. These data suggest that the broad beneficial effects of EGCG could be mediated, at least in part, through its membranotropic effects and therefore membrane functioning. Communicated by Ramaswamy H. Sarma.

Testing for physical validity in molecular simulations

Advances in recent years have made molecular dynamics (MD) and Monte Carlo (MC) simulations powerful tools in molecular-level research, allowing the prediction of experimental observables in the study of systems such as proteins, membranes, and polymeric materials. However, the quality of any prediction based on molecular dynamics results will strongly depend on the validity of underlying physical assumptions. Unphysical behavior of simulations can have significant influence on the results and reproducibility of these simulations, such as folding of proteins and DNA or properties of lipid bilayers determined by cutoff treatment, dynamics of peptides and polymers affected by the choice of thermostat, or liquid properties depending on the simulation time step. Motivated by such examples, we propose a two-fold approach to increase the robustness of molecular simulations. The first part of this approach involves tests which can be performed by the users of MD programs on their respective systems and setups. We present a number of tests of different complexity, ranging from simple post-processing analysis to more involved tests requiring additional simulations. These tests are shown to significantly increase the reliability of MD simulations by catching a number of common simulation errors violating physical assumptions, such as non-conservative integrators, deviations from the Boltzmann ensemble, and lack of ergodicity between degrees of freedom. To make the usage as easy as possible, we have developed an open-source and platform-independent Python library (https://physical-validation.readthedocs.io) implementing these tests. The second part of the approach involves testing for code correctness. While unphysical behavior can be due to poor or incompatible choices of parameters by the user, it can just as well originate in coding errors within the program. We therefore propose to include physical validation tests in the code-checking mechanism of MD software packages. We have implemented such a validation for the GROMACS software package, ensuring that every major release passes a number of physical sanity checks performed on selected representative systems before shipping. It is, to our knowledge, the first major molecular mechanics software package to run such validation routinely. The tests are, as the rest of the package, open source software, and can be adapted for other software packages.

Atomistic-level study of the interactions between hIAPP protofibrils and membranes: Influence of pH and lipid composition

The pathology of type 2 diabetes mellitus is associated with the aggregation of human islet amyloid polypeptide (hIAPP) and aggregation-mediated membrane disruption. The interactions of hIAPP aggregates with lipid membrane, as well as the effects of pH and lipid composition at the atomic level, remain elusive. Herein, using molecular dynamics simulations, we investigate the interactions of hIAPP protofibrillar oligomers with lipids, and the membrane perturbation that they induce, when they are partially inserted in an anionic dipalmitoyl-phosphatidylglycerol (DPPG) membrane or a mixed dipalmitoyl-phosphatidylcholine (DPPC)/DPPG (7:3) lipid bilayer under acidic/neutral pH conditions. We observed that the tilt angles and insertion depths of the hIAPP protofibril are strongly correlated with the pH and lipid composition. At neutral pH, the tilt angle and insertion depth of hIAPP protofibrils at a DPPG bilayer reach ~52° and ~1.62 nm with respect to the membrane surface, while they become ~77° and ~1.75 nm at a mixed DPPC/DPPG membrane. The calculated tilt angle of hIAPP at DPPG membrane is consistent with a recent chiral sum frequency generation spectroscopic study. The acidic pH induces a smaller tilt angle of ~40° and a shallower insertion depth (~1.24 nm) of hIAPP at the DPPG membrane surface, mainly due to protonation of His18 near the turn region. These differences mainly result from a combination of distinct electrostatic, van der Waals, hydrogen bonding and salt-bridge interactions between hIAPP and lipid bilayers. The hIAPP-membrane interaction energy analysis reveals that besides charged residues K1, R11 and H18, aromatic residues Phe15 and Phe23 also exhibit strong interactions with lipid bilayers, revealing the crucial role of aromatic residues in stabilizing the membrane-bound hIAPP protofibrils. hIAPP-membrane interactions disturb the lipid ordering and the local bilayer thickness around the peptides. Our results provide atomic-level information of membrane interaction of hIAPP protofibrils, revealing pH-dependent and membrane-modulated hIAPP aggregation at the early stage.

Cell Membrane Disruption by Vertical Micro-/Nanopillars: Role of Membrane Bending and Traction Forces

Gaining access to the cell interior is fundamental for many applications, such as electrical recording and drug and biomolecular delivery. A very promising technique consists of culturing cells on micro-/nanopillars. The tight adhesion and high local deformation of cells in contact with nanostructures can promote the permeabilization of lipids at the plasma membrane, providing access to the internal compartment. However, there is still much experimental controversy regarding when and how the intracellular environment is targeted and the role of the geometry and interactions with surfaces. Consequently, we investigated, by coarse-grained molecular dynamics simulations of the cell membrane, the mechanical properties of the lipid bilayer under high strain and bending conditions. We found out that a high curvature of the lipid bilayer dramatically lowers the traction force necessary to achieve membrane rupture. Afterward, we experimentally studied the permeabilization rate of the cell membrane by pillars with comparable aspect ratios but different sharpness values at the edges. The experimental data support the simulation results: even pillars with diameters in the micron range may cause local membrane disruption when their edges are sufficiently sharp. Therefore, the permeabilization likelihood is connected to the local geometric features of the pillars rather than diameter or aspect ratio. The present study can also provide significant contributions to the design of three-dimensional biointerfaces for tissue engineering and cellular growth.

A fast and accurate computational method for the linear-combination-based isotropic periodic sum

An isotropic periodic sum (IPS) is a powerful technique to reasonably calculate intermolecular interactions for wide range of molecular systems under periodic boundary conditions. A linear-combination-based IPS (LIPS) has been developed to attain computational accuracy close to an exact lattice sum, such as the Ewald sum. The algorithm of the original LIPS method has a high computational cost because it needs long-range interaction calculations in real space. This becomes a performance bottleneck for long-time molecular simulations. In this work, the combination of an LIPS and fast Fourier transform (FFT) was developed, and evaluated on homogeneous and heterogeneous molecular systems. This combinational approach of LIPS/FFT attained computational efficiency close to that of a smooth particle mesh Ewald while maintaining the same high accuracy as the original LIPS. We concluded that LIPS/FFT has great potential to extend the capability of IPS techniques for the fast and accurate computation of many types of molecular systems.

Influence of Cholesterol on the Dynamics of Hydration in Phospholipid Bilayers

We investigate the dynamics of interfacial waters in dipalmitoylphosphatidylcholine (DPPC) bilayers upon the addition of cholesterol, by molecular dynamics simulations. Our data reveal that the inclusion of cholesterol modifies the membrane aqueous interfacial dynamics: waters diffuse faster, their rotational decay time is shorter, and the DPPC/water hydrogen bond dynamics relaxes faster than in the pure DPPC membrane. The observed acceleration of the translational water dynamics agrees with recent experimental results, in which, by means of NMR techniques, an increment of the surface water diffusivity is measured upon the addition of cholesterol. A microscopic analysis of the lipid/water hydrogen bond network at the interfacial region suggests that the mechanism underlying the observed water mobility enhancement is given by the rupture of a fraction of interlipid water bridge hydrogen bonds connecting two different DPPC molecules, concomitant to the formation of new lipid/solvent bonds, whose dynamics is faster than that of the former. The consideration of a simple two-state model for the decay of the hydrogen bond correlation function yielded excellent results, obtaining two well-separated characteristic time scales: a slow one (∼250 ps) associated with bonds linking two DPPC molecules, and a fast one (∼15 ps), related to DPPC/solvent bonds.