CHARMM at 45: Enhancements in Accessibility, Functionality, and Speed.

Since its inception nearly a half century ago, CHARMM has been playing a central role in computational biochemistry and biophysics. Commensurate with the developments in experimental research and advances in computer hardware, the range of methods and applicability of CHARMM have also grown. This review summarizes major developments that occurred after 2009 when the last review of CHARMM was published. They include the following: new faster simulation engines, accessible user interfaces for convenient workflows, and a vast array of simulation and analysis methods that encompass quantum mechanical, atomistic, and coarse-grained levels, as well as extensive coverage of force fields. In addition to providing the current snapshot of the CHARMM development, this review may serve as a starting point for exploring relevant theories and computational methods for tackling contemporary and emerging problems in biomolecular systems. CHARMM is freely available for academic and nonprofit research at https://academiccharmm.org/program.

Metadynamics Study of Lipid-Mediated Antibacterial Toxin Binding to the EmrE Multiefflux Protein.

EmrE is a bacterial efflux protein in the small multidrug-resistant (SMR) family present in Escherichia coli. Due to its small size, 110 residues in each dimer subunit, it is an ideal model system to study ligand-protein-membrane interactions. Here in our work, we have calculated the free energy landscape of benzyltrimetylammonium (BTMA) and tetraphenyl phosphonium (TPP) binding to EmrE using the enhanced sampling method-multiple walker metadynamics. We estimate that the free energy of BTMA binding to EmrE is -21.2 ± 3.3 kJ/mol and for TPP is -43.6 ± 3.8 kJ/mol. BTMA passes through two metastable states to reach the binding pocket, while TPP has a more complex binding landscape with four metastable states and one main binding site. Our simulations show that the ligands interact with the membrane lipids at a distance 1 nm away from the binding site which forms a broad local minimum, consistent for both BTMA and TPP. This site can be an alternate entry point for ligands to partition from the membrane into the protein, especially for bulky and/or branched ligands. We also observed the membrane lipid and C-terminal 110HisA form salt-bridge interactions with the helix-1 residue 22LysB. Our free energy estimates and clusters are in close agreement with experimental data and give us an atomistic view of the ligand-protein-lipid interactions. Understanding the binding pathway of these ligands can guide us in future design of ligands that can alter or halt the function of EmrE.

Investigation of the Effect of Peptide p5 Targeting CDK5-p25 Hyperactivity on Munc18-1 (P67) Regulating Neuronal Exocytosis Using Molecular Simulations.

Munc18-1 is an SM (sec1/munc-like) family protein involved in vesicle fusion and neuronal exocytosis. Munc18-1 is known to regulate the exocytosis process by binding with closed- and open-state conformations of Syntaxin1, a protein belonging to the SNARE family established to be central to the exocytosis process. Our previous work studied peptide p5 as a promising drug candidate for CDK5-p25 complex, an Alzheimer's disease (AD) pathological target. Experimental in vivo and in vitro studies suggest that Munc18-1 promotes p5 to selectively inhibit the CDK5-p25 complex without affecting the endogenous CDK5 activity, a characteristic of remarkable therapeutic implications. In this paper, we identify several binding modes of p5 with Munc18-1 that could potentially affect the Munc18-1 binding with SNARE proteins and lead to off-target effects on neuronal communication using molecular dynamics simulations. Recent studies indicate that disruption of Munc18-1 function not only disrupts neurotransmitter release but also results in neurodegeneration, exhibiting clinical resemblance to other neurodegenerative conditions such as AD, causing diagnostic and treatment challenges. We characterize such interactions between p5 and Munc18-1, define the corresponding pharmacophores, and provide guidance for the in vitro validation of our findings to improve therapeutic efficacy and safety of p5.

Modeling asymmetric cell membranes at all-atom resolution.

Molecular dynamics (MD) simulations are a useful tool when studying the properties of membranes as they allow for a molecular view of lipid interactions with proteins, nucleic acids, or small molecules. While model membranes are usually symmetric in their lipid composition between leaflets and include a small number of lipid components, physiological membranes are highly complex and vary in the level of asymmetry. Simulation studies have shown that changes in leaflet asymmetry can alter the properties of a membrane. It is therefore necessary to carefully build asymmetric membranes to accurately simulate membranes. This chapter carefully describes the different methods for building asymmetric membranes and the advantages/disadvantages of each method. The simplest methods involve building a membrane with either an equal number of lipids per leaflet or an equal initial surface area (SA) estimated by the area per lipid. More detailed methods include combining two symmetric membranes of equal SA or altering an asymmetric membrane and adjusting the number of lipids after equilibration to minimize an observable such as differential stress (0-DS). More complex methods that require specific simulation software are also briefly described. The challenges and assumptions are listed for each method which should help guide the researcher to choose the best method for their unique MD simulation of an asymmetric membrane.

Proposed dual membrane contact with full-length Osh4.

Membrane contacts sites (MCSs) play important roles in lipid trafficking across cellular compartments and maintain the widespread structural diversity of organelles. We have utilized microsecond long all-atom (AA) molecular dynamics (MD) simulations and enhanced sampling techniques to unravel the MCS structure targeting by yeast oxysterol binding protein (Osh4) in an environment that mimics the interface of membranes with an increased proportion of anionic lipids using CHARMM36m forcefield with additional CUFIX parameters for lipid-protein electrostatic interactions. In a dual-membrane environment, unbiased MD simulations show that Osh4 briefly interacts with both membranes, before aligning itself with a single membrane, adopting a ß-crease-bound conformation similar to observations in a single-membrane scenario. Targeted molecular dynamics simulations followed by microsecond-long AA MD simulations have revealed a distinctive dual-membrane bound state of Osh4 at MCS, wherein the protein interacts with the lower membrane via the ß-crease surface, featuring its PHE-239 residue positioned below the phosphate plane of membrane, while concurrently establishing contact with the opposite membrane through the extended α6-α7 region. Osh4 maintains these dual membrane contacts simultaneously over the course of microsecond-long MD simulations. Moreover, binding energy calculations highlighted the essential roles played by the phenylalanine loop and the α6 helix in dynamically stabilizing dual-membrane bound state of Osh4 at MCS. Our computational findings were corroborated through frequency of contact analysis, showcasing excellent agreement with past experimental cross-linking data. Our computational study reveals a dual-membrane bound conformation of Osh4, providing insights into protein-membrane interactions at membrane contact sites and their relevance to lipid transfer processes.

Expanding the CHARMM36 United Atom Chain Model for the Inclusion of Sphingolipids.

The inclusion of accurate yet computationally inexpensive lipid force fields (FF) is pertinent for the study of lipids and lipid-containing systems using molecular dynamics (MD). Within the past decade, the implementation and further expansion of a united atom (UA) FF for lipids have been developed in the CHARMM family of FFs. The most recent version of the UA presented more accurate descriptions of lipid properties for several phospholipids with saturated and monounsaturated chains, termed C36UAr. However, the original C36UAr model lacks parameters for an important class of lipids, such as sphingolipids. The focus of this article is to broaden the scope of the C36UAr chain model to incorporate these lipids. In this study, two common sphingolipids, N-palmitoyl sphingomyelin and N-stearoyl sphingomyelin are converted to a UA-chain representation and simulated to investigate the accuracy and speed over the all-atom FF model for sphingolipids. Improvements were found among multiple parameters, for example, in the surface area per lipid (SA/lip) and hydrogen order parameters, over the all-atom simulations of these sphingomyelins in C36, while as much as halving the simulation time for simulations of the same setup otherwise. Thus, the accuracy and efficiency found in this study are consistent with those found in the C36UAr model for phospholipids and expand the application of C36UAr to a wider array of membrane models to better match that in vivo.

Chlorpromazine inhibits EAG1 channels by altering the coupling between the PAS, CNBH and pore domains.

EAG1 depolarization-activated potassium selective channels are important targets for treatment of cancer and neurological disorders. EAG1 channels are formed by a tetrameric subunit assembly with each subunit containing an N-terminal Per-Arnt-Sim (PAS) domain and C-terminal cyclic nucleotide-binding homology (CNBH) domain. The PAS and CNBH domains from adjacent subunits interact and form an intracellular tetrameric ring that regulates the EAG1 channel gating, including the movement of the voltage sensor domain (VSD) from closed to open states. Small molecule ligands can inhibit EAG1 channels by binding to their PAS domains. However, the allosteric pathways of this inhibition are not known. Here we show that chlorpromazine, a PAS domain small molecule binder, alters interactions between the PAS and CNBH domains and decreases the coupling between the intracellular tetrameric ring and the pore of the channel, while having little effect on the coupling between the PAS and VSD domains. In addition, chlorpromazine binding to the PAS domain did not alter Cole-Moore shift characteristic of EAG1 channels, further indicating that chlorpromazine has no effect on VSD movement from the deep closed to opened states. Our study provides a framework for understanding global pathways of EAG1 channel regulation by small molecule PAS domain binders.

Multiscale Molecular Dynamics Simulations of the Homodimer Accessory Protein ORF7b of SARS-CoV-2.

The SARS-CoV-2 ORF7b protein has drawn attention for its potential role in viral pathogenesis, but its structural details and lateral membrane associations remain elusive. In this study, we conducted multiscale molecular dynamics simulations to provide detailed molecular insights of the protein's dimerization, which is crucial for unraveling its structural model of protein-protein interface important to regulating cellular immune response. To gain a deeper understanding of homodimer configurations, we employed a machine learning algorithm for structural-based clustering. Clusters were categorized into three distinct groups for both parallel and antiparallel orientations, highlighting the influence of the initial monomer conformation on dimer configurations. Analysis of hydrogen bonding and π-π and π-cation stacking interactions within clusters revealed variations in interactions between clusters. In parallel dimers, weak stacking interactions in the transmembrane (TM) region were observed. In contrast, antiparallel dimers exhibited strong hydrogen bonding and stacking interactions contributing to tight dimeric packing, both within and outside the TM domain. Overall, our study provides a comprehensive view of the structural dynamics of ORF7b homodimerization in both parallel and antiparallel orientations. These findings shed light on the molecular interactions involved in ORF7b dimerization, which are crucial for understanding its potential roles in SARS-CoV-2 pathogenesis. This knowledge could inform future research and therapeutic strategies targeting this viral protein.

Molecular mechanism of EAG1 channel inhibition by imipramine binding to the PAS domain.

Ether-a-go-go (EAG) channels are key regulators of neuronal excitability and tumorigenesis. EAG channels contain an N-terminal Per-Arnt-Sim (PAS) domain that can regulate currents from EAG channels by binding small molecules. The molecular mechanism of this regulation is not clear. Using surface plasmon resonance and electrophysiology we show that a small molecule ligand imipramine can bind to the PAS domain of EAG1 channels and inhibit EAG1 currents via this binding. We further used a combination of molecular dynamics (MD) simulations, electrophysiology, and mutagenesis to investigate the molecular mechanism of EAG1 current inhibition by imipramine binding to the PAS domain. We found that Tyr71, located at the entrance to the PAS domain cavity, serves as a "gatekeeper" limiting access of imipramine to the cavity. MD simulations indicate that the hydrophobic electrostatic profile of the cavity facilitates imipramine binding and in silico mutations of hydrophobic cavity-lining residues to negatively charged glutamates decreased imipramine binding. Probing the PAS domain cavity-lining residues with site-directed mutagenesis, guided by MD simulations, identified D39 and R84 as residues essential for the EAG1 channel inhibition by imipramine binding to the PAS domain. Taken together, our study identified specific residues in the PAS domain that could increase or decrease EAG1 current inhibition by imipramine binding to the PAS domain. These findings should further the understanding of molecular mechanisms of EAG1 channel regulation by ligands and facilitate the development of therapeutic agents targeting these channels.

Conformational Fluctuations in ß2-Microglubulin Using Markov State Modeling and Molecular Dynamics.

Conformational dynamics in proteins can give rise to aggregation prone states during folding, and these kinetically stable states could form oligomers and aggregates. In this study, we investigate the intermediate states and near-folded states of ß2-microglobulin and their physico-chemical properties using molecular dynamics and Markov state modeling. Analysis of hundreds of microseconds simulation show the importance of the edge strands in the misfolded states that give rise to a high exposure of hydrophobic residues in the core of the protein that could initiate oligomerization and aggregate formation. Our study sheds light on the first step of aggregation of ß2m monomers and gave a better picture of the landscape of protein misfolding and aggregation.

Understanding how transmembrane domains regulate interactions between human BST-2 and the SARS-CoV-2 accessory protein ORF7a.

Severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), the causative agent of COVID, replicates at intracellular membranes. Bone marrow stromal antigen 2 (BST-2; tetherin) is an antiviral response protein that inhibits transport of viral particles after budding within infected cells. RNA viruses such as SARS-CoV-2 use various strategies to disable BST-2, including use of transmembrane 'accessory' proteins that interfere with BST-2 oligomerization. ORF7a is a small, transmembrane protein present in SARS-CoV-2 shown previously to alter BST-2 glycosylation and function. In this study, we investigated the structural basis for BST-2 ORF7a interactions, with a particular focus on transmembrane and juxtamembrane interactions. Our results indicate that transmembrane domains play an important role in BST-2 ORF7a interactions and mutations to the transmembrane domain of BST-2 can alter these interactions, particularly single-nucleotide polymorphisms in BST-2 that result in mutations such as I28S. Using molecular dynamics simulations, we identified specific interfaces and interactions between BST-2 and ORF7a to develop a structural basis for the transmembrane interactions. Differences in glycosylation are observed for BST-2 transmembrane mutants interacting with ORF7a, consistent with the idea that transmembrane domains play a key role in their heterooligomerization. Overall, our results indicate that ORF7a transmembrane domain interactions play a key role along with extracellular and juxtamembrane domains in modulating BST-2 function.

Drude Polarizable Lipid Force Field with Explicit Treatment of Long-Range Dispersion: Parametrization and Validation for Saturated and Monounsaturated Zwitterionic Lipids.

Accurate empirical force fields of lipid molecules are a critical component of molecular dynamics simulation studies aimed at investigating properties of monolayers, bilayers, micelles, vesicles, and liposomes, as well as heterogeneous systems, such as protein-membrane complexes, bacterial cell walls, and more. While the majority of lipid force field-based simulations have been performed using pairwise-additive nonpolarizable models, advances have been made in the development of the polarizable force field based on the classical Drude oscillator model. In the present study, we undertake further optimization of the Drude lipid force field, termed Drude2023, including improved treatment of the phosphate and glycerol linker region of PC and PE headgroups, additional optimization of the alkene group in monounsaturated lipids, and inclusion of long-range Lennard-Jones interactions using the particle-mesh Ewald method. Initial optimization targeted quantum mechanical (QM) data on small model compounds representative of the linker region. Subsequent optimization targeted QM data on larger model compounds, experimental data, and dihedral potentials of mean force from the CHARMM36 additive lipid force field using a parameter reweighting protocol. The use of both experimental and QM target data during the reweighting protocol is shown to produce physically reasonable parameters that reproduce a collection of experimental observables. Target data for optimization included surface area/lipid for DPPC, DSPC, DMPC, and DLPC bilayers and nuclear magnetic resonance (NMR) order parameters for DPPC bilayers. Validation data include prediction of membrane thickness, scattering form factors, electrostatic potential profiles, compressibility moduli, surface area per lipid, water permeability, NMR T1 relaxation times, diffusion constants, and monolayer surface tensions for a variety of saturated and unsaturated lipid mono- and bilayers. Overall, the agreement with experimental data is quite good, though the results are less satisfactory for the NMR T1 relaxation times for carbons near the ester groups. Notable improvements compared to the additive C36 force field were obtained for membrane dipole potentials, lipid diffusion coefficients, and water permeability with the exception of monounsaturated lipid bilayers. It is anticipated that the optimized polarizable Drude2023 force field will help generate more accurate molecular simulations of pure bilayers and heterogeneous systems containing membranes, advancing our understanding of the role of electronic polarization in these systems.

Antimicrobial Peptide Mechanism Studied by Scattering-Guided Molecular Dynamics Simulation.

In an effort to combat rising antimicrobial resistance, our labs have rationally designed cationic, helical, amphipathic antimicrobial peptides (AMPs) as alternatives to traditional antibiotics since AMPs incur bacterial resistance in weeks, rather than days. One highly positively charged AMP, WLBU2 (+13e), (RRWV RRVR RWVR RVVR VVRR WVRR), has been shown to be effective in killing both Gram-negative (G(-)) and Gram-positive (G(+)) bacteria by directly perturbing the bacterial membrane nonspecifically. Previously, we used two equilibrium experimental methods: synchrotron X-ray diffuse scattering (XDS) providing lipid membrane thickness and neutron reflectometry (NR) providing WLBU2 depth of penetration into three lipid model membranes (LMMs). The purpose of the present study is to use the results from the scattering experiments to guide molecular dynamics (MD) simulations to investigate the detailed biophysics of the interactions of WLBU2 with LMMs of Gram-negative outer and inner membranes, and Gram-positive cell membranes, to elucidate the mechanisms of bacterial killing. Instead of coarse-graining, backmapping, or simulating without bias for several microseconds, all-atom (AA) simulations were guided by the experimental results and then equilibrated for ∼0.5 µs. Multiple replicas of the inserted peptide were run to probe stability and reach a combined time of at least 1.2 µs for G(-) and also 2.0 µs for G(+). The simulations with experimental comparisons help rule out certain structures and orientations and propose the most likely set of structures, orientations, and effects on the membrane. The simulations revealed that water, phosphates, and ions enter the hydrocarbon core when WLBU2 is positioned there. For an inserted peptide, the three types of amino acids, arginine, tryptophan, and valine (R, W, V), are arranged with the 13 Rs extending from the hydrocarbon core to the phosphate group, Ws are located at the interface, and Vs are more centrally located. For a surface state, R, W, and V are positioned relative to the bilayer interface as expected from their hydrophobicities, with Rs closest to the phosphate group, Ws close to the interface, and Vs in between. G(-) and G(+) LMMs are thinned ∼1 Å by the addition of WLBU2. Our results suggest a dual anchoring mechanism for WLBU2 both in the headgroup and in the hydrocarbon region that promotes a defect region where water and ions can flow across the slightly thinned bacterial cell membrane.

Molecular dynamics simulations of the human ocular lens with age and cataract.

The human ocular lens consists primarily of elongated, static fibers characterized by high stability and low turnover, which differ dramatically in their composition and properties from other biological membranes. Cholesterol (Chol) and sphingolipids (SL) are present at high concentrations, including saturated SLs, such as dihyrosphingomyelin (DHSM). Past molecular dynamics simulations demonstrated that the presence of DHSM and high Chol concentration contributes to higher order in lipid membranes. This current study simulated more complex models of human lens membranes. Models were developed representing physiological compositions in cataractous lenses aged 74 ± 6 years and in healthy lenses aged 22 ± 4, 41 ± 6, and 69 ± 3 years. With older age, Chol and ceramide concentrations increase and glycerophospholipid concentration decreases. With cataract, ceramide concentration increases and Chol and glycerophospholipid concentrations decrease. Surface area per lipid, deuterium order parameters (SCD), sterol tilt angle, electron density profiles, bilayer thickness, chain interdigitation, two-dimensional radial distribution functions (2D-RDF), lipid clustering, and hydrogen bonding were calculated for all simulations. All systems exhibited low surface area per lipid and high bilayer thickness, indicative of strong vertical packing. SCD parameters suggest similarly, with saturated tails in the hydrophobic core of the membrane having elevated order. Vertical packing and acyl tail order increased with both age and cataract condition. Lateral diffusion decreased with age and cataracts, with the older and cataractous models demonstrating increased long-range structure by the 2D-RDF analysis. In future work examining the membrane proteins of the lens, these models can serve as a physiologically accurate representation of the lens lipidome.

Computational Study of the Allosteric Effects of p5 on CDK5-p25 Hyperactivity as Alternative Inhibitory Mechanisms in Neurodegeneration.

The cyclin-dependent kinase (CDK5) forms a stable complex with its activator p25, leading to the hyperphosphorylation of tau proteins and to the formation of plaques and tangles that are considered to be one of the typical causes of Alzheimer's disease (AD). Hence, the pathological CDK5-p25 complex is a promising therapeutic target for AD. Small peptides, obtained from the truncation of CDK5 physiological activator p35, have shown promise in inhibiting the pathological complex effectively while also crossing the blood-brain barrier. One such small 24-residue peptide, p5, has shown selective inhibition toward the pathological complex in vivo. Our previous research focused on the characterization of a computationally predicted CDK5-p5 binding mode and of its pharmacophore, which was consistent with competitive inhibition. In continuation of our previous work, herein, we investigate four additional binding modes to explore other possible mechanisms of interaction between CDK5 and p5. The quantitative description of the pharmacophore is consistent with both competitive and allosteric p5-induced inhibition mechanisms of CDK5-p25 pathology. The gained insights can direct further in vivo/in vitro tests and help design small peptides, linear or cyclic, or peptidomimetic compounds as adjuvants of orthosteric inhibitors or as part of a cocktail of drugs with enhanced effectiveness and lower side effects.

A reengineered common chain cytokine augments CD8+ T cell-dependent immunotherapy.

Cytokine therapy is limited by undesirable off-target side effects as well as terminal differentiation and exhaustion of chronically stimulated T cells. Here, we describe the signaling properties of a potentially unique cytokine by design, where T cell surface binding and signaling are separated between 2 different families of receptors. This fusion protein cytokine, called OMCPmutIL-2, bound with high affinity to the cytotoxic lymphocyte-defining immunoreceptor NKG2D but signaled through the common γ chain cytokine receptor. In addition to precise activation of cytotoxic T cells due to redirected binding, OMCPmutIL-2 resulted in superior activation of both human and murine CD8+ T cells by improving their survival and memory cell generation and decreasing exhaustion. This functional improvement was the direct result of altered signal transduction based on the reorganization of surface membrane lipid rafts that led to Janus kinase-3-mediated phosphorylation of the T cell receptor rather than STAT/AKT signaling intermediates. This potentially novel signaling pathway increased CD8+ T cell response to low-affinity antigens, activated nuclear factor of activated T cells transcription factors, and promoted mitochondrial biogenesis. OMCPmutIL-2 thus outperformed other common γ chain cytokines as a catalyst for in vitro CD8+ T cell expansion and in vivo CD8+ T cell-based immunotherapy.

GraphVAMPNet, using graph neural networks and variational approach to Markov processes for dynamical modeling of biomolecules.

Finding a low dimensional representation of data from long-timescale trajectories of biomolecular processes, such as protein folding or ligand-receptor binding, is of fundamental importance, and kinetic models, such as Markov modeling, have proven useful in describing the kinetics of these systems. Recently, an unsupervised machine learning technique called VAMPNet was introduced to learn the low dimensional representation and the linear dynamical model in an end-to-end manner. VAMPNet is based on the variational approach for Markov processes and relies on neural networks to learn the coarse-grained dynamics. In this paper, we combine VAMPNet and graph neural networks to generate an end-to-end framework to efficiently learn high-level dynamics and metastable states from the long-timescale molecular dynamics trajectories. This method bears the advantages of graph representation learning and uses graph message passing operations to generate an embedding for each datapoint, which is used in the VAMPNet to generate a coarse-grained dynamical model. This type of molecular representation results in a higher resolution and a more interpretable Markov model than the standard VAMPNet, enabling a more detailed kinetic study of the biomolecular processes. Our GraphVAMPNet approach is also enhanced with an attention mechanism to find the important residues for classification into different metastable states.

Modeling the membrane binding mechanism of a lipid transport protein Osh4 to single membranes.

All-atom (AA) molecular dynamics simulations are used to unravel the binding mechanism of yeast oxysterol binding protein (Osh4) to model membranes with varying anionic lipid concentration using AA and the highly mobile membrane mimetic (HMMM) representations. For certain protein-lipid interactions, an improved forcefield description is used (CUFIX) to accurately describe lipid-protein electrostatic interactions. Our detailed computational studies have identified a single, ß-crease oriented, membrane-bound conformation of Osh4 for all anionic membranes. The penetration of the PHE-239 residue below the membrane phosphate plane is the characteristic signature of the membrane-bound state of Osh4. As the phenylalanine loop anchors itself deeply in the membrane; the other regions of the Osh4, namely, ALPS motif, ß6- ß7 loop, ß14- ß15 loop, and ß16- ß17 loop, maximize their contact with the membrane. Furthermore, loose lipid packing and higher mobility of HMMM enable stronger association of the ALPS motif with the membrane lipids through its hydrophobic surface. After the HMMM is converted to AA and equilibrated, the binding is two to three times stronger compared with simulations started with the AA representation, yielding the major importance of the ALPS motif to binding. Quantitative estimation of binding energy revealed that the phenylalanine loop plays a crucial role in stable membrane attachment of Osh4 and contributes significantly toward overall binding process. The CUFIX parameters provide a more balanced picture of hydrophobic and electrostatic interactions between the protein and the membrane, which differs from our past work that showed salt bridges alone stabilized Osh4-membrane contact. Our study provides a comprehensive picture of the binding mechanism of Osh4 with model single membranes and, thus, understanding of the initial interactions is important for elucidating the biological function of this protein to shuttle lipids between organelles.

Leaflet Asymmetry Modeling in the Lipid Composition of Escherichia coli Cytoplasmic Membranes.

Lipid composition asymmetry between leaflets is important to cell function and plays a key role in the "positive inside" rule in transmembrane proteins. In this work, Escherichia coli inner plasma membrane models reflecting this asymmetry have been investigated at the early-log and stationary stages during the bacterial lifecycle using all-atom molecular dynamics simulations. The CHARMM36 lipid force field is used, and selected membrane properties are tested for variations between two leaflets and whole membranes. Our models include bacterial lipids with a cyclopropane moiety on the sn-2 acyl chain in the stationary membrane model. The PE/PG ratio for two leaflets reflects the "positive inside" rule of membrane proteins, set to 6.8 and 2.8 for the inner and outer leaflets of the two models, respectively. We are the first to model leaflet asymmetry in the lipid composition of E. coli cytoplasmic membranes and observe the effect on membrane properties in leaflets and whole membranes. Specifically, our results show that for the stationary phase bilayer, the surface area per lipid (SA/lipid) is larger, the thickness (2DC and DB) is smaller, the tilt angle is larger, the tilt modulus is smaller, and the deuterium order parameters (SCD) of sn-1 and sn-2 tails are lower, compared to the early-log stage. Moreover, the stationary stage bilayer has a positive spontaneous curvature, while the early-log stage has a near flat spontaneous curvature. For leaflet asymmetry, the inner leaflet has a larger SA/lipid, a smaller thickness, a smaller elastic tilt modulus (a larger tilt angle), and lower SCD, compared to the outer leaflet in both stages. Moreover, an asymmetric membrane involves a lipid tilt and a lateral extension, varying from a reference state of a pre-equilibrium membrane. This work encourages a more profound exploration of leaflet asymmetry in various other membrane models and how this might affect the structure and function of membrane-associated peptides and proteins.

All-Atom Modeling of Complex Cellular Membranes.

Cell membranes are composed of a variety of lipids and proteins where they interact with each other to fulfill their roles. The first step in modeling these interactions in molecular simulations is to have reliable mimetics of the membrane's lipid environment. This Feature Article presents our recent efforts to model complex cellular membranes using all-atom force fields. A short review of the CHARMM36 (C36) lipid force field and its recent update to incorporate the long-range dispersion is presented. Key examples of model membranes mimicking various species and organelles are given. These include single-celled organisms such as bacteria (E. coli., chlamydia, and P. aeruginosa) and yeast (plasma membrane, endoplasmic reticulum, and trans-Golgi network) and more advanced ones such as plants (soybean and Arabidopsis thaliana) and mammals (ocular lens, stratum corneum, and peripheral nerve myelin). Leaflet asymmetry in composition has also been applied to some of these models. With the increased lipid diversity in the C36 lipid FF, these complex models can better reflect the structural, mechanical, and dynamic properties of realistic membranes and open an opportunity to study biological processes involving other molecules.