Nonimmobilizers are agents that lack anesthetic properties, although their chemical structure is very similar to known anesthetics. The primary action site of both agents, whether at the membrane or target protein level, is still a matter of debate. However, increasing evidence points to the distinct modifications of the membrane physical properties that such agents induce. Such modification may play a role in the mechanism of anesthesia, and may therefore be related to the differences in their clinical behavior.
Molecular dynamics (MD) computer simulations have been used to investigate the distribution of a nonimmobilizer, hexafluroethane (HFE, C(2)F(6)), in a lipid membrane. The biologically relevant liquid-crystal phase of a hydrated dimyristoyl phosphatidyl choline (DMPC) bilayer was used as a membrane model. Two MD simulations corresponding to HFE mole fractions of 6% and 25% have been performed at room temperature and constant ambient pressure, for a duration of 2 nanoseconds each.
The equilibrium configurations of HFE in the bilayer show that the nonimmobilizer molecules are evenly distributed along the lipid hydrocarbon chains with a slight preference for the bilayer center. This partitioning induces an expansion of the bilayer thickness and a lateral contraction of the membrane (decrease of the area per lipid). The presence of HFE has essentially no effect on the lipid acyl chain conformations in agreement with nuclear magnetic resonance (NMR) measurements of the chain order parameters. The partitioning of the nonimmobilizer does not influence the orientation of the lipid head-group dipole moment.
The modifications induced by the presence of the nonimmobilizer HFE on a model membrane are distinct from those previously found for halothane (CF(3)CHBrCl), its anesthetic analogue, and appear to result from different distributions in the lipid bilayer. The results of the MD simulations show that (1) the changes in the average area per lipid and in the membrane thickness are opposite for the two agents and (2) HFE induces no change in the lipid head-group orientation, in contrast to halothane. These different effects (1) on the physical properties of the lipid bilayer and (2) on the electrostatic properties of the membrane-water interface may be linked to different clinical effects, and thus might contribute to the mechanism of general anesthesia.
MANY volatile compounds when inhaled induce transient analgesia, amnesia, and immobility in response to a noxious stimulus. 1This large and chemically diverse group of compounds is known as inhaled general anesthetics (IAs). A remarkable correlation of IA potency with oil solubility, known as the Meyer-Overton rule, 2has served for many years to predict anesthetic potency. Recently, a number of compounds with structures and lipophilicities similar to those of known IAs, and therefore predicted to be good anesthetics, were found to lack anesthetic properties. 3–5These molecules, known as nonimmobilizers (also called nonanesthetics) do not suppress movement, but cause amnesia. Because of the chemical similarity between IAs and nonimmobilizers, it was suggested that different mechanisms may be responsible for anesthesia and immobility. 6–9This has triggered an interest in comparative studies of the effects of nonimmobilizers versus anesthetics, aimed at better understanding the mechanism of IA action.
The site of anesthetic action is still uncertain and very little is known at the molecular level about the mechanism of IA action. While the lipid theory of narcosis that suggests that the site of action is the lipid membrane is supported by numerous observations, several other experimental studies shifted the focus from the lipids, suggesting that IA action may involve direct binding to proteins and perturbation of their function. For instance, transmembrane ligand–gated ion–channel proteins have been implicated in IA action. 10–12Most recently, renewed interest in the role of membrane lipids was ignited by a hypothesis from Cantor regarding modifications in membrane structure induced by the presence of IAs, which in turn may indirectly alter membrane protein function. 13–15
It is along these lines that recent NMR experiments were performed for various IAs and nonimmobilizers in diverse model membranes. 16These experiments measure the order parameters along the acyl chains of the lipids forming the bilayer, thus allowing characterization at the molecular level of the effects of these agents on the membrane. The overall picture emerging from the measurements indicates that anesthetics do perturb the structure of the membrane while nonimmobilizers have no effect on the acyl chain order parameters. This contrasting behavior is attributed to differences in the location of the agents in the membrane. The results of 19F NMR experiments 17indeed show that anesthetics reside preferentially near the lipid head-group while nonimmobilizers partition within the bilayer interior.
Computer simulation techniques have been a useful tool in supplying a molecular level description of the interactions of anesthetics and nonimmobilizers with lipids comprising model membranes. For instance, the difference in affinities of the agents for the water–membrane interface versus the lipid core was predicted by Pohorille et al. 18–20from energy profile distributions. Further details on the molecular distribution of anesthetics in lipid membranes and on their effects have been provided by molecular dynamics (MD) computer simulation studies. In an earlier work we reported the results obtained for halothane in a model membrane of hydrated dipalmitoylphosphatidylcholine (DPPC) lipid bilayer. 21–22We found that the molecules preferentially distribute and segregate below the lipid glycerol backbone, in the upper part of the lipid acyl chains. The presence of halothane induces a lateral expansion of the bilayer and a concomitant contraction in its thickness and in agreement with experiments, a decrease of the lipid chain order parameters. The aggregation of halothane in the upper region of the acyl chain also induces significant changes in the orientation of the head-group, implying measurable modifications in the electrostatic properties of the membrane–water interface.
In this paper, using MD simulations, we investigate the distribution and the effects of hexafluroethane (HFE), a nonimmobilizer analogue of halothane, on a dimyristoyl phosphatidyl choline (DMPC) model membrane. The aim of the study is to compare the effects of the nonimmobilizer and the anesthetic on the membrane properties.
Before proceeding with a detailed description of the methodologies used, we briefly review some basic concept and terminology of the MD simulation technique used in the present paper and refer the reader to the cited references for details. 23–26
Molecular dynamic simulations involve solving Newton's equations of motion for all the molecules in the system under consideration (i.e. , lipids, water, anesthetic), using forces derived from an intermolecular potential energy function, V . 23In the classic approach adopted here, the function V is represented in terms of interactions between pairs of atoms and is typically divided into “bonded” and “nonbonded” terms. The bonded interactions consist of harmonic bond stretching, angle bending, and dihedral angle deformation of individual molecules. Nonbonded interactions include electrostatic and van der Waals terms. The former are calculated from partial charges assigned to each atom, whereas the latter describe short-range repulsion and long-range attractions between pairs of atoms. A typical potential function V is of the form:27
Here, rijis the distance between atoms i and j ;q iis the partial charge on atom i ; ϵijand ςijare Lennard-Jones parameters for van der Waals interactions; kb, kθand kφand b 0, θ0and φ0are the force constants and equilibrium values for the bond streach, angle bend, and dihedral torsion deformations.
Molecular dynamic simulations generate a set of atomic positions and velocities as a function of time that evolve deterministically from an initial configuration according to the interaction potential V . Using statistical mechanics, the simulation results may be used to calculate observable quantities if the trajectory is long enough to yield satisfactory time averages. In practice, trajectories for complex systems such as lipid membranes span a few nanoseconds (10−9s). This is not long enough, in general, to permit complete relaxation of the system and determine local structural properties. The quality of the proper sampling of the trajectory may be assessed, in part, by comparing calculated quantities to experiment. 26
For the computation to be tractable, the number of lipids in the system is typically less than 100. To reduce the effects caused by the relatively small system size, the sample is placed in a central cell, called the simulations box, which is replicated infinitely in three-dimensions using periodic boundary conditions (PBCs). Every molecule in the box interacts therefore with its neighbors and their replicas. The long-range van der Waals interactions are typically turned off beyond a certain cutoff distance (generally, 10–15 Å, which is smaller than half the box size to avoid interactions between a molecule and its own image. Through the use of PBCs a simulation box containing a small patch of a hydrated lipid bilayer effectively corresponds to an infinite multilamellar system. 23
A fully hydrated DMPC lipid bilayer was used as the model membrane. This choice of the lipid is mainly motivated by the fact that the system is in the biologically relevant liquid crystal phase at room temperature. 28The initial configuration for the present simulations was taken from a previously well equilibrated pure DMPC bilayer containing 64 DMPC molecules (32 per leaflet) and 1645 water molecules. The simulated DMPC system is characterized by a surface area of A = 59.2Å2per lipid head-group, and a lamellar spacing d = 61.5Å, which is within 2% the experimental values. 28This excellent agreement, between the MD results and experimental structural data is taken as justification for the simulation protocols and the quality of the forcefield (potential functions, V ).
The initial configuration for the two MD simulations of HFE in DMPC at low and high concentrations were generated as follows: Four and sixteen HFE molecules (6.5% and 25% mole fractions, respectively), were incorporated in the lipid regions of DMPC. For the low concentration system, the molecules were initially placed near the lipid head group glycerol region. For the higher concentration, the HFE molecules were randomly distributed in the lipid core. These initial configurations were similar to those chosen previously in the halothane–DPPC simulations 21–22to allow comparison between the two systems. The HFE molecules were initially treated as small pointlike molecules with zero atomic charges and van der Waals parameters, which mediates possible unfavorable repulsive contacts with neighboring lipid molecules. Using short consecutive MD runs, the molecules were “grown” in size by extending their intramolecular bonds, and their interactions with the surroundings. 29The systems were then equilibrated for 200 ps (picoseconds) at fixed volume and temperature (T = 30°C) using an orthorhombic simulation cell, with three-dimensional PBCs. This was followed by approximately 2 nanoseconds trajectories at constant pressure (1 atmosphere [atm]) and temperature (T = 30°C), which allowed the system to adjust its desired structure. For the analysis of membrane properties, averaged quantities were calculated over the last 1.2 nanoseconds of the MD run. This averaging is taken over every 20 time steps and every lipid in the system.
The CHARMM 27force field parameters (functional form similar to in equation 1) were employed for the fully atomistic representation of DMPC lipid molecules and TIP3 parameters for water. 30A set of parameters for HFE were fitted based on our previous work on halothane. 31In accord with Tang et al. 32, the C and F atoms have relatively large partial charges (tables 1 and 2). MD simulations were performed for pure liquid HFE, considering a system containing 125 molecules at T =−90°C and constant ambient pressure. The potential parameters were fitted in a series of MD runs by comparing the simulation results to a number of physical properties of HFE such as density and heat of vaporization. The liquid state structural parameters calculated from an average over the last 500 ps of a 750 ps MD run are within 2% of the experimental results. 33,34
To carry out the MD simulations described herein we employed the recently developed simulation package (PINY-MD). 29This code uses the so-called Nosé-Hoover chain thermostat extended system method 35to generate the constant temperature and pressure trajectories. A recently developed reversible multiple time-step algorithm was used to speed up the calculation (4 fs MD time step). Electrostatic interactions between the atomic charges were calculated using the Particle-Mesh-Ewald (PME) method. 36The minimum image convention 23was employed to calculate the Lennard-Jones interactions and the real-space part of the Ewald sum, using a spherical truncation of 7Å and 10Å, for the short- and long-range parts of the force decomposition, respectively. Shake–Roll and Rattle–Roll methods were implemented to constrain all the bonds involving hydrogen atoms to their equilibrium values. 35The calculations were run using 16 processors in parallel at the Pittsburgh Supercomputer Center, Pittsburgh, PA, and required about 2.5 cpu hours per picosecond.
Distribution of Agents in the Membrane
Snapshots of the simulation system taken near the beginning and end of the MD simulations are displayed in figure 1. The figure shows that most of the HFE molecules have migrated significantly from their initial positions during the 2 nanosecond run. At low concentration, three HFE molecules moved toward the bilayer center. At high concentration, most nonimmobilizer molecules seem to partition along the lipid hydrocarbon acyl chains except two molecules, which migrated to the water region.
The average distributions of the HFE molecules in the lipid bilayer, at high concentrations, are reported in terms of average electron density profiles along the direction normal to the bilayer surface, Z (fig. 2A). These profiles describe the distributions of the HFE molecules in the bilayer (e.g. , as seen by x-ray studies) near the beginning and end of the simulation run. The results indicate that HFE molecules have a preference for the core of the bilayer, away from the interface. In contrast, similar simulations for halothane in DPPC, reported elsewhere, 22indicate their preference for the upper part of the lipid acyl chains (Fig. 2B).
Effects on the Lipid Bilayer
To monitor the changes of the lipid bilayer properties induced by the presence of HFE molecules we have calculated the membrane thickness and the area per lipid. The results reported in figure 3show an expansion of the former as measured by the increase of the d-spacing (the repeat distance of the bilayer, quantifying the lipid + water thickness as measured by x-rays or neutron scattering studies), and a decrease of the average area per lipid. These changes are more pronounced for the higher concentration, and clearly opposite to those induced by the presence of halothane. 21,22
The overall effect of HFE molecules on the membrane lipid packing can be deduced from changes in the deuterium NMR order parameters along the acyl chains. 37Experimentally, the NMR order parameters, SCD, are derived from the measured residual quadrupole splitting of deuterated lipids. 38They are usually expressed in terms of molecular order parameters Smol(SCD=−½Smol), which can be calculated from a MD trajectory using the identity: Smol=½< 3 cos2θ−1 >. Here, θ is the angle between the normal to the bilayer and the normal to the plane of the two C–D bonds of the lipid chain methylene groups, and the brackets denote a time average over a MD trajectory. In practice, we assume that the order parameters are not sensitive to the H/D substitution and therefore, compare directly the estimates from simulations (where the lipids were not deuterated) to experiment.
It is useful to consider relative changes in the order parameters. We define R , as the ratio of Smolfor the DMPC/HFE system and that of the pure lipid bilayer. Figure 4reports the calculated R values for low and high concentrations of HFE in DMPC, as a function of the carbon number along the acyl chain. In both cases, the effects of the nonimmobilizer are moderate, showing only a slight increase of the order parameters for the highest concentration. In contrast, results for halothane, reported on the same figure, show a drastic decrease of the order parameters when the molecules are present at high concentration in the membrane.
To investigate HFE effects on the membrane water interface, we report data on the average lipid head-group orientation. Specifically, we estimate from the simulations the probability distribution of the direction of the local dipole connecting the P and N atoms of the phosphatidylcholine head-group, with respect to the normal to the bilayer. The distributions of the P-N vector orientation for the systems with and without HFE are displayed in figure 5. The results show that the presence of the nonimmobilizer essentially has no effect on the head-group orientations. Similar plots for the DPPC–halothane systems are reported in figure 5B. They show that upon addition of the anesthetic, the P-N orientational distribution changes as the concentration of agent increases. It is worth noting here that the dipole distributions for both the pure DMPC and DPPC are not similar. This may be caused by the fact that these systems are simulated at different temperatures (T = 30°C and T = 50°C, respectively), which will significantly affect the hydrogen bonding capability of water and its interaction with the headgroups. These temperatures were chosen to yield the biologically relevant liquid crystal phases for each bilayer. Whether temperature is the main factor causing these changes in the head-group interfacial organizations of the lipid is currently under study. However, the present focus is on the relative changes of this organization when guest molecules are present in the system.
We have performed a series of MD simulations of a model membrane with low and high concentrations of hexafluroethane (HFE), a nonimmobilizer, and compared the results to those obtained previously for halothane, its anesthetic analogue. The results focus on the membrane response to the presence of the guest agents.
The membrane is modeled by a hydrated DMPC lipid bilayer, in its bio-relevant liquid crystal phase, at constant pressure and T = 30°C. Unlike biologic membranes, this system is composed of a single lipid, with saturated acyl-chains. The choice of the temperature and the hydration level ensures however that its characteristics are similar to those of the fluid membrane phase, despite the fact that T = 30°C used herein is below physiologic temperature.
The simulation protocols carried out here are similar to those previously used for the DPPC–halothane system. The initial HFE–model membrane system consisted of hydrated DMPC lipid bilayer, in which guest molecules were randomly incorporated. The calculations converged after 2 nanoseconds to equilibrium states that are shown to be in good agreement with most available experimental data.
The use of two concentrations, one around typical IA clinical levels and the other well above, is motivated by two reasons. First, most available experimental data on these systems concern high concentrations. Comparison of our results to this data are crucial since it ensures that the models and methodologies used are valid. Second, because of the limitations in time and length scales (i.e. , size of the system and number of agents considered), MD simulations containing only a few guest molecules are less reliable. In this study, the changes induced in the bilayer structural properties were monitored at low and high concentration HFE. The overall picture that emerges from the simulations is that the changes induced in the lipid bilayer, for both halothane and HFE are enhanced as the concentration of agent increases.
The main finding is that during the MD simulation run, HFE molecules migrate from their initial position to an equilibrium configuration where they distribute along the lipid acyl chains and in the methyl trough, with no apparent preference for the head-group region. This is in contrast to halothane, which was shown to prefer to populate regions near the upper part of the lipid membrane. This indicates that the calculation was long enough to allow the inserted molecules to respond to their immediate environment, and to partition in the membrane according to their affinity for the lipid–water interface vs. the lipid hydrophobic core. The present results agree well with 19F NMR experiments and give support to theoretical investigations based on energy profile calculations. 18–20The results are also in agreement with the interpretations of NMR study of anesthetics and nonimmobilizers at various concentrations, 17which attributed distinct locations to anesthetics and nonimmobilizer compounds in the membrane. As suggested by Pohorille et al., 18–20the affinity of HFE for the lipid acyl chain region can be interpreted in terms of a predominance of hydrophobic interactions with the membrane. The tendency of halothane to aggregate near the interface, exhibiting a nonuniform distribution across the membrane, suggests a tendency to maximize the electrostatic interactions with the interface at the expense of hydrophobic interactions with the lipid core.
The presence of guest agents in the bilayer modifies the overall membrane structure. Indeed, HFE molecules induce an expansion of the bilayer thickness and a lateral contraction manifested by a decrease of the area occupied per lipid. These changes are enhanced when passing from low to high concentration of nonimmobilizer and are opposite to those induced by the presence of halothane. Different effects on the lipid chain onformations are observed. Upon insertion of HFE, the lipid chain order parameters appear to change only slightly, while for halothane a significant decrease toward the end of the chains is observed.
The difference in the partitioning of HFE and halothane in the membrane not only affects the lipid core characteristics in opposite ways, but also influences differently the electrostatic properties of the membrane– water interface. Indeed, the results show that the distribution of the lipid head-group dipole (Phosphate-Choline, P–N pair), is not perturbed by the nonimmobilizer presence, even at high concentrations. In the case of halothane, our previous MD results show the opposite effect, namely a modification of the average orientation of the P-N dipole. The phosphatidylcholine lipid head group dipole orientation makes a significant contribution to the electrostatic potential across the membrane–water interface. The changes induced by the presence of halothane produce a measurable increase of this potential. From the lack of effect on the head groups on insertion of HFE, it is expected that the interfacial electrostatic properties remain relatively unchanged.
The effects on the membrane characteristics may all be attributed to the distinct distributions of these guest molecules within the lipid region. Indeed, aggregation of the molecules near the lipid head group region (anesthetic case), is expected to increase the head-group to head-group spacing, and therefore increase the free volume in the hydrophobic lipid core toward the end of the acyl chains. To reduce the excess of free volume, the acyl chains adopt new conformations that result in a decrease in the overall thickness of the membrane. On the contrary, when the guest molecules are homogeneously distributed along the lipid acyl chains (nonimmobilizer case), they are likely to make the membrane core more rigid, increasing, therefore, the order of the acyl chains, which results in an increase of the membrane thickness. The partitioning of the anesthetic near the lipid head-group influences also the organization of the zwitterionic phosphatidylcholine moiety. These changes do not occur when the guest molecules partition in the lipid hydrophobic core away from the interface.
The current simulation results show that membranes exhibit different responses to the introduction of guest anesthetic and nonimmobilizer molecules. While most of the underlined differences are more pronounced for the high concentration samples, i.e. , well above typical IA clinical concentrations, the tendency is also likely to be present for the low concentration systems. Such is the case for the effects on the membrane thickness and lateral lipids organization, and more importantly, for the effects on the interfacial electrostatic properties. Both predictions need to be confirmed by experimental measurements.
Whether the changes in membrane properties induced at low guest-molecules concentrations are enough to explain the distinct clinical behavior of the two agents remains unclear. Do anesthetics exert their in vivo effects by direct interaction with protein targets or indirectly by modifying the membrane properties? In light of experimental studies one cannot rule out a combination of the two. Our study does not address the possibility of direct binding of the agents to membrane proteins. However, the function of membrane proteins is known to be modulated by the properties of the supporting membrane, crucial among which are the membrane packing and the interfacial electrostatic potential. We have provided evidence that in this regard, the nonimmobilizers and the anesthetics induce opposite changes in these properties.
Historically, nonimmobilizers received attention because they were not anesthetic at concentrations where the Overton Meyer correlation said they should be—i.e. , there was the same amount in the membrane, but this did not cause anesthesia. This fact has been viewed as evidence against lipids being an important target. However, the present results clearly show that despite similar amounts dissolved in the membrane, some physical parameters are controlled by the molecular-level distribution. Moreover, since this distribution is different for halothane and HFE, and since membrane structural parameters could control protein activity, one cannot exclude the possibility that lipid membranes are an important direct target. The MD simulations demonstrate that mole fraction alone is not predictive of all lipid effects; essentially because lipid membranes are not isotropic solvents like olive oil. A systematic study of diverse model membranes and a range of guest molecules would be useful for a better understanding of this compelling problem.
The authors thank Professor Roderick Eckenhoff, Ph.D., and Assistant Professor Jonas S. Johansson, Ph.D. (both from the Department of Anesthesia, University of Pennsylvania, Philadelphia, PA), and Research Fellow Carlos F. Lopez, B.Sc. (Center for Molecular Modeling, Department of Chemistry, Philadelphia, PA) for contributions that have greatly improved this manuscript.