Electrostatics calculations: recent methodological advances and applications to membranes

In vitro and in silico perspectives on the activation of antioxidant responsive element by citrus-derived flavonoids

Introduction

Oxidative stress plays an essential role in the pathogenesis of chronic diseases. Disrupting the Keap1-Nrf2 pathway by binding Keap1 is identified as a potential strategy to prevent oxidative stress-related chronic diseases. Therefore, of special interest is the utilization of dietary antioxidations from citrus, including narirutin, naringenin, hesperetin, hesperidin, naringin, neohesperidin dihydrochalcone, neohesperidin, and nobiletin, has been exploited as a prospective way to treat or prevent several human pathologies as Keap1-Nrf2 inhibitors for modulation of antioxidant properties.

Methods

To probe into the structural foundation of the molecular identification of citrus-derived antioxidations, we calculated the antioxidant responsive element activation ability of citrus-derived flavonoids after binding with Keap1. Also, the quantum chemistry properties and binding mode were performed theoretically with frontier molecular orbitals, molecular electrostatic potential analysis, molecular docking, and absorption, distribution, metabolism, excretion (ADME) calculation.

Results and discussion

Experimental findings combining computational assays revealed that the tested citrus-derived flavonoids can be grouped into strong agonists and weak agonists. The citrus-derived antioxidations were well housed in the bound zone of Keap1 via stable hydrogen bonding and hydrophobic interaction. Eventually, three of eight antioxidations were identified after ADME and physicochemical evaluations. The citrus-derived flavonoids were identified as potential dietary antioxidants of the Keap1-Nrf2 interaction, and can be used to improve oxidative stress-related chronic diseases.

Heterogeneous dielectric generalized Born model with a van der Waals term provides improved association energetics of membrane-embedded transmembrane helices

The heterogeneous dielectric generalized Born (HDGB) implicit membrane formalism is extended by the addition of a van der Waals dispersion term to better describe the nonpolar components of the free energy of solvation. The new model, termed HDGBvdW, improves the energy estimates in the hydrophobic interior of the membrane, where polar and charged species are rarely found and nonpolar interactions become significant. The implicit van der Waals term for the membrane environment extends the model from Gallicchio et al. (J. Comput. Chem. 2004, 25, 479) by combining separate contributions from each of the membrane components. The HDGBvdW model is validated with a series of test cases ranging from membrane insertion and pair association profiles of amino acid side chain analogs and transmembrane helices. Overall, the HDGBvdW model leads to increased agreement with explicit membrane simulation results and experimental data. © 2016 Wiley Periodicals, Inc.

Efficient Multistate Reactive Molecular Dynamics Approach Based on Short-Range Effective Potentials

Nonbonded interactions between molecules usually include the van der Waals force and computationally expensive long-range electrostatic interactions. This article develops a more efficient approach: the effective-interaction multistate empirical-valence-bond (EI-MS-EVB) model. The EI-MS-EVB method relies on a mapping of all interactions onto a short-range and thus, computationally efficient effective potential. The effective potential is tabulated by matching its force to known trajectories obtained from the full-potential empirical multistate empirical-valence-bond (MS-EVB) model. The effective pairwise interaction depends on and is uniquely determined by the atomic configuration of the system, varying only with respect to the hydrogen-bonding topology. By comparing the EI-MS-EVB and full MS-EVB calculations of several equilibrium and dynamic properties important to hydrated excess proton solvation and transport, we show that the EI-MS-EVB model produces very accurate results for the specific system in which the tabulated potentials were generated. The EI-MS-EVB potential also transfers reasonably well to similar systems with different temperatures and box sizes. The EI-MS-EVB method also reduces the computational cost of the nonbonded interactions by about 1 order of magnitude in comparison with the full algorithm.

Isotropic Periodic Sum Treatment of Long-Range Electrostatic Interactions in Combined Quantum Mechanical and Molecular Mechanical Calculations

The isotropic periodic sum (IPS) method was extended to describe long-range electrostatic interactions in combined quantum mechanical and molecular mechanical (QM/MM) calculations. The resulting method, designated QM/MM-IPS, was tested for two ion association processes and a model SN2 reaction in aqueous solution. Potential of mean force (PMF) profiles and radial distribution functions computed from the QM/MM-IPS simulations were compared with those obtained by using the existing QM/MM-Ewald sum and cutoff (QM/MM-Cutoff) methods. In contrast to the QM/MM-Cutoff method, with which PMFs of ion separation tend to display a spurious linear drift, the QM/MM-IPS method successfully eliminates such artifacts, in excellent agreement with the QM/MM-Ewald results. The PMF obtained with the QM/MM-IPS method for the SN2 reaction that transfers an NH3 group between two chloride anions closely resembles that from the QM/MM-Ewald simulations. Compared with QM/MM-Ewald, the QM/MM-IPS method reduces the computational cost by 60-70% when a local region of 12 to 14 Å is used. These results suggest that the QM/MM-IPS method can be used as a reliable and efficient alternative to the QM/MM-Ewald method to incorporate long-range electrostatic effects in simulating solution-phase chemical reactions.

Simple and accurate scheme to compute electrostatic interaction: zero-dipole summation technique for molecular system and application to bulk water

The zero-dipole summation method was extended to general molecular systems, and then applied to molecular dynamics simulations of an isotropic water system. In our previous paper [I. Fukuda, Y. Yonezawa, and H. Nakamura, J. Chem. Phys. 134, 164107 (2011)], for evaluating the electrostatic energy of a classical particle system, we proposed the zero-dipole summation method, which conceptually prevents the nonzero-charge and nonzero-dipole states artificially generated by a simple cutoff truncation. Here, we consider the application of this scheme to molecular systems, as well as some fundamental aspects of general cutoff truncation protocols. Introducing an idea to harmonize the bonding interactions and the electrostatic interactions in the scheme, we develop a specific algorithm. As in the previous study, the resulting energy formula is represented by a simple pairwise function sum, enabling facile applications to high-performance computation. The accuracy of the electrostatic energies calculated by the zero-dipole summation method with the atom-based cutoff was numerically investigated, by comparison with those generated by the Ewald method. We obtained an electrostatic energy error of less than 0.01% at a cutoff length longer than 13 Å for a TIP3P isotropic water system, and the errors were quite small, as compared to those obtained by conventional truncation methods. The static property and the stability in an MD simulation were also satisfactory. In addition, the dielectric constants and the distance-dependent Kirkwood factors were measured, and their coincidences with those calculated by the particle mesh Ewald method were confirmed, although such coincidences are not easily attained by truncation methods. We found that the zero damping-factor gave the best results in a practical cutoff distance region. In fact, in contrast to the zero-charge scheme, the damping effect was insensitive in the zero-charge and zero-dipole scheme, in the molecular system we treated. We discussed the origin of this difference between the two schemes and the dependence of this fact on the physical system. The use of the zero damping-factor will enhance the efficiency of practical computations, since the complementary error function is not employed. In addition, utilizing the zero damping-factor provides freedom from the parameter choice, which is not trivial in the zero-charge scheme, and eliminates the error function term, which corresponds to the time-consuming Fourier part under the periodic boundary conditions.

Assessment of Common Simulation Protocols for Simulations of Nanopores, Membrane Proteins, and Channels

Molecular dynamics (MD) simulation has become a common technique to study biological systems. Transport of small molecules through carbon nanotubes and membrane proteins has been an intensely studied topic, and MD simulations have been able to provide valuable predictions, many of which have later been experimentally proven. Simulations of such systems pose challenges, and unexpected problems in commonly used protocols and methods have been found in the past few years. The two main reasons why some were not found before are that most of these newly discovered errors do not lead to unstable simulations. Furthermore, some of them manifest themselves only after relatively long simulation times. We assessed the reliability of the most common simulations protocols by MD and stochastic dynamics (SD) or Langevin dynamics, simulations of an alpha hemolysin nanochannel embedded in a palmitoyloleoylphosphatidylcholine (POPC) lipid bilayer. Our findings are that (a) reaction field electrostatics should not be used in simulations of such systems, (b) local thermostats should be preferred over global ones since the latter may lead to an unphysical temperature distribution,

The ELBA force field for coarse-grain modeling of lipid membranes

A new coarse-grain model for molecular dynamics simulation of lipid membranes is presented. Following a simple and conventional approach, lipid molecules are modeled by spherical sites, each representing a group of several atoms. In contrast to common coarse-grain methods, two original (interdependent) features are here adopted. First, the main electrostatics are modeled explicitly by charges and dipoles, which interact realistically through a relative dielectric constant of unity (ε(r) = 1). Second, water molecules are represented individually through a new parametrization of the simple Stockmayer potential for polar fluids; each water molecule is therefore described by a single spherical site embedded with a point dipole. The force field is shown to accurately reproduce the main physical properties of single-species phospholipid bilayers comprising dioleoylphosphatidylcholine (DOPC) and dioleoylphosphatidylethanolamine (DOPE) in the liquid crystal phase, as well as distearoylphosphatidylcholine (DSPC) in the liquid crystal and gel phases. Insights are presented into fundamental properties and phenomena that can be difficult or impossible to study with alternative computational or experimental methods. For example, we investigate the internal pressure distribution, dipole potential, lipid diffusion, and spontaneous self-assembly. Simulations lasting up to 1.5 microseconds were conducted for systems of different sizes (128, 512 and 1058 lipids); this also allowed us to identify size-dependent artifacts that are expected to affect membrane simulations in general. Future extensions and applications are discussed, particularly in relation to the methodology's inherent multiscale capabilities.

Molecular dynamics scheme for precise estimation of electrostatic interaction via zero-dipole summation principle

We propose a novel idea, zero-dipole summation, for evaluating the electrostatic energy of a classical particle system, and have composed an algorithm for effectively utilizing the idea for molecular dynamics. It conceptually prevents the nonzero-charge and nonzero-dipole states artificially generated by a simple cutoff truncation. The resulting energy formula is nevertheless represented by a simple pairwise function sum, which enables facile application to high-performance computation. By following a heuristic approach to derive the current electrostatic energy formula, we developed an axiomatic approach to construct the method consistently. Explorations of the theoretical details of our method revealed the structure of the generated error, and we analyzed it by comparisons with other methods. A numerical simulation using liquid sodium chloride confirmed that the current method with a small damping factor yielded sufficient accuracy with a practical cutoff distance region. The current energy function also conducts stable numerical integration in a liquid MD simulation. Our method is an extension of the charge neutralized summation developed by Wolf et al. [J. Chem. Phys. 110, 8254 (1999)]. Furthermore, we found that the current method becomes a generalization of the preaveraged potential method proposed by Yakub and Ronchi [J. Chem. Phys. 119, 11556 (2003)], which is based on a viewpoint different from the neutrality. The current study presents these relationships and suggests possibilities for their further applications.

Reactivity of a phospholipid monolayer model under periodic boundary conditions: a density functional theory study of the Schiff base formation between phosphatidylethanolamine and acetaldehyde

A mechanism for the formation of the Schiff base between an acetaldehyde and an amine-phospholipid monolayer model based on Dmol3/density functional theory calculations under periodic boundary conditions was constructed. This is the first time such a system has been modeled to examine its chemical reactivity at this computation level. Each unit cell contains two phospholipid molecules, one acetaldehyde molecule, and nine water molecules. One of the amine-phospholipid molecules in the cell possesses a neutral amino group that is used to model the nucleophilic attack on the carboxyl group of acetaldehyde, whereas the other has a charged amino group acting as a proton donor. The nine water molecules form a hydrogen bond network along the polar heads of the phospholipids that facilitates very fast proton conduction at the interface. Using periodic boundary conditions afforded proton transfer between different cells. The reaction takes place in two steps, namely, (1) formation of a carbinolamine and (2) its dehydration to the Schiff base. The carbinolamine is the primary reaction intermediate, and dehydration is the rate-determining step of the process, consistent with available experimental evidence for similar reactions. On the basis of the results, the cell membrane surface environment may boost phospholipid glycation via a neighboring catalyst effect.

Coarse-grained ions without charges: reproducing the solvation structure of NaCl in water using short-ranged potentials

A coarse-grained model of NaCl in water is presented where the ions are modeled without charge to avoid computationally challenging electrostatics. A monatomic model of water [V. Molinero and E. B. Moore, J. Phys. Chem. B 113, 4008 (2009)] is used as the basis for this coarse-grain approach. The ability of Na(+) to disrupt the native tetrahedral arrangement of water molecules, and of Cl(-) to integrate within this organization, is preserved in this mW-ion model through parametrization focused on water's solvation of these ions. This model successfully reproduces the structural effect of ions on water, referenced to observations from experiments and atomistic molecular dynamics simulations, while using extremely short-ranged potentials. Without Coulomb interactions the model replicates details of the ion-water structure such as distinguishing contact and solvent-separated ion pairs and the free energy barriers between them. The approach of mimicking ionic effects with short-ranged interactions results in performance gains of two orders of magnitude compared to Ewald methods. Explored over a broad range of salt concentration, the model reproduces the solvation structure and trends of diffusion relative to atomistic simulations and experimental results. The functional form of the mW-ion model can be parametrized to represent other electrolytes. With increased computational efficiency and reliable structural fidelity, this model promises to be an asset for accessing significantly longer simulation time scales with an explicit solvent in a coarse-grained system involving, for example, polyelectrolytes such as proteins, nucleic acids, and fuel-cell membranes.

Using DelPhi to compute electrostatic potentials and assess their contribution to interactions

There is a general agreement that electrostatic interactions play a significant role in the structure and function of biological molecules. However, obtaining quantitative estimation of the electrostatic energy requires computational models that capture the microscopic nature of the heterogeneous environment of macromolecules. This protocol offers elaboration on one of the common methods to calculate the electrostatic energetic contributions using continuum electrostatics. The method involves solving the Poisson-Boltzmann (PB) equation numerically and regarding the solute as having a homogenous dielectric constant. In order to apply this method, a three dimensional structure of the molecule derived from experimental data (crystallography, NMR) or modeling techniques is required. The protocol will focus on the DelPhi program (Accelrys Inc. San Diego), which is one of the most common programs used for the estimation of electrostatic free energy contribution. A simple procedure of assigning criteria and parameters (charge distribution, solvent and solute dielectric constants, iterations, grid resolution, etc) enables one to illustrate an electrostatic potential map and estimate the electrostatic free energy, although with limited accuracy.

A quantitative coarse-grain model for lipid bilayers

A simplified particle-based computer model for hydrated phospholipid bilayers has been developed and applied to quantitatively predict the major physical features of fluid-phase biomembranes. Compared with available coarse-grain methods, three novel aspects are introduced. First, the main electrostatic features of the system are incorporated explicitly via charges and dipoles. Second, water is accurately (yet efficiently) described, on an individual level, by the soft sticky dipole model. Third, hydrocarbon tails are modeled using the anisotropic Gay-Berne potential. Simulations are conducted by rigid-body molecular dynamics. Our technique proves 2 orders of magnitude less demanding of computational resources than traditional atomic-level methodology. Self-assembled bilayers quantitatively reproduce experimental observables such as electron density, compressibility moduli, dipole potential, lipid diffusion, and water permeability. The lateral pressure profile has been calculated, along with the elastic curvature constants of the Helfrich expression for the membrane bending energy; results are consistent with experimental estimates and atomic-level simulation data. Several of the results presented have been obtained for the first time using a coarse-grain method. Our model is also directly compatible with atomic-level force fields, allowing mixed systems to be simulated in a multiscale fashion.

Is the Ewald summation still necessary? Pairwise alternatives to the accepted standard for long-range electrostatics

We investigate pairwise electrostatic interaction methods and show that there are viable computationally efficient (O(N)) alternatives to the Ewald summation for typical modern molecular simulations. These methods are extended from the damped and cutoff-neutralized Coulombic sum originally proposed by Wolf et al. [J. Chem. Phys. 110, 8255 (1999)]. One of these, the damped shifted force method, shows a remarkable ability to reproduce the energetic and dynamic characteristics exhibited by simulations employing lattice summation techniques. Comparisons were performed with this and other pairwise methods against the smooth particle-mesh Ewald summation to see how well they reproduce the energetics and dynamics of a variety of molecular simulations.

Acyl chain order parameter profiles in phospholipid bilayers: computation from molecular dynamics simulations and comparison with 2H NMR experiments

Order parameters from deuterium NMR are often used to validate or calibrate molecular dynamics simulations. This paper gives a short overview of the literature in which experimental order parameters from (2)H NMR are compared to those calculated from MD simulations. The different ways in which order parameters from experiment are used to calibrate and validate simulations are reviewed. In the second part of this review, a case study of cholesterol in a DMPC bilayer is presented. It is concluded that the agreement between experimental data and simulation is favorable in the hydrophobic region of the membrane, for both the phospholipids and cholesterol. In the interfacial region the agreement is less satisfactory, probably because of the high polarity of this region which makes the correct computation of the electrostatics more complex.

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

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

Three hydrolases and a transferase: comparative analysis of active-site dynamics via the BioSimGrid database

Comparative molecular dynamics (MD) simulations enable us to explore the conformational dynamics of the active sites of distantly related enzymes. We have used the BioSimGrid (http://www.biosimgrid.org) database to facilitate such a comparison. Simulations of four enzymes were analyzed. These included three hydrolases and a transferase, namely acetylcholinesterase, outer-membrane phospholipase A, outer-membrane protease T, and PagP (an outer-membrane enzyme which transfers a palmitate chain from a phospholipid to lipid A). A set of 17 simulations were analyzed corresponding to a total of approximately 0.1 micros simulation time. A simple metric for active-site integrity was used to demonstrate the existence of clusters of dynamic conformational behaviour of the active sites. Small (i.e. within a cluster) fluctuations appear to be related to the function of an enzymatically active site. Larger fluctuations (i.e. between clusters) correlate with transitions between catalytically active and inactive states. Overall, these results demonstrate the potential of a comparative MD approach to analysis of enzyme function. This approach could be extended to a wider range of enzymes using current high throughput MD simulation and database methods.

Electrostatics calculations: latest methodological advances

Electrostatics plays a major role in the stabilization and function of biomolecules; as such, it remains a major focus of theoretical and computational studies of macromolecules. Electrostatic interactions are long range, and strongly dependent on the solvent and ions surrounding the biomolecule under study. During the past year, progress has been reported in the treatment of electrostatics in explicit and implicit solvent models. Interesting new developments of explicit solvent models include more efficient Ewald summation methods, as well as alternative approaches based on reaction field theory, periodic images and Euler summations. Implicit solvent models remain divided into those that solve the Poisson-Boltzmann equation numerically and those based on the generalized Born formalism. Both approaches are now included in molecular dynamics simulations and their accuracies may be assessed by direct comparison against experimental data. It is worth mentioning the recent development of web interfaces that facilitate access to and usage of existing tools for computing electrostatic interactions.

Molecular simulation study of structural and dynamic properties of mixed DPPC/DPPE bilayers

Molecular dynamics simulations have been used to study structural and dynamic properties of fully hydrated mixed 1,2-dipalmitoyl-sn-glycero-3-phosphocholine (DPPC) and 1,2-dipalmitoyl-sn-glycero-3-phosphoethanolamine (DPPE) bilayers at 0, 25, 50, 75, and 100 mol % DPPE. Simulations were performed for 50 ns at 350 K and 1 bar for the liquid-crystalline state of the mixtures. Results show that the average area per headgroup reduces from 0.65 +/- 0.01 nm(2) in pure DPPC to 0.52 +/- 0.01 nm(2) in pure DPPE systems. The lipid tails become more ordered with increasing DPPE concentration, resulting in a slight increase in membrane thickness (3.43 +/- 0.01 nm in pure DPPC to 4.00 +/- 0.01 nm in pure DPPE). The calculated area per headgroup and order parameter for pure DPPE deviates significantly from available experimental measurements, suggesting that the force field employed requires further refinement. In-depth analysis of the hydrogen-bond distribution in DPPE molecules shows that the amine groups strongly interact with the phosphate and carbonyl groups through inter/intramolecular hydrogen bonds. This yields a bilayer structure with DPPE headgroups preferentially located near the lipid phosphate and ester oxygens. It is observed that increasing DPPE concentrations causes competitive hydrogen bonding between the amine groups (hydrogen-donor) and the phosphate/carbonyl groups or water (hydrogen-acceptor). Due to the increasing number of hydrogen-donors from DPPE molecules with increasing concentration, DPPE becomes more hydrated. Trajectory analysis shows that DPPE molecules in the lipid mixtures move laterally and randomly around the membrane surface and the movement becomes more localized with increasing DPPE concentrations. For the conditions and simulation time considered, no aggregation or phase separation was observed between DPPC and DPPE.

Interface connections of a transmembrane voltage sensor

Voltage-sensitive ion channels open and close in response to changes in transmembrane (TM) potential caused by the motion of the S4 voltage sensors. These sensors are alpha-helices that include four or more positively charged amino acids, most commonly arginine. The so-called paddle model, based on the high-resolution structure of the KvAP K+ channel [Jiang, et al. (2003) Nature 423, 33-41], posits that the S4 sensors move within the membrane bilayer in response to TM voltage changes. Direct exposure of S4 sensors to lipid is contrary to the classical expectation that the dielectric contrast between the membrane hydrocarbon core and water presents an insurmountable energetic penalty to burial of electric charges. Nevertheless, recent experiments have shown that a helix with the sequence of KvAP S4 can be inserted across the endoplasmic reticulum membrane. To reconcile this result with the classical energetics argument, we have carried out a molecular dynamics simulation of an isolated TM S4 helix in a lipid bilayer. The simulation reveals a stabilizing hydrogen-bonded network of water and lipid phosphates around the arginines that reduces the effective thickness of the bilayer hydrocarbon core to approximately 10 A in the vicinity of the helix. It suggests that bilayer phospholipids can adapt locally to strongly perturbing protein elements, causing the phospholipids to become a structural extension of the protein.

Molecular dynamics simulation approaches to K channels: conformational flexibility and physiological function

Molecular modeling and simulations enable extrapolation for the structure of bacterial potassium channels to the function of their mammalian homologues. Molecular dynamics simulations have revealed the concerted single-file motion of potassium ions and water molecules through the selectivity filter of K channels and the role of filter flexibility in ion permeation and in "fast gating." Principal components analysis of extended K channel simulations suggests that hinge-bending of pore-lining M2 (or S6) helices plays a key role in K channel gating. Based on these and other simulations, a molecular model for gating of inward rectifier K channel gating is presented.

Asymmetry of lipid bilayers induced by monovalent salt: Atomistic molecular-dynamics study

Interactions between salt ions and lipid components of biological membranes are essential for the structure, stability, and functions of the membranes. The specific ionic composition of aqueous buffers inside and outside of the cell is known to differ considerably. To model such a situation we perform atomistic molecular-dynamics (MD) simulations of a single-component phosphatidylcholine lipid bilayer which separates two aqueous reservoirs with and without NaCl salt. To implement the difference in electrolyte composition near two membrane sides, a double bilayer setup (i.e., two bilayers in a simulation box) is employed. It turns out that monovalent salt, being in contact with one leaflet only, induces a pronounced asymmetry in the structural, electrostatic, and dynamical properties of bilayer leaflets after 50 ns of MD simulations. Binding of sodium ions to the carbonyl region of the leaflet which is in contact with salt results in the formation of "Na-lipids" complexes and, correspondingly, reduces mobility of lipids of this leaflet. In turn, attractive interactions of chloride ions (mainly located in the aqueous phase close to the water-lipid interface) with choline lipid groups lead to a substantial (more vertical) reorientation of postphatidylcholine headgroups of the leaflet adjoined to salt. The difference in headgroup orientation on two sides of a bilayer, being coupled with salt-induced reorientation of water dipoles, leads to a notable asymmetry in the charge-density profiles and electrostatic potentials of bilayer constitutes of the two leaflets. Although the overall charge density of the bilayer is found to be almost insensitive to the presence of salt, a slight asymmetry in the charge distribution between the two bilayer leaflets results in a nonzero potential difference of about 85 mV between the two water phases. Thus, a transmembrane potential of the order of the membrane potential in a cell can arise without ionic charge imbalance between two aqueous compartments.

Molecular Dynamics Simulations of Aqueous Pullulan Oligomers

Molecular dynamics (MD) simulations have been used to model small-angle X-ray scattering (SAXS) data on aqueous solutions of four oligomeric segments of the glucan pullulan: the trimer G(3) (comprising one polymer repeating unit), the hexamer (G(3))(2), the nonamer (G(3))(3), and the dodecamer (G(3))(4). The AMBER force field was used in conjunction with the GB/SA continuum solvation model to calculate both the mean global dimensions of the oligomers from the limiting small angle scattering behavior and the shorter range structural information implicit in the Debye scattering function at larger scattering angles. This same force field and solvation treatment were employed earlier by Liu et al. (Macromolecules 1999, 32, 8611-8620) with apparent success in a rotational isomeric state (RIS) treatment of the same experimental data. The present work discloses that, despite numerical success in modeling the SAXS data, the RIS treatment, which includes only the interactions within dimeric segments of the polymer chain, fails to account accurately for excluded volume effects at the range of 3-12 sugar residues in the polymer backbone. It is suggested that MD simulations using continuum solvation models can be used to circumvent errors inherent in the computationally efficient RIS treatments of polymer nano- and picosecond dynamics while at the same time avoiding the heavy computational requirements of all-atom methods.

Physical landscapes in biological membranes: physico-chemical terrains for spatio-temporal control of biomolecular interactions and behaviour

The evolving complexities of biological membranes are discussed from the point of view of potential roles of the physical constitution of the membrane. These include features of the surface and dipole potentials and membrane 'rafts'. These properties are outlined; they emphasize that protein-lipid and specific lipid environments are influential parameters in how biomolecular interactions may take place with and within membranes. Several fluorescence detection technologies directed towards measurement of these properties are also outlined that permit high-resolution experimental determination of intermolecular interactions with membranes by measuring small changes of these potentials. These point to the possibility that the membrane dipole potential in particular is enormously influential in determining the behaviour of receptor and signalling systems within membrane rafts, and offers the means of a novel mechanism for biological control.

New Method for Direct Linear-Scaling Calculation of Electron Density of Proteins

A new scheme for direct linear-scaling quantum mechanical calculation of electron density of protein systems is developed. The new scheme gives much improved accuracy of electron density for proteins than the original MFCC (molecular fractionation with conjugate caps) approach in efficient linear-scaling calculation for protein systems. In this new approach, the error associated with each cut in the MFCC approach is estimated by computing the two neighboring amino acids in both cut and uncut calculations and is corrected. Numerical tests are performed on six oligopeptide taken from PDB (protein data bank), and the results show that the new scheme is efficient and accurate.

The transmembrane domain of the acetylcholine receptor: insights from simulations on synthetic peptide models

We have studied the structure and properties of a bundle of alpha-helical peptides embedded in a 1,2-dimyristoyl-3-phosphatidylcholine phospholipid bilayer by molecular dynamics simulations. The bundle of five transmembrane deltaM2 segments constitutes the model for the pore region of the nicotinic acetylcholine receptor, which is the neurotransmitter-gated ion-channel responsible for the fast propagation of electrical signals between cells at the nerve-muscle synapse. The deltaM2 segments were shown to oligomerize in biomembranes resulting in ion-channel activity with characteristics similar to the native protein, and the structure of the isolated peptides was studied in 1,2-dimyristoyl-3-phosphatidylcholine bilayers and micelles by NMR experiments (Opella, S. J., et al. 1999. Nat. Struct. Biol. 6:374-379). Our analyses indicate that the structure, helix tilt, and the overall shape of the channel are in good agreement with the NMR experiments and the proposed model for the channel, which we show is formed by rings of functional residues. The studied geometry resulted in a closed pore state, where the channel is partially dehydrated at the hydrophobic extracellular half and the extracellular mouth of the channel blocked by the hydrocarbon chains of Arg+ residues. The arginine amino acids form intermolecular salt-bridges with the C-terminus, which contribute as well to the bundle stabilization.

Lipid-protein interactions of integral membrane proteins: a comparative simulation study

The interactions between membrane proteins and their lipid bilayer environment play important roles in the stability and function of such proteins. Extended (15-20 ns) molecular dynamics simulations have been used to explore the interactions of two membrane proteins with phosphatidylcholine bilayers. One protein (KcsA) is an alpha-helix bundle and embedded in a palmitoyl oleoyl phosphatidylcholine bilayer; the other (OmpA) is a beta-barrel outer-membrane protein and is in a dimyristoyl phosphatidylcholine bilayer. The simulations enable analysis in detail of a number of aspects of lipid-protein interactions. In particular, the interactions of aromatic amphipathic side chains (i.e., Trp, Tyr) with lipid headgroups, and "snorkeling" interactions of basic side chains (i.e., Lys, Arg) with phosphate groups are explored. Analysis of the number of contacts and of H-bonds reveal fluctuations on an approximately 1- to 5-ns timescale. There are two clear bands of interacting residues on the surface of KcsA, whereas there are three such bands on OmpA. A large number of Arg-phosphate interactions are seen for KcsA; for OmpA, the number of basic-phosphate interactions is smaller and shows more marked fluctuations with respect to time. Both classes of interaction occur in clearly defined interfacial regions of width approximately 1 nm. Analysis of lateral diffusion of lipid molecules reveals that "boundary" lipid molecules diffuse at about half the rate of bulk lipid. Overall, these simulations present a dynamic picture of lipid-protein interactions: there are a number of more specific interactions but even these fluctuate on an approximately 1- to 5-ns timescale.

The range and shielding of dipole-dipole interactions in phospholipid bilayers

In molecular dynamics simulations of lipid bilayers, the structure is sensitive to the precise treatment of electrostatics. The dipole-dipole interactions between headgroup dipoles are not long-ranged, but the area per lipid and, through it, other properties of the bilayer are very sensitive to the detailed balance between the perpendicular and in-plane components of the headgroup dipoles. This is affected by the detailed properties of the cutoff scheme or if long-range interactions are included by Ewald or particle-mesh Ewald techniques. Interaction between the in-plane components of the headgroup dipoles is attractive and decays as the inverse sixth power of distance. The interaction is screened by the square of a dielectric permittivity close to the value for water. Interaction between the components perpendicular to the membrane plane is repulsive and decays as the inverse third power of distance. These interactions are screened by a dielectric permittivity of the order 10. Thus, despite the perpendicular components being much smaller in magnitude than the in-plane components, they will dominate the interaction energies at large distances.

Molecular simulation study on the influence of dimethylsulfoxide on the structure of phospholipid bilayers

Molecular dynamics simulations of dipalmitoylphosphatidylcholine (DPPC) lipid bilayer/water systems were performed in the presence of dimethylsulfoxide (DMSO) at 2, 5, 10, and 100 mol % DMSO (lipid-free basis). The equilibrium structure and several dynamic properties were determined for these systems. Results show that DMSO penetrates much deeper into the bilayer than water does. It is also found that DMSO molecules do not interact with the polar groups of the lipid headgroup, but exhibit a preference to remain either directly below the headgroup or in the aqueous phase, which is a consequence of the chemical characteristics of DMSO. As the temperature increases, a higher DMSO concentration is observed in the bilayer side of the interface. The area per headgroup in the presence of DMSO is significantly increased from 66.8 A2 for the pure bilayer to as high as 87.0 A2 at 10 mol % DMSO at 350 K. DMSO hydrogen-bonds strongly with water and exhibits unfavorable interactions with the polar headgroups of DPPC, thereby inducing a dehydration of the headgroups.

Molecular simulation study of phospholipid bilayers and insights of the interactions with disaccharides

Molecular simulations of hydrated dipalmitoylphosphatidylcholine lipid bilayers have been performed for temperatures in the range of 250-450 K. The area per headgroup increases with temperature from 58 to 77 A(2). Other properties such as hydration number, alkyl tail order parameter, diffusion coefficients, and radial distribution functions exhibit a clear dependence on temperature. Simulations of bilayers have also been performed in the presence of two disaccharides, namely trehalose and sucrose, at concentrations of up to 18 wt % (lipid-free basis). The simulated area per headgroup of the bilayer is not affected by the presence of the disaccharides, suggesting that the overall structure of the bilayer remains undisturbed. The results of simulations reveal that the interaction of disaccharide molecules with the bilayer occurs at the surface of the bilayer, and it is governed by the formation of multiple hydrogen bonds to specific groups of the lipid. Disaccharide molecules are observed to adopt specific conformations to fit onto the surface topology of the bilayer, often interacting with up to three different lipids simultaneously. At high disaccharide concentrations, the results of simulations indicate that disaccharides can serve as an effective replacement for water under anhydrous conditions, which helps explain their effectiveness as lyophilization agents for liposomes and cells.

Cationic DMPC/DMTAP lipid bilayers: molecular dynamics study

Cationic lipid membranes are known to form compact complexes with DNA and to be effective as gene delivery agents both in vitro and in vivo. Here we employ molecular dynamics simulations for a detailed atomistic study of lipid bilayers consisting of a mixture of cationic dimyristoyltrimethylammonium propane (DMTAP) and zwitterionic dimyristoylphosphatidylcholine (DMPC). Our main objective is to examine how the composition of the DMPC/DMTAP bilayers affects their structural and electrostatic properties in the liquid-crystalline phase. By varying the mole fraction of DMTAP, we have found that the area per lipid has a pronounced nonmonotonic dependence on the DMTAP concentration, with a minimum around the point of equimolar DMPC/DMTAP mixture. We show that this behavior has an electrostatic origin and is driven by the interplay between positively charged TAP headgroups and the zwitterionic phosphatidylcholine (PC) heads. This interplay leads to considerable reorientation of PC headgroups for an increasing DMTAP concentration, and gives rise to major changes in the electrostatic properties of the lipid bilayer, including a significant increase of total dipole potential across the bilayer and prominent changes in the ordering of water in the vicinity of the membrane. Moreover, chloride counterions are bound mostly to PC nitrogens implying stronger screening of PC heads by Cl ions compared to TAP headgroups. The implications of these findings are briefly discussed.

Homology modelling and molecular dynamics simulations: comparative studies of human aquaporin-1

The structures of the mammalian water transport protein Aqp1 and of its bacterial homologue GlpF enables us to test whether homology models can be used to explore relationships between structure, dynamics and function in mammalian transport proteins. Molecular dynamics simulations (totalling almost 40 ns) were performed starting from: the X-ray structure of Aqp1; a homology model of Aqp1 based on the GlpF structure; and intermediate resolution structures of Aqp1 derived from electron microscopy. Comparisons of protein RMSDs vs. time suggest that the homology models are of comparable conformational stability to the X-ray structure, whereas the intermediate resolution structures exhibit significant conformation drift. For simulations based on the X-ray structure and on homology models, the flexibility profile vs. residue number correlates well with the crystallographic B-values for each residue. In the simulations based on intermediate resolution structures, mobility of the highly conserved NPA loops is substantially higher than in the simulations based on the X-ray structure or the homology models. Pore radius profiles remained relatively constant in the X-ray and homology model simulations but showed substantial fluctuations (reflecting the higher NPA loop mobility) in the intermediate resolution simulations. The orientation of the dipoles of water molecules within the pore is of key importance in maintaining low proton permeability through Aqp1. This property seems to be quite robust to the starting model used in the simulation. These simulations suggest that homology models based on bacterial homologues may be used to derive functionally relevant information on the structural dynamics of mammalian transport proteins.

A computational model for the electrostatic sequestration of PI(4,5)P2 by membrane-adsorbed basic peptides

The multivalent acidic phospholipid phosphatidylinositol 4,5-bisphosphate (PI(4,5)P2) plays a key role in many biological processes. Recent studies show that unstructured clusters of basic residues from a number of peripheral proteins can laterally sequester PI(4,5)P2 in membranes. Specifically, experiments suggest that the basic effector domain of the myristoylated alanine-rich C kinase substrate (MARCKS), or a peptide corresponding to this domain, MARCKS(151-175), sequesters several PI(4,5)P2 and that this sequestration is due to nonspecific electrostatic interactions. Here, we use the finite difference Poisson-Boltzmann method to test this hypothesis by calculating the electrostatic free energy of lateral sequestration of PI(4,5)P2 by membrane-adsorbed basic peptides: Lys-7, Lys-13, and FA-MARCKS(151-175), a peptide based on MARCKS(151-175). In agreement with experiments, we find that the electrostatic free energy becomes more favorable when: 1), Lys-13 and FA-MARCKS(151-175) sequester several PI(4,5)P2; 2), the linear charge density of the basic peptide increases; 3), the mol percent monovalent acidic lipid in the membrane decreases; and 4), the ionic strength of the solution decreases. In addition, the electrostatic sequestration free energy is in excess of the entropic penalty associated with localizing PI(4,5)P2. Our calculations, thus, provide a structural and quantitative description of the observed interaction of PI(4,5)P2 with membrane-adsorbed basic sequences.

Rearrangement of charge-charge interactions in variant ubiquitins as detected by double-mutant cycles and NMR

Previous studies of ubiquitin disclosed numerous charge-charge interactions on the protein's surface. To investigate how neighboring residues influence the strength of these interactions, double-mutant cycles are combined with pK(a) determinations by 2D NMR. More specifically, the environment around the Asp21-Lys29 ion pair has been altered through mutations at position 25, which is an asparagine in mammalian ubiquitin and a positively-charged residue in many other ubiquitin-like proteins. The pK(a) value of Asp21 decreases by 0.4 to 0.7 pH unit when Asn25 is substituted with a positively charged residue, suggesting a new and favorable ion pair interaction between positions 21 and 25. However, analysis of double mutants reveals that the favorable interaction between Asp21 and Lys29 is weakened when position 25 is a positively charged residue. Interestingly, while the pK(a) value of His25 in the N25H variant agrees with model compound values, additional mutants reveal that this agreement is fortuitous, resulting from a balance of favorable and unfavorable interactions; similar results were observed previously for Glu34 in ubiquitin and His8 in staphylococcal nuclease. Ionizable groups may thus have pK(a) values similar to model compound values and yet still be involved in significant interactions with other protein groups. One surprising result of introducing positively charged residues at position 25 is a new interaction between Lys29 and Glu18, an interaction not present in wild-type ubiquitin. This unanticipated result illustrates a key advantage of using NMR to determine pK(a) values for many residues simultaneously in the variant proteins. Overall, the strength of an interaction between two residues at the surface of ubiquitin is sensitive to the identity of neighboring residues. The results also demonstrate that relatively conservative and common point mutations such as substitutions of polar with charged residues and vice versa can have effects on interactions beyond the site of mutation per se.