A molecular dynamics study of gating in dioxolane-linked gramicidin A channels

Calcium-gated potassium channel blockade via membrane-facing fenestrations

Quaternary ammonium blockers were previously shown to bind in the pore to block both open and closed conformations of large-conductance calcium-activated potassium (BK and MthK) channels. Because blocker entry was assumed through the intracellular entryway (bundle crossing), closed-pore access suggested that the gate was not at the bundle crossing. Structures of closed MthK, a Methanobacterium thermoautotrophicum homolog of BK channels, revealed a tightly constricted intracellular gate, leading us to investigate the membrane-facing fenestrations as alternative pathways for blocker access directly from the membrane. Atomistic free energy simulations showed that intracellular blockers indeed access the pore through the fenestrations, and a mutant channel with narrower fenestrations displayed no closed-state TPeA block at concentrations that blocked the wild-type channel. Apo BK channels display similar fenestrations, suggesting that blockers may use them as access paths into closed channels. Thus, membrane fenestrations represent a non-canonical pathway for selective targeting of specific channel conformations, opening novel ways to selectively drug BK channels.

Mechanism of C-type inactivation in the hERG potassium channel

The fast C-type inactivation displayed by the voltage-activated potassium channel hERG plays a critical role in the repolarization of cardiac cells, and malfunction caused by nonspecific binding of drugs or naturally occurring missense mutations affecting inactivation can lead to pathologies. Because of its impact on human health, understanding the molecular mechanism of C-type inactivation in hERG represents an advance of paramount importance. Here, long-time scale molecular dynamics simulations, free energy landscape calculations, and electrophysiological experiments are combined to address the structural and functional impacts of several disease-associated mutations. Results suggest that C-type inactivation in hERG is associated with an asymmetrical constricted-like conformation of the selectivity filter, identifying F627 side-chain rotation and the hydrogen bond between Y616 and N629 as key determinants. Comparison of hERG with other K⁺ channels suggests that C-type inactivation depends on the degree of opening of the intracellular gate via the filter-gate allosteric coupling.

A Computational Pipeline to Predict Cardiotoxicity: From the Atom to the Rhythm

RATIONALE

Drug-induced proarrhythmia is so tightly associated with prolongation of the QT interval that QT prolongation is an accepted surrogate marker for arrhythmia. But QT interval is too sensitive a marker and not selective, resulting in many useful drugs eliminated in drug discovery.

OBJECTIVE

To predict the impact of a drug from the drug chemistry on the cardiac rhythm.

METHODS AND RESULTS

In a new linkage, we connected atomistic scale information to protein, cell, and tissue scales by predicting drug-binding affinities and rates from simulation of ion channel and drug structure interactions and then used these values to model drug effects on the hERG channel. Model components were integrated into predictive models at the cell and tissue scales to expose fundamental arrhythmia vulnerability mechanisms and complex interactions underlying emergent behaviors. Human clinical data were used for model framework validation and showed excellent agreement, demonstrating feasibility of a new approach for cardiotoxicity prediction.

CONCLUSIONS

We present a multiscale model framework to predict electrotoxicity in the heart from the atom to the rhythm. Novel mechanistic insights emerged at all scales of the system, from the specific nature of proarrhythmic drug interaction with the hERG channel, to the fundamental cellular and tissue-level arrhythmia mechanisms. Applications of machine learning indicate necessary and sufficient parameters that predict arrhythmia vulnerability. We expect that the model framework may be expanded to make an impact in drug discovery, drug safety screening for a variety of compounds and targets, and in a variety of regulatory processes.

Simulating ion channel activation mechanisms using swarms of trajectories

Atomic-level studies of protein activity represent a significant challenge as a result of the complexity of conformational changes occurring on wide-ranging timescales, often greatly exceeding that of even the longest simulations. A prime example is the elucidation of protein allosteric mechanisms, where localized perturbations transmit throughout a large macromolecule to generate a response signal. For example, the conversion of chemical to electrical signals during synaptic neurotransmission in the brain is achieved by specialized membrane proteins called pentameric ligand-gated ion channels. Here, the binding of a neurotransmitter results in a global conformational change to open an ion-conducting pore across the nerve cell membrane. X-ray crystallography has produced static structures of the open and closed states of the proton-gated GLIC pentameric ligand-gated ion channel protein, allowing for atomistic simulations that can uncover changes related to activation. We discuss a range of enhanced sampling approaches that could be used to explore activation mechanisms. In particular, we describe recent application of an atomistic string method, based on Roux's "swarms of trajectories" approach, to elucidate the sequence and interdependence of conformational changes during activation. We illustrate how this can be combined with transition analysis and Brownian dynamics to extract thermodynamic and kinetic information, leading to understanding of what controls ion channel function. © 2019 Wiley Periodicals, Inc.

Accelerating the weighted histogram analysis method by direct inversion in the iterative subspace

The weighted histogram analysis method (WHAM) for free energy calculations is a valuable tool to produce free energy differences with the minimal errors. Given multiple simulations, WHAM obtains from the distribution overlaps the optimal statistical estimator of the density of states, from which the free energy differences can be computed. The WHAM equations are often solved by an iterative procedure. In this work, we use a well-known linear algebra algorithm which allows for more rapid convergence to the solution. We find that the computational complexity of the iterative solution to WHAM and the closely-related multiple Bennett acceptance ratio (MBAR) method can be improved by using the method of direct inversion in the iterative subspace. We give examples from a lattice model, a simple liquid and an aqueous protein solution.

Molecular simulation of nonfacilitated membrane permeation

This is a review. Non-electrolytic compounds typically cross cell membranes by passive diffusion. The rate of permeation is dependent on the chemical properties of the solute and the composition of the lipid bilayer membrane. Predicting the permeability coefficient of a solute is important in pharmaceutical chemistry and toxicology. Molecular simulation has proven to be a valuable tool for modeling permeation of solutes through a lipid bilayer. In particular, the solubility-diffusion model has allowed for the quantitative calculation of permeability coefficients. The underlying theory and computational methods used to calculate membrane permeability are reviewed. We also discuss applications of these methods to examine the permeability of solutes and the effect of membrane composition on permeability. The application of coarse grain and polarizable models is discussed. This article is part of a Special Issue entitled: Membrane Proteins edited by J.C. Gumbart and Sergei Noskov.

Modeling ion channels: past, present, and future

Ion channels are membrane-bound enzymes whose catalytic sites are ion-conducting pores that open and close (gate) in response to specific environmental stimuli. Ion channels are important contributors to cell signaling and homeostasis. Our current understanding of gating is the product of 60 plus years of voltage-clamp recording augmented by intervention in the form of environmental, chemical, and mutational perturbations. The need for good phenomenological models of gating has evolved in parallel with the sophistication of experimental technique. The goal of modeling is to develop realistic schemes that not only describe data, but also accurately reflect mechanisms of action. This review covers three areas that have contributed to the understanding of ion channels: traditional Eyring kinetic theory, molecular dynamics analysis, and statistical thermodynamics. Although the primary emphasis is on voltage-dependent channels, the methods discussed here are easily generalized to other stimuli and could be applied to any ion channel and indeed any macromolecule.

Molecular strategies to achieve selective conductance in NaK channel variants

A recent crystallization of several ion channels has provided strong impetus for efforts aimed at understanding the different strategies employed by nature for selective ion transport. In this work, we used two variants of the selectivity filter of NaK channel to explore molecular mechanisms that give rise to K(+)-selectivity. We computed one-dimensional (1D) and two-dimensional (2D) potentials of mean force (PMFs) for ion permeation across the channel. The results indicate that the energies for Na(+) and K(+) permeation across the selectivity filter display significant differences in positions of the binding sites and barriers. One characteristic signature of a K(+)-selective channel is the apparent preservation of the site analogous to that of S2 in KcsA. The S2-bound ion can be almost ideally dehydrated and coordinated by 6 to 8 carbonyls. In a striking contrast, the PMFs controlling transport of ions in a nonselective variant show almost identical profiles for either K(+) or Na(+) and significant involvement of water molecules in ion coordination across the entire selectivity filter. An analysis of differences in 1D PMFs for Na(+) and K(+) suggests that coordination number alone is an insufficient predictor of site selectivity, while chemical composition (ratio of carbonyls and water molecules) correlates well with preference for K(+). Multi-ion effects such as dependence of the barriers and wells for permeant ion on the type of copermeant ion were found to play a significant role in the selectivity signature of the channel as well.

Internal Diffusion-Controlled Enzyme Reaction: The Acetylcholinesterase Kinetics

Acetylcholinesterase is an enzyme with a very high turnover rate; it quenches the neurotransmitter, acetylcholine, at the synapse. We have investigated the kinetics of the enzyme reaction by calculating the diffusion rate of the substrate molecule along an active site channel inside the enzyme from atomic-level molecular dynamics simulations. In contrast to the previous works, we have found that the internal substrate diffusion is the determinant of the acetylcholinesterase kinetics in the low substrate concentration limit. Our estimate of the overall bimolecular reaction rate constant for the enzyme is in good agreement with the experimental data. In addition, the present calculation provides a reasonable explanation for the effects of the ionic strength of solution and the mutation of surface residues of the enzyme. The study suggests that internal diffusion of the substrate could be a key factor in understanding the kinetics of enzymes of similar characteristics.

CHARMM: the biomolecular simulation program

CHARMM (Chemistry at HARvard Molecular Mechanics) is a highly versatile and widely used molecular simulation program. It has been developed over the last three decades with a primary focus on molecules of biological interest, including proteins, peptides, lipids, nucleic acids, carbohydrates, and small molecule ligands, as they occur in solution, crystals, and membrane environments. For the study of such systems, the program provides a large suite of computational tools that include numerous conformational and path sampling methods, free energy estimators, molecular minimization, dynamics, and analysis techniques, and model-building capabilities. The CHARMM program is applicable to problems involving a much broader class of many-particle systems. Calculations with CHARMM can be performed using a number of different energy functions and models, from mixed quantum mechanical-molecular mechanical force fields, to all-atom classical potential energy functions with explicit solvent and various boundary conditions, to implicit solvent and membrane models. The program has been ported to numerous platforms in both serial and parallel architectures. This article provides an overview of the program as it exists today with an emphasis on developments since the publication of the original CHARMM article in 1983.

Molecular dynamics simulation of water in cytochrome c oxidase reveals two water exit pathways and the mechanism of transport

We have examined the network of connected internal cavities in cytochrome c oxidase along which water produced at the catalytic center is removed from the enzyme. Using combination of structural analysis, molecular dynamics simulations, and free energy calculations we have identified two exit pathways that connect the Mg2+ ion cavity to the outside of the enzyme. Each pathway has a well-defined bottleneck, which determines the overall rate of water traffic along the exit pathway, and a specific cooperative mechanism of passing it. One of the pathways is going via Arg438/439 (in bovine numbering) toward the CuA center, approaching closely its His204B ligand and Lys171B residue; and the other is going toward Asp364 and Thr294. Comparison of the pathways among different aa3-type enzymes shows that they are well conserved. Possible connections of the finding to redox-coupled proton pumping mechanism are discussed. We propose specific mutations near the bottlenecks of the exit pathways that can test some of our hypotheses.

Using stochastic models calibrated from nanosecond nonequilibrium simulations to approximate mesoscale information

We demonstrate how the surrogate process approximation (SPA) method can be used to compute both the potential of mean force along a reaction coordinate and the associated diffusion coefficient using a relatively small number (10-20) of bidirectional nonequilibrium trajectories coming from a complex system. Our method provides confidence bands which take the variability of the initial configuration of the high-dimensional system, continuous nature of the work paths, and thermal fluctuations into account. Maximum-likelihood-type methods are used to estimate a stochastic differential equation (SDE) approximating the dynamics. For each observed time series, we estimate a new SDE resulting in a collection of SPA models. The physical significance of the collection of SPA models is discussed and methods for exploiting information in the population of estimated SPA models are demonstrated and suggested. Molecular dynamics simulations of potassium ion dynamics inside a gramicidin A channel are used to demonstrate the methodology, although SPA-type modeling has also proven useful in analyzing single-molecule experimental time series [J. Phys. Chem. B 113, 118 (2009)].

The gramicidin ion channel: A model membrane protein

The linear peptide gramicidin forms prototypical ion channels specific for monovalent cations and has been extensively used to study the organization, dynamics and function of membrane-spanning channels. In recent times, the availability of crystal structures of complex ion channels has challenged the role of gramicidin as a model membrane protein and ion channel. This review focuses on the suitability of gramicidin as a model membrane protein in general, and the information gained from gramicidin to understand lipid-protein interactions in particular. Special emphasis is given to the role and orientation of tryptophan residues in channel structure and function and recent spectroscopic approaches that have highlighted the organization and dynamics of the channel in membrane and membrane-mimetic media.

Conduction properties of KcsA measured using brownian dynamics with flexible carbonyl groups in the selectivity filter

In the narrow segment of an ion conducting pathway, it is likely that a permeating ion influences the positions of the nearby atoms that carry partial or full electronic charges. Here we introduce a method of incorporating the motion of charged atoms lining the pore into Brownian dynamics simulations of ion conduction. The movements of the carbonyl groups in the selectivity filter of the KcsA channel are calculated explicitly, allowing their bond lengths, bond angles, and dihedral angels to change in response to the forces acting upon them. By systematically changing the coefficients of bond stretching and of angle bending, the carbon and oxygen atoms can be made to fluctuate from their fixed positions by varying mean distances. We show that incorporating carbonyl motion in this way does not alter the mechanism of ion conduction and only has a small influence on the computed current. The slope conductance of the channel increases by approximately 25% when the root mean-square fluctuations of the carbonyl groups are increased from 0.01 to 0.61 A. The energy profiles and the number of resident ions in the channel remain unchanged. The method we utilized here can be extended to allow the movement of glutamate or aspartate side chains lining the selectivity filters of other ionic channels.

Proton transfer in water wires in proteins: modulation by local constraint and polarity in gramicidin a channels

The transfer of protons in membrane proteins is an essential phenomenon in biology. However, the basic rules by which H(+) transfer occurs in water wires inside proteins are not well characterized. In particular, the effects of specific atoms and small groups of atoms on the rate of H(+) transfer in water wires are not known. In this study, new covalently linked gramicidin-A (gA) peptides were synthesized, and the effects of specific atoms and peptide constraints on the rate of H(+) transfer were measured in single molecules. The N-termini of two gA peptides were linked to various molecules: S,S-cyclopentane diacid, R,R-cyclopentane diacid, and succinic acid. Single-channel proton conductances (g(H)) were measured at various proton concentrations ([H(+)]) and compared to previous measurements obtained in the S,S- and R,R-dioxolane-linked as well as in native gA channels. Replacing the S,S-dioxolane by an S,S-cyclopentane had no effects on the g(H)-[H(+)] relationships, suggesting that the constrained and continuous transition between the two gA peptides via these S,S linkers is ultimately responsible for the two- to fourfold increase in g(H) relative to native gA channels. It is likely that constraining a continuous transition between the two gA peptides enhances the rate of H(+) transfer in water wires by decreasing the number of water wire configurations that do not transfer H(+) at higher rates as in native gA channels (a decrease in the activation entropy of the system). On the other hand, g(H) values in the R,R-cyclopentane are considerably larger than those in R,R-dioxolane-linked gA channels. One explanation would be that the electrostatic interactions between the oxygens in the dioxolane and adjacent carbonyls in the R,R-dioxolane-linked gA channel attenuate the rate of H(+) transfer in the middle of the pore. Interestingly, g(H)-[H(+)] relationships in the R,R-cyclopentane-linked gA channel are quite similar to those in native gA channels. g(H) values in succinyl-linked gA channels display a wide distribution of values that is well represented by a bigaussian. The larger peaks of these distributions are similar to g(H) values measured in native gA channel. This observation is also consistent with the notion that constraining the transition between the two beta-helical gA peptides enhances the rate of H(+) transfer in water wires by decreasing the activation entropy of the system.

On the thermodynamic stability of a charged arginine side chain in a transmembrane helix

Biological membranes consist of bilayer arrangements of lipids forming a hydrophobic core that presents a physical barrier to all polar and charged molecules. This long-held notion has recently been challenged by biological translocon-based experiments that report small apparent free energies to insert charged side chains near the center of a transmembrane (TM) helix. We have carried out fully atomistic simulations to provide the free-energy profile for moving a TM helix containing a protonated Arg side chain across a lipid bilayer. Our results reveal the fundamental thermodynamics governing the stability of charged side chains in membranes and the microscopic interactions involved. Despite local membrane deformations, where large amounts of water and lipid head groups are pulled into the bilayer to interact with Arg, the free-energy barrier is 17 kcal/mol. We provide a rationale for the differences in our microscopic free energies and cell biological experiments using free-energy calculations that indicate that a protonated Arg at the central residue of a TM helix of the Leader peptidase might reside close to the interface and not at the membrane center. Our findings have implications for the gating mechanisms of voltage-gated ion channels, suggesting that movements of protonated Arg residues through the membrane will be prohibited.

Framework models of ion permeation through membrane channels and the generalized King-Altman method

A modern approach to studying the detailed dynamics of biomolecules is to simulate them on computers. Framework models have been developed to incorporate information from these simulations in order to calculate properties of the biomolecules on much longer time scales than can be achieved by the simulations. They also provide a simple way to think about the simulated dynamics. This article develops a method for the solution of framework models, which generalizes the King-Altman method of enzyme kinetics. The generalized method is used to construct solutions of two framework models which have been introduced previously, the single-particle and Grotthuss (proton conduction) models. The solution of the Grotthuss model is greatly simplified in comparison with direct integration. In addition, a new framework model is introduced, generalizing the shaking stack model of ion conduction through the potassium channel.

Ion permeation through a narrow channel: using gramicidin to ascertain all-atom molecular dynamics potential of mean force methodology and biomolecular force fields

We investigate methods for extracting the potential of mean force (PMF) governing ion permeation from molecular dynamics simulations (MD) using gramicidin A as a prototypical narrow ion channel. It is possible to obtain well-converged meaningful PMFs using all-atom MD, which predict experimental observables within order-of-magnitude agreement with experimental results. This was possible by careful attention to issues of statistical convergence of the PMF, finite size effects, and lipid hydrocarbon chain polarizability. When comparing the modern all-atom force fields of CHARMM27 and AMBER94, we found that a fairly consistent picture emerges, and that both AMBER94 and CHARMM27 predict observables that are in semiquantitative agreement with both the experimental conductance and dissociation coefficient. Even small changes in the force field, however, result in significant changes in permeation energetics. Furthermore, the full two-dimensional free-energy surface describing permeation reveals the location and magnitude of the central barrier and the location of two binding sites for K(+) ion permeation near the channel entrance--i.e., an inner site on-axis and an outer site off-axis. We conclude that the MD-PMF approach is a powerful tool for understanding and predicting the function of narrow ion channels in a manner that is consistent with the atomic and thermally fluctuating nature of proteins.

Molecular dynamics - potential of mean force calculations as a tool for understanding ion permeation and selectivity in narrow channels

Ion channels catalyze the permeation of charged molecules across cell membranes and are essential for many vital physiological functions, including nerve and muscle activity. To understand better the mechanisms underlying ion conduction and valence selectivity of narrow ion channels, we have employed free energy techniques to calculate the potential of mean force (PMF) for ion movement through the prototypical gramicidin A channel. Employing modern all-atom molecular dynamics (MD) force fields with umbrella sampling methods that incorporate one hundred 1-2 ns trajectories, we find that it is possible to achieve semi-quantitative agreement with experimental binding and conductance measurements. We also examine the sensitivity of the MD-PMF results to the choice of MD force field and compare PMFs for potassium, calcium and chloride ions to explore the basis for the valence selectivity of this narrow and uncharged ion channel. A large central barrier is observed for both anions and divalent ions, consistent with lack of experimental conductance. Neither anion or divalent cation is seen to be stabilized inside the channel relative to the bulk electrolyte and each leads to large disruptions to the protein and membrane structure when held deep inside the channel. Weak binding of calcium ions outside the channel corresponds to a free energy well that is too shallow to demonstrate channel blocking. Our findings emphasize the success of the MD-PMF approach and the sensitivity of ion energetics to the choice of biomolecular force field.

Diffusion constant of K+ inside Gramicidin A: a comparative study of four computational methods

The local diffusion constant of K(+) inside the Gramicidin A (GA) channel has been calculated using four computational methods based on molecular dynamics (MD) simulations, specifically: Mean Square Displacement (MSD), Velocity Autocorrelation Function (VACF), Second Fluctuation Dissipation Theorem (SFDT) and analysis of the Generalized Langevin Equation for a Harmonic Oscillator (GLE-HO). All methods were first tested and compared for K(+) in bulk water--all predicted the correct diffusion constant. Inside GA, MSD and VACF methods were found to be unreliable because they are biased by the systematic force exerted by the membrane-channel system on the ion. SFDT and GLE-HO techniques properly unbias the influence of the systematic force on the diffusion properties and predicted a similar diffusion constant of K(+) inside GA, namely, ca. 10 times smaller than in the bulk. It was found that both SFDT and GLE-HO methods require extensive MD sampling on the order of tens of nanoseconds to predict a reliable diffusion constant of K(+) inside GA.

Et tu, Grotthuss! and other unfinished stories

This review article is divided into three sections. In Section 1, a short biographical note on Freiherr von Grotthuss is followed by a detailed summary of the main findings and ideas present in his 1806 paper. Attempts to place Grotthuss contribution in the context of the science done at his time were also made. In Section 2, the modern version of the Grotthuss mechanism is reviewed. The classical Grotthuss model has been recently questioned and new mechanisms and ideas regarding proton transfer are briefly discussed. The last section discusses the significance of a classical Grotthuss mechanism for proton transfer in water chains inside protein cavities. This has been an interesting new twist in the ongoing history of the Grotthuss mechanism. A summary and discussion of what was learned from probably the simplest currently available experimental models of proton transfer in water wires in semi-synthetic ion channels are critically presented. This review ends discussing some of the questions that need to be addressed in the near future.

A gate in the selectivity filter of potassium channels

The selectivity filter of potassium channels is the structural element directly responsible for the selective and rapid conduction of K+, whereas other parts of the protein are thought to function as a molecular gate that either permits or blocks the passage of ions. However, whether the selectivity filter itself also possesses the ability to play the role of a gate is an unresolved question. Using free energy molecular dynamics simulations, it is shown that the reorientation of two peptide linkages in the selectivity filter of the KcsA K+ channel can lead to a stable nonconducting conformational state. Two microscopic factors influence the transition toward such a conformational state: the occupancy of one specific cation binding site in the selectivity filter (S2), and the strength of intersubunit interactions involving the GYG signature sequence. These results suggest that such conformational transitions occurring in the selectivity filter might be related to different K+ channel gating events, including C-type (slow) inactivation.

A framework model based on the Smoluchowski equation in two reaction coordinates

The general form of the Smoluchowski equation in two reaction coordinates is obtained as the diffusion limit of a random walk on an infinite square grid using transition probabilities that satisfy detailed balance at thermodynamic equilibrium. The diffusion limit is then used to construct a generalization of the single-particle model to two reaction coordinates. The state space includes a square on which diffusion takes place and an isolated empty state. Boundary conditions on opposite sides of the square correspond to transitions between the empty state and the square. The two-dimensional (2D) model can be reduced to a 1D single-particle model by adiabatic elimination. A finite element solution of the 2D boundary value problem is described. The method used to construct the 2D model can be adapted to state spaces that have been constructed by other authors to model K+ conduction through gramicidin, proton conduction through dioxolane-linked gramicidin, and chloride conduction through the bacterial H(+)-Cl- antiporter.

Structure of gramicidin a in a lipid bilayer environment determined using molecular dynamics simulations and solid-state NMR data

Two different high-resolution structures recently have been proposed for the membrane-spanning gramicidin A channel: one based on solid-state NMR experiments in oriented phospholipid bilayers (Ketchem, R. R.; Roux, B.; Cross, T. A. Structure 1997, 5, 1655-1669; Protein Data Bank, PDB:1MAG); and one based on two-dimensional NMR in detergent micelles (Townsley, L. E.; Tucker, W. A.; Sham, S.; Hinton, J. F. Biochemistry 2001, 40, 11676-11686; PDB:1JNO). Despite overall agreement, the two structures differ in peptide backbone pitch and the orientation of several side chains; in particular that of the Trp at position 9. Given the importance of the peptide backbone and Trp side chains for ion permeation, we undertook an investigation of the two structures using molecular dynamics simulation with an explicit lipid bilayer membrane, similar to the system used for the solid-state NMR experiments. Based on 0.1 micros of simulation, both backbone structures converge to a structure with 6.25 residues per turn, in agreement with X-ray scattering, and broad agreement with SS backbone NMR observables. The side chain of Trp 9 is mobile, more so than Trp 11, 13, and 15, and undergoes spontaneous transitions between the orientations in 1JNO and 1MAG. Based on empirical fitting to the NMR results, and umbrella sampling calculations, we conclude that Trp 9 spends 80% of the time in the 1JNO orientation and 20% in the 1MAG orientation. These results underscore the utility of molecular dynamics simulations in the analysis and interpretation of structural information from solid-state NMR.

A microscopic view of ion conduction through the K+ channel

Recent results from x-ray crystallography and molecular dynamics free-energy simulations have revealed the existence of a number of specific cation-binding sites disposed along the narrow pore of the K+ channel from Streptomyces lividans (KcsA), suggesting that K+ ions might literally "hop" in single file from one binding site to the next as permeation proceeds. In support of this view, it was found that the ion configurations correspond to energy wells of similar depth and that ion translocation is opposed only by small energy barriers. Although such features of the multiion potential energy surface are certainly essential for achieving a high throughput rate, diffusional and dissipative dynamical factors must also be taken into consideration to understand how rapid conduction of K+ is possible. To elucidate the mechanism of ion conduction, we established a framework theory enabling the direct simulation of nonequilibrium fluxes by extending the results of molecular dynamics over macroscopically long times. In good accord with experimental measurements, the simulated maximum conductance of the channel at saturating concentration is on the order of 550 and 360 pS for outward and inward ions flux, respectively, with a unidirectional flux-ratio exponent of 3. Analysis of the ion-conduction process reveals a lack of equivalence between the cation-binding sites in the selectivity filter.

Theoretical study of the structure and dynamic fluctuations of dioxolane-linked gramicidin channels

Gramicidin is a hydrophobic peptide that assembles as a head-to-head dimer in lipid membranes to form water-filled channels selective to small monovalent cations. Two diastereoisomeric forms, respectively SS and RR, of chemically modified channels in which a dioxolane ring links the formylated N-termini of two gramicidin monomers, were shown to form ion channels. To investigate the structural basis underlying experimentally measured differences in proton conductance in the RR and SS channels, we construct atomic-resolution models of dioxolane-linked gramicidin dimers by analogy with the native dimer. A parametric description of the linker compatible with the CHARMM force field used for the peptide is derived by fitting geometry, vibrational frequencies, and energy to the results of ab initio calculations. The linker region of the modified gramicidin dimers is subjected to an extensive conformational search using high-temperature simulated annealing, and free-energy surfaces underlying the structural fluctuations of the channel backbone at 298K are computed from molecular dynamics simulations. The overall secondary structure of the beta-helical gramicidin pore is retained in both linked channels. The SS channel is found in a single conformation resembling that of the native dimer, with its peptide bonds undergoing rapid librations with respect to the channel axis. By contrast, its RR counterpart is characterized by local backbone distortions in which the two peptide bonds flanking the linker are markedly tilted in order to satisfy the pitch of the helix. In these distorted structures, each of the two carbonyl groups points either in or out of the lumen. Flipping these two peptides in and out involves thermally activated transitions, which results in four distinct conformational states at equilibrium with one another on a nanosecond time scale. This work opens the way to detailed comparative studies of structure-function relationships in biological proton ducts.

Computational studies of the gramicidin channel

Ion channels are highly specific membrane-spanning protein structures which serve to facilitate the passage of selected ions across the lipid barrier. In the past decade, molecular dynamics simulations based on atomic models and realistic microscopic interactions with explicit solvent and membrane lipids have been used to gain insight into the function of these complex systems. These calculations have considerably expanded our view of ion permeation at the microscopic level. This Account will mainly focus on computational studies of the gramicidin A channel, one of the simplest and best characterized molecular pore.

On the origin of closing flickers in gramicidin channels: a new hypothesis

The submillisecond closing events (flickers) and the single channel conductances to protons (g(H)) were studied in native gramicidin A (gA) and in the SS and RR diastereoisomers of dioxolane-linked gA channels in planar bilayers. Bilayers were formed from glycerylmonooleate (GMO) in various solvents. In GMO/decane (thick) bilayers, the largest flicker frequency occurred in the SS channel (39 s(-1)), followed by the RR (4 s(-1)) and native gA channels (3 s(-1)). These frequencies were attenuated in GMO/squalene (thin) bilayers by 100-, 30-, and 70-fold in the SS, RR, and native gA channels, respectively. In thin bilayers, the average burst duration of native gA channels was 30-fold longer than in thick bilayers. The RR dioxolane-linked gA dimer "inactivated" in GMO/decane but not in squalene-containing bilayers. The mean closed time of flickers (approximately 0.12 ms) was essentially the same in various gA channels. In thin bilayers, g(H) values were larger by approximately 10% (SS), 30% (RR), and 20% (native gA) in relation to thick bilayers. It is concluded that flickers are not related to pre-dissociation or dissociation states of gA monomers, and do not seem to be caused by intrinsic conformational changes of channel proteins. It is proposed that flickers are caused by undulations of the bilayer that obliterate the openings of gA channels. Differences between flicker frequencies in various gA channels are likely to result from variations in channel geometries at the bilayer/channel interface. The smaller g(H) in thick bilayers suggests that the deformation of these bilayers around the gA channel creates a diffusional pathway next to the mouths of the channel that is longer and more restrictive than in thin GMO bilayers. A possible molecular interpretation for these effects is attempted.

Covalently linked gramicidin channels: effects of linker hydrophobicity and alkaline metals on different stereoisomers

The direct role of the dioxolane group on the gating and single-channel conductance of different stereoisomers of the dioxolane-linked gramicidin A (gA) channels reconstituted in planar lipid bilayers was investigated. Four different covalently linked gA dimers were synthesized. In two of them, the linker was the conventional dioxolane described previously (SS and RR channels). Two gAs were covalently linked with a novel modified dioxolane group containing a retinal attachment (ret-SS and ret-RR gA dimers). These proteins also formed ion channels in lipid bilayers and were selective for monovalent cations. The presence of the bulky and hydrophobic retinal group immobilizes the dioxolane linker in the bilayer core preventing its rotation into the hydrophilic lumen of the pore. In 1 M HCl the gating kinetics of the SS or RR dimers were indistinguishable from their retinal counterparts; the dwell-time distributions of the open and closed states in the SS and ret-SS were basically the same. In particular, the inactivation of the RR was not prevented by the presence of the retinal group. It is concluded that neither the fast closing events in the SS or RR dimers nor the inactivation of the RR are likely to be a functional consequence of the flipping of the dioxolane inside the pore of the channel. On the other hand, the inactivation of the RR dimer was entirely eliminated when alkaline metals (Cs(+) or K(+)) were the permeating cations in the channel. In fact, the open state of the RR channel became extremely stable, and the gating characteristics of both the SS and RR channels were different from what was seen before with permeating protons. As in HCl, the presence of a retinal in the dioxolane linker did not affect the gating behavior of the SS and RR in Cs(+)- or K(+)-containing solutions. Alternative hypotheses concerning the gating of linked gA dimers are discussed.

A combined molecular dynamics and diffusion model of single proton conduction through gramicidin

We develop a model for proton conduction through gramicidin based on the molecular dynamics simulations of Pomès and Roux (Biophys. J. 72:A246, 1997). The transport of a single proton through the gramicidin pore is described by a potential of mean force and diffusion coefficient obtained from the molecular dynamics. In addition, the model incorporates the dynamics of a defect in the hydrogen bonding structure of pore waters without an excess proton. Proton entrance and exit were not simulated by the molecular dynamics. The single proton conduction model includes a simple representation of these processes that involves three free parameters. A reasonable value can be chosen for one of these, and the other two can be optimized to yield a good fit to the proton conductance data of, Ann. N.Y. Acad. Sci. 339:8-20) for pH > or = 1.7. A sensitivity analysis shows the significance of this fit.

Proton mobilities in water and in different stereoisomers of covalently linked gramicidin A channels

Proton conductivities in bulk solution (lambda(H)) and single-channel proton conductances (g(H)) in two different stereoisomers of the dioxolane-linked gramicidin A channel (the SS and RR dimers) were measured in a wide range of bulk proton concentrations ([H], 0.1-8000 mM). Proton mobilities (micro(H)) in water as well as in the SS and RR dimers were calculated from the conductivity data. In the concentration range of 0.1-2000 mM, a straight line with a slope of 0.75 describes the log (g(H))-log ([H]) relationship in the SS dimer. At [H] > 2000 mM, saturation is followed by a decline in g(H). The g(H)-[H] relationship in the SS dimer is qualitatively similar to the [H] dependence of lambda(H). However, the slope of the straight line in the log(lambda(H))-log([H]) plot is 0.96, indicating that the rate-limiting step for proton conduction through the SS dimer is not the diffusion of protons in bulk solution. The significant difference between the slopes of those linear relationships accounts for the faster decline of micro(H) as a function of [H] in the SS dimer in relation to bulk solution. In the high range of [H], saturation and decline of g(H) in the SS dimer can be accounted for by the significant decrease of micro(H) in bulk solution. At any given [H], g(H) in the RR dimer is significantly smaller than in the SS. Moreover, the g(H)-[H] relationship in the RR stereoisomer is qualitatively different from that in the SS. Between 1 and 50 mM [H], g(H) can be fitted with an adsorption isotherm, suggesting the presence of a proton-binding site inside the pore (pK(a) approximately 2), which limits proton exit from the channel. At 100 mM < [H] < 3000 mM, g(H) increases linearly with [H]. The distinctive shape of the g(H)-[H] relationship in the RR dimer suggests that the channel can be occupied simultaneously by more than one proton. At higher [H], the saturation and decline of g(H) in the RR dimer reflect the properties of micro(H) in bulk solution. In the entire range of [H], protons seem to cross the SS and RR channels via a Grotthuss-like mechanism. The rate-limiting step for proton transfer in the SS dimer is probably the membrane-channel/bulk solution interface. It is also proposed that the smaller g(H) in the RR dimer is the consequence of a different organization and dynamics of the H-bonded network of water molecules inside the pore of the channel, resulting in a slower proton transfer and multiple pore occupancy by protons.

The conduction of protons in different stereoisomers of dioxolane-linked gramicidin A channels

Two different stereoisomers of the dioxolane-linked gramicidin A (gA) channels were individually synthesized (the SS and RR dimers;. Science. 244:813-817). The structural differences between these dimers arise from different chiralities within the dioxolane linker. The SS dimer mimics the helicity and the inter- and intramolecular hydrogen bonding of the monomer-monomer association of gA's. In contrast, there is a significant disruption of the helicity and hydrogen bonding pattern of the ion channel in the RR dimer. Single ion channels formed by the SS and RR dimers in planar lipid bilayers have different proton transport properties. The lipid environment in which the different dimers are reconstituted also has significant effects on single-channel proton conductance (g(H)). g(H) in the SS dimer is about 2-4 times as large as in the RR. In phospholipid bilayers with 1 M [H(+)](bulk), the current-voltage (I-V) relationship of the SS dimer is sublinear. Under identical experimental conditions, the I-V plot of the RR dimer is supralinear (S-shaped). In glycerylmonooleate bilayers with 1 M [H(+)](bulk), both the SS and RR dimers have a supralinear I-V plot. Consistent with results previously published (. Biophys. J. 73:2489-2502), the SS dimer is stable in lipid bilayers and has fast closures. In contrast, the open state of the RR channel has closed states that can last a few seconds, and the channel eventually inactivates into a closed state in either phospholipid or glycerylmonooleate bilayers. It is concluded that the water dynamics inside the pore as related to proton wire transfer is significantly different in the RR and SS dimers. Different physical mechanisms that could account for this hypothesis are discussed. The gating of the synthetic gA dimers seems to depend on the conformation of the dioxolane link between gA's. The experimental results provide an important framework for a detailed investigation at the atomic level of proton conduction in different and relatively simple ion channel structures.

Proton conduction in gramicidin A and in its dioxolane-linked dimer in different lipid bilayers

Gramicidin A (gA) molecules were covalently linked with a dioxolane ring. Dioxolane-linked gA dimers formed ion channels, selective for monovalent cations, in planar lipid bilayers. The main goal of this study was to compare the functional single ion channel properties of natural gA and its covalently linked dimer in two different lipid bilayers and HCl concentrations (10-8000 mM). Two ion channels with different gating and conductance properties were identified in bilayers from the product of dimerization reaction. The most commonly observed and most stable gramicidin A dimer is the main object of this study. This gramicidin dimer remained in the open state most of the time, with brief closing flickers (tau(closed) approximately 30 micros). The frequency of closing flickers increased with transmembrane potential, making the mean open time moderately voltage dependent (tau(open) changed approximately 1.43-fold/100 mV). Such gating behavior is markedly different from what is seen in natural gA channels. In PEPC (phosphatidylethanolamine-phosphatidylcholine) bilayers, single-channel current-voltage relationships had an ohmic behavior at low voltages, and a marked sublinearity at relatively higher voltages. This behavior contrasts with what was previously described in GMO (glycerylmonooleate) bilayers. In PEPC bilayers, the linear conductance of single-channel proton currents at different proton concentrations was essentially the same for both natural and gA dimers. g(max) and K(D), obtained from fitting experimental points to a Langmuir adsorption isotherm, were approximately 1500 pS and 300 mM, respectively, for both the natural gA and its dimer. In GMO bilayers, however, proton affinities of gA and the dioxolane-dimer were significantly lower (K(D) of approximately 1 and 1.5 M, respectively), and the g(max) higher (approximately 1750 and 2150 pS, respectively) than in PEPC bilayers. Furthermore, the relationship between single-channel conductance and proton concentration was linear at low bulk concentrations of H+ (0.01-2 M) and saturated at concentrations of more than 3 M. It is concluded that 1) The mobility of protons in gramicidin A channels in different lipid bilayers is remarkably similar to proton mobilities in aqueous solutions. In particular, at high concentrations of HCl, proton mobilities in gramicidin A channel and in solution differ by only 25%. 2) Differences between proton conductances in gramicidin A channels in GMO and PEPC cannot be explained by surface charge effects on PEPC membranes. It is proposed that protonated phospholipids adjacent to the mouth of the pore act as an additional source of protons for conduction through gA channels in relation to GMO bilayers. 3) Some experimental results cannot be reconciled with simple alterations in access resistance to proton flow in gA channels. Said differences could be explained if the structure and/or dynamics of water molecules inside gramicidin A channels is modulated by the lipid environment and by modifications in the structure of gA channels. 4) The dioxolane ring is probably responsible for the closing flickers seen in the dimer channel. However, other factors can also influence closing flickers.

Molecular dynamics study of free energy profiles for organic cations in gramicidin A channels

The free energy profiles for four organic cations in right-handed single-helix gramicidin A dimers were computed by using umbrella sampling molecular dynamics with CHARMM. Ion-water column translocations were facilitated by using a novel "water-tunnel" approach. The overlapping pieces of free energy profile for adjacent windows were selected from three trajectories that differed in initial ion rotation and were aligned by the method of umbrella potential differences. Neglected long-range electrostatic energies from the bulk water and the bilayer were computed with DelPhi and added to the profile. The approach was corroborated for the formamidinium-guanidinium pair by using perturbation dynamics at axial positions 0, 6, 12, and 15 A from the channel center. The barrier to ethylammonium entry was prohibitive at 21 kcal/mol, whereas for methylammonium it was 5.5 kcal/mol, and the profile was quite flat through the channel, roughly consistent with conductance measurements. The profile for formamidinium was very similar to that of methylammonium. Guanidinium had a high entry barrier (deltaF = +8.6 kcal/mol) and a narrow deep central well (deltaF = -2.6 kcal/mol), qualitatively consistent with predictions from voltage-dependent potassium current block measurements. Its deep central well, contrasting with the flat profile for formamidinium, was verified with perturbation dynamics and was correlated with its high propensity to form hydrogen bonds with the channel at the dimer junction (not shared by the other three cations). Analysis of the ensemble average radial forces on the ions demonstrates that all four ions undergo compressive forces in the channel that are at maximum at the center of the monomer and relieved at the dimer junction, illustrating increased flexibility of the channel walls in the center of the channel.

Permeation through the calcium release channel of cardiac muscle

Current voltage (I-V) relations were measured from the calcium release channel (CRC) of the sarcoplasmic reticulum of cardiac muscle in 12 KCl solutions, symmetrical and asymmetrical, from 25 mM to 2 M. I-V curves are nearly linear, in the voltage range +/- 150 mV approximately 12kT/e, even in asymmetrical solutions, e.g., 2 M // 100 mM. It is awkward to describe straight lines as sums of exponentials in a wide range of solutions and potentials, and so traditional barrier models have difficulty fitting this data. Diffusion theories with constant fields predict curvilinear I-V relations, and so they are also unsatisfactory. The Poisson and Nernst-Planck equations (PNP) form a diffusion theory with variable fields. They fit the data by using adjustable parameters for the diffusion constant of each ion and for the effective density of fixed (i.e., permanent) charge P(x) along the channel's "filter" (7-A diameter, 10 A long). If P(x) is described by just one parameter, independent of x (i.e., P(x) = P0 = -4.2 M), the fits are satisfactory (RMS error/RMS current = 6.4/67), and the estimates of diffusion coefficients are reasonable D(K) = 1.3 x 10(-6) cm2/s, D(Cl) = 3.9 x 10(-6) cm2/s. The CRC seems to have a small selectivity filter with a very high density of permanent charge. This may be a design principle of channels specialized for large flux. The Appendix derives barrier models, and their prefactor, from diffusion theories (with variable fields) and argues that barrier models are poor descriptions of CRCs in particular and open channels in general.

Molecular dynamics study of an alpha-cyclodextrin-phosphatidylinositol inclusion complex

Cyclodextrins are cyclic oligosaccharides known for their ability to include substrate molecules in their hydrophobic cavity. Moreover, cyclodextrins show a hemolytic activity when mM concentrations are added to blood. This hemolysis is commonly interpreted as a massive dissociation of phospholipids from the cell membrane due to the formation of complexes with the cyclodextrins. In the literature, a complexation between alpha-cyclodextrin (alpha CD) and phosphatidylinositol (PI) specific to the inositol headgroup has been proposed. But the need for the detailed interaction mechanism between the two molecules motivated the present work based on molecular dynamics simulations. Investigation of long range electrostatic interactions shows that a mutual approach of the molecules is only possible when the primary hydroxyl side of alpha CD faces the inositol headgroup of PI. This orientation is also the most favourable from adiabatic-and free-energy profiles calculated along a reaction coordinate that leads to an inclusion of PI into alpha CD. For free energy simulations, partial hydration of the model has been used. A study of glycosidic bond dihedral angles in alpha CD shows an increase in dihedral fluctuations before complexation and a dihedral "freezing" once the complex is formed.