Igor Vorobyov

CHARMM general force field: A force field for drug‐like molecules compatible with the CHARMM all‐atom additive biological force fields

The widely used CHARMM additive all‐atom force field includes parameters for proteins, nucleic acids, lipids, and carbohydrates. In the present article, an extension of the CHARMM force field to drug‐like molecules is presented. The resulting CHARMM General Force Field (CGenFF) covers a wide range of chemical groups present in biomolecules and drug‐like molecules, including a large number of heterocyclic scaffolds. The parametrization philosophy behind the force field focuses on quality at the expense of transferability, with the implementation concentrating on an extensible force field. Statistics related to the quality of the parametrization with a focus on experimental validation are presented. Additionally, the parametrization procedure, described fully in the present article in the context of the model systems, pyrrolidine, and 3‐phenoxymethylpyrrolidine will allow users to readily extend the force field to chemical groups that are not explicitly covered in the force field as well as add functional groups to and link together molecules already available in the force field. CGenFF thus makes it possible to perform “all‐CHARMM” simulations on drug‐target interactions thereby extending the utility of CHARMM force fields to medicinally relevant systems.

Update of the CHARMM all-atom additive force field for lipids: Validation on six lipid types

A significant modification to the additive all-atom CHARMM lipid force field (FF) is developed and applied to phospholipid bilayers with both choline and ethanolamine containing head groups and with both saturated and unsaturated aliphatic chains. Motivated by the current CHARMM lipid FF (C27 and C27r) systematically yielding values of the surface area per lipid that are smaller than experimental estimates and gel-like structures of bilayers well above the gel transition temperature, selected torsional, Lennard-Jones and partial atomic charge parameters were modified by targeting both quantum mechanical (QM) and experimental data. QM calculations ranging from high-level ab initio calculations on small molecules to semiempirical QM studies on a 1,2-dipalmitoyl-sn-phosphatidylcholine (DPPC) bilayer in combination with experimental thermodynamic data were used as target data for parameter optimization. These changes were tested with simulations of pure bilayers at high hydration of the following six lipids: DPPC, 1,2-dimyristoyl-sn-phosphatidylcholine (DMPC), 1,2-dilauroyl-sn-phosphatidylcholine (DLPC), 1-palmitoyl-2-oleoyl-sn-phosphatidylcholine (POPC), 1,2-dioleoyl-sn-phosphatidylcholine (DOPC), and 1-palmitoyl-2-oleoyl-sn-phosphatidylethanolamine (POPE); simulations of a low hydration DOPC bilayer were also performed. Agreement with experimental surface area is on average within 2%, and the density profiles agree well with neutron and X-ray diffraction experiments. NMR deuterium order parameters (S(CD)) are well predicted with the new FF, including proper splitting of the S(CD) for the aliphatic carbon adjacent to the carbonyl for DPPC, POPE, and POPC bilayers. The area compressibility modulus and frequency dependence of (13)C NMR relaxation rates of DPPC and the water distribution of low hydration DOPC bilayers also agree well with experiment. Accordingly, the presented lipid FF, referred to as C36, allows for molecular dynamics simulations to be run in the tensionless ensemble (NPT), and is anticipated to be of utility for simulations of pure lipid systems as well as heterogeneous systems including membrane proteins.

Determination of Electrostatic Parameters for a Polarizable Force Field Based on the Classical Drude Oscillator

A procedure to determine the electrostatic parameters has been developed for a polarizable empirical force field based on the classical Drude oscillator model. Atomic charges and polarizabilities for a given molecule of interest were derived from restrained fitting to quantum-mechanical electrostatic potentials (ESP) calculated at the B3LYP/ cc-pVDZ or B3LYP/aug-cc-pVDZ levels on grid points located on concentric Connolly surfaces. The determination of the atomic polarizabilities requires a series of perturbed ESP maps, each one representing the electronic response of the molecule in the presence of a background charge placed on Connolly surfaces primarily along chemical bonds and lone pairs. Reference values for the partial atomic charges were taken from the CHARMM27 additive all-atom force field, and those for the polarizabilities were based on adjusted Millers ahp atomic polarizability values. The fitted values of atomic polarizabilities were scaled to reflect the reduced polarization expected for the condensed media and/or to correct for the systematic underestimation of experimental molecular polarizabilities by B3LYP calculations. Following correction of the polarizabilities, the atomic charges were adjusted to reproduce gas-phase dipole moments. The developed scheme has been tested on a set of small molecules representing functional moieties of nucleic acids. The derived electrostatic parameters have been successfully applied in a preliminary polarizable molecular dynamics simulation of a DNA octamer in a box of water with sodium counterions. Thus, this study confirms the feasibility of the use of a polarizable force field based on a classical Drude model for simulations of biomolecules in the condensed phase.

Additive and classical drude polarizable force fields for linear and cyclic ethers

Empirical force field parameters consistent with the CHARMM additive and classical Drude based polarizable force fields are presented for linear and cyclic ethers. Initiation of the optimization process involved validation of the aliphatic parameters based on linear alkanes and cyclic alkanes. Results showed the transfer to cyclohexane to yield satisfactory agreement with target data; however, in the case of cyclopentane direct transfer of the Lennard-Jones parameters was not sufficient due to ring strain, requiring additional optimization of these parameters for this molecule. Parameters for the ethers were then developed starting with the available aliphatic parameters, with the nonbond parameters for the oxygens optimized to reproduce both gas- and condensed-phase properties. Nonbond parameters for the polarizable model include the use of an anisotropic electrostatic model on the oxygens. Parameter optimization emphasized the development of transferable parameters between the ethers of a given class. The ether models are shown to be in satisfactory agreement with both pure solvent and aqueous solvation properties, and the resulting parameters are transferable to test molecules. The presented force field will allow for simulation studies of ethers in condensed phase and provides a basis for ongoing developments in both additive and polarizable force fields for biological molecules.

Potential of mean force and pKa profile calculation for a lipid membrane-exposed arginine side chain

The issue of ionizable protein side chains interacting with lipid membranes has been the focus of much attention since the proposal of the paddle model of voltage-gated ion channels, which suggested multiple arginine (Arg) side chains may move through the hydrocarbon core of a lipid membrane. Recent cell biology experiments have also been interpreted to suggest that these side chains would face only small free energy penalties to cross membranes, challenging a long-standing view in membrane biophysics. Here, we employ side chain analog and transmembrane helix models to determine the free energy of an Arg side chain, as a function of protonation state, across a membrane. We observe high free energy barriers for both the charged and neutral states that would prohibit lipid-exposed movement. The mechanisms for charged and neutral Arg transport are, however, very different, with the neutral state experiencing simple dehydration, whereas the charged state experiences a complex mechanism involving connections to the bilayer interfaces that deform the local membrane structure. We employ special methods to ensure sampling of these interfacial connections and decompose the free energy to shed light on the mechanisms. These deformations are found to preferentially stabilize the protonated form, such that the Arg side chain remains almost exclusively charged inside the membrane, with a pKa shift of <or=4.5 units. In contrast, the analog models are found to exaggerate the variations in energetics across the membrane and have larger pKa shifts. These results have implications for models of voltage gated ion channels, suggesting that although Arg side chains are ideally suited for carrying charge, the thermodynamics dictate that they must remain sequestered from the lipid bilayer environment.

Assessing atomistic and coarse-grained force fields for protein-lipid interactions: the formidable challenge of an ionizable side chain in a membrane.

Ionizable amino acid side chains play important roles in membrane protein structure and function, including the activation of voltage-gated ion channels, where it has been previously suggested that charged side chains may move through the hydrocarbon core of the membrane. However, all-atom molecular dynamics simulations have demonstrated large free energy barriers for such lipid-exposed motions. These simulations have also revealed that the membrane will deform due to the presence of a charged side chain, leading to a complex solvation microenvironment for which empirical force fields were not specifically parametrized. We have tested the ability of the all-atom CHARMM, Drude polarizable CHARMM, and a recent implementation of a coarse-grained force field to measure the thermodynamics of arginine-membrane interactions as a function of protonation state. We have employed model systems to attempt to match experimental bulk partitioning and quantum mechanical interactions within the membrane and found that free energy profiles from nonpolarizable and polarizable CHARMM simulations are accurate to within 1-2 kcal/mol. In contrast, the coarse-grained simulations failed to reproduce the same membrane deformations, exhibit interactions that are an order of magnitude too small, and thus, have incorrect free energy profiles. These results illustrate the need for careful parametrization of coarse-grained force fields and demonstrate the utility of atomistic molecular dynamics for providing quantitative thermodynamic and mechanistic analysis of protein-lipid interactions.

The different interactions of lysine and arginine side chains with lipid membranes.

The basic amino acids lysine (Lys) and arginine (Arg) play important roles in membrane protein activity, the sensing of membrane voltages, and the actions of antimicrobial, toxin, and cell-penetrating peptides. These roles are thought to stem from the strong interactions and disruptive influences of these amino acids on lipid membranes. In this study, we employ fully atomistic molecular dynamics simulations to observe, quantify, and compare the interactions of Lys and Arg with saturated phosphatidylcholine membranes of different thickness. We make use of both charged (methylammonium and methylguanidinium) and neutral (methylamine and methylguanidine) analogue molecules, as well as Lys and Arg side chains on transmembrane helix models. We find that the free energy barrier experienced by a charged Lys crossing the membrane is strikingly similar to that of a charged Arg (to within 2 kcal/mol), despite the two having different chemistries, H-bonding capability, and hydration free energies that differ by ∼10 kcal/mol. In comparison, the barrier for neutral Arg is higher than that for neutral Lys by around 5 kcal/mol, being more selective than that for the charged species. This can be explained by the different transport mechanisms for charged or neutral amino acid side chains in the membrane, involving membrane deformations or simple dehydration, respectively. As a consequence, we demonstrate that Lys would be deprotonated in the membrane, whereas Arg would maintain its charge. Our simulations also reveal that Arg attracts more phosphate and water in the membrane, and can form extensive H-bonding with its five H-bond donors to stabilize Arg-phosphate clusters. This leads to enhanced interfacial binding and membrane perturbations, including the appearance of a trans-membrane pore in a thinner membrane. These results highlight the special role played by Arg as an amino acid to bind to, disrupt, and permeabilize lipid membranes, as well as to sense voltages for a range of peptide and protein activities in nature and in engineered bionanodevices.

Conformational studies of sphingolipids by NMR spectroscopy. I. Dihydrosphingomyelin

The conformational features of dihydrosphingomyelin (DHSM), the major phospholipid of human lens membranes, were investigated by 1H and 31P nuclear magnetic resonance spectroscopy. Several postulates emerge from the observed trends: (a) in partially hydrated samples of DHSM in CDCl3 above 13 mM, at which lipid-lipid interactions prevail, the amide proton is mostly involved in intermolecular H-bonds that link neighboring phospholipids through bridging water molecules. In the absence of water, the NH group is involved in an intramolecular H-bond that restricts the mobility of the phosphate group. (b) In the monomeric form of the lipid molecule, the amide proton of the major conformer is bound intramolecularly with one of the anionic and/or ester oxygens of the phosphate group. A minor conformer may also be present in which the NH proton participates in an intramolecular H-bond linking to the OH group of the sphingoid base. (c) Complete hydration leads to an extension of the head group as water molecules bind to the phosphate and NH groups via H-bonds, thus disrupting the intramolecular H-bonds prevalent at low concentrations.

Glycero- versus sphingo-phospholipids: correlations with human and non-human mammalian lens growth

The human lens differs from other mammalian lenses in its very slow growth and unusual phospholipid composition of its cell membranes. Dihydrosphingomyelins (DHSMs) make up about half of all phospholipids in adult human fiber membranes. In all other membranes, sphingomyelins(SMs) with a trans double bond in their backbone, are prevalent. In our quest to understand the biological implications of such elevated DHSM levels, we analyzed membranes from various regions of human, elephant, giraffe, polar bear, pig and cow lenses. The levels of DHSMs were minor in non-human lens membranes. A strong correlation was observed between growth rate and relative contents of phosphatidylcholines(PCs) in epithelia and outer cortical fibers. Sphingomyelins became increasingly predominant in differentiated fibers and this increase was age dependent. Indeed, nuclear fiber membranes of aged non-human mammals were composed, almost exclusively, of (SMs). Although human lens membranes followed comparable compositional trends, the magnitude of the changes was much smaller. We postulate that the high relative contents of DHSMs provide a biochemically inert matrix in which only small amounts of PCs and SMs and their metabolites, known to promote and arrest growth, respectively, are present. This compositional difference is proposed to contribute to the slow multiplication and elongation of human lens cells.

Ion conduction and conformational flexibility of a bacterial voltage-gated sodium channel

Significance Voltage-gated sodium channels are one of the most fundamental electrical components in the nervous system and are key targets for local anesthesia and therapeutics for neurological and cardiac disorders. We have used multimicrosecond simulations to provide molecular-level descriptions of sodium channel function. We describe an almost barrier-less three-ion conduction mechanism involving competing knock-on and “pass-by” processes, intimately linked to signature glutamate ring protonation and structural isomerizations. These simulations have uncovered a high degree of protein flexibility, with conformational fluctuations in the pore domain involving residues central to slow-type inactivation, leading to gate collapse, helix bending, filter disruption, and changes in lipid-facing fenestrations linked to Nav drug pathways. Voltage-gated Na+ channels play an essential role in electrical signaling in the nervous system and are key pharmacological targets for a range of disorders. The recent solution of X-ray structures for the bacterial channel NavAb has provided an opportunity to study functional mechanisms at the atomic level. This channel’s selectivity filter exhibits an EEEE ring sequence, characteristic of mammalian Ca2+, not Na+, channels. This raises the fundamentally important question: just what makes a Na+ channel conduct Na+ ions? Here we explore ion permeation on multimicrosecond timescales using the purpose-built Anton supercomputer. We isolate the likely protonation states of the EEEE ring and observe a striking flexibility of the filter that demonstrates the necessity for extended simulations to study conduction in this channel. We construct free energy maps to reveal complex multi-ion conduction via knock-on and “pass-by” mechanisms, involving concerted ion and glutamate side chain movements. Simulations in mixed ionic solutions reveal relative energetics for Na+, K+, and Ca2+ within the pore that are consistent with the modest selectivity seen experimentally. We have observed conformational changes in the pore domain leading to asymmetrical collapses of the activation gate, similar to proposed inactivated structures of NavAb, with helix bending involving conserved residues that are critical for slow inactivation. These structural changes are shown to regulate access to fenestrations suggested to be pathways for lipophilic drugs and provide deeper insight into the molecular mechanisms connecting drug activity and slow inactivation.