| HOME | HELP | FEEDBACK | SUBSCRIPTIONS | ARCHIVE | SEARCH | TABLE OF CONTENTS |
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Department of Physics, University of Illinois at Urbana-Champaign, and Beckman Institute for Advanced Science and Technology, Urbana, Illinois
Correspondence: Address reprint requests to Klaus Schulten, E-mail: kschulte{at}ks.uiuc.edu.
 |
ABSTRACT |
|---|
 TOP ABSTRACT INTRODUCTIONMETHODSRESULTS AND DISCUSSIONCONCLUSIONSAPPENDIXSUPPLEMENTARY MATERIALACKNOWLEDGEMENTSREFERENCES |
|---|
|
INTRODUCTION |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTS AND DISCUSSIONCONCLUSIONSAPPENDIXSUPPLEMENTARY MATERIALACKNOWLEDGEMENTSREFERENCES |
|---|
; Sukharev and Corey, 2004). These membrane proteins open and close upon application of mechanical stimuli, thus generating an ionic current that eventually becomes an electrical signal (mechanotransduction) or that simply helps, for instance, in the regulation of cell volume or the adaptation to osmotic pressures.
In bacteria, different mechanosensitive ion channels (MSCs) gated by mechanical strain in the cellular membrane have been identified (Martinac et al., 1987; Berrier et al., 1989; Sukharev et al., 1994). In particular, it has been shown that the mechanosensitive channel of large conductance, MscL, and the mechanosensitive channel of small conductance, MscS, both play an important role in osmoregulation and prevention of cell lysis in Escherichia coli (Levina et al., 1999; Blount and Moe, 1999; Sukharev, 2002).
MscL, with a conductance of 
3.2 nS and an activation midpoint of 11 dyn/cm, has been extensively studied (see Hamill and Martinac, 2001; Sukharev and Corey, 2004, and references therein) and its crystallographic structure was the first one solved for a mechanosensitive ion channel (Chang et al., 1998). Gating mechanisms for this nonselective MSC have been explored and modeled using the available closed structure as a starting point (Sukharev et al., 2001; Betanzos et al., 2002; Perozo et al., 2002a,b; Gullingsrud et al., 2001; Gullingsrud and Schulten, 2003; Kong et al., 2002; Colombo et al., 2003; Markin and Sachs, 2004).
On the other hand, MscS, with a conductance of 1 nS and an activation midpoint of 5 dyn/cm (Sukharev, 2002; Bezanilla and Perozo, 2002), is now becoming the subject of intense scrutiny due to its voltage sensitivity and the recent release of its x-ray crystallographic structure at 3.9 Å resolution (Bass et al., 2002), the first structure available for a voltage-modulated ion channel.
The overall architecture of MscS is shown in Fig. 1. The protein is formed by seven identical subunits of 286 amino acids each; residues 27–280 are present in the crystal structure. The transmembrane domain has three 
-helices per subunit as shown in Fig. 1 A. The first two transmembrane helices (TM1 and TM2) are tilted 27–35° with respect to the normal, and have charged residues (arginines) that may act as voltage sensors (Bass et al., 2002; Bezanilla and Perozo, 2002). The domain formed by these two transmembrane helices is similar to the controversial structure of the voltage sensor in the voltage-gated potassium channel KvAP (Jiang et al., 2003a,b; Perozo and Rees, 2003). The third transmembrane helix is divided by a pronounced kink into two subunits: TM3A and TM3B (Fig. 1 B). The division is observed at residue Gly113 and marks the end of the pore's cytoplasmic side. The transmembrane and highly hydrophobic pore, formed by TM3A, appears to be in an open state (Fig. 1 D) with a diameter of 6.5 Å in its narrowest section, located at the seven Leu105 side chains (Fig. 1 C). The second part of the third transmembrane helix, TM3B, is almost parallel to the membrane plane (Fig. 1 B), and the end of TM3B of one subunit slips below the loop between TM1 and TM2 of the next subunit (Fig. 1 D).
|
|
; B. Martinac, personal communication). In addition, MscS shows slight anionic preference and is voltage-modulated (Martinac et al., 1987; Sukharev, 2002). The estimated energy of opening is 
G = 11.4 ± 0.5 kT and the transition-related area change of the pore is A = 8.4 ± 0.4 nm2, values that have been obtained using a two-state Boltzmann model and curvature data from video imaging (Sukharev, 2002).
Experiments have also shown that small deletions in the C-terminal region (bottom end of the cytoplasmic domain shown in Fig. 1 A, B, and E) failed to produce functional channels. However, the extracellular N-terminal segment, denoted at the top of Fig. 1 B, seems to tolerate significant (residues 8–21) deletions (Miller et al., 2003a). Moreover, cross-linking of C termini suggested that they move apart upon channel opening and are involved in channel activation (Koprowski and Kubalski, 2003). In addition, cross-linking of site-specific cysteine substitutions suggest that residues in position 58 (situated at the loop between helices TM1 and TM2, Fig. 1 B) can be located within 3 Å of each other compared to being 33 Å apart as in the crystal structure (Miller et al., 2003b). Altogether, these experimental results suggest that the protein is very flexible, its C-terminus is structurally highly relevant, and there are significant structural differences between the closed and open states.
Random mutagenesis experiments have identified some functionally relevant residues. The double mutation I48D/S49P (both residues located at the cytoplasmic end of TM1, Fig. 1 B) inactivates the channel, whereas T93R, L109S, A102P, V40K, and V40D have been identified as gain-of-function phenotypes (Miller et al., 2003a; Okada et al., 2002). It is not surprising that residues at positions 93, 102, and 109 affect the conduction properties of the channel, since they are located on the pore-lining helices (TM3A in Fig. 1 B). The role of positions 40, 48, and 49 is less clear (all located throughout TM1, Fig. 1 B).
Simulations of membrane proteins are becoming a standard tool to gain insight into the mechanisms of ion channels (Roux, 2002; Ash et al., 2004; Roux and Schulten, 2004). Although most of the simulations have focused on ion or water permeation, increasing attention has been paid to the channel gating, especially in MscL (Elmore and Dougherty, 2003; Gullingsrud and Schulten, 2003; Colombo et al., 2003) and in the KcsA potassium channel (Biggin and Sansom, 2002).
A molecular dynamics study of water permeation and ion conduction in MscS (Anishkin and Sukharev, 2004) suggested that the crystal structure of MscS does not represent an open state and that MscS gating involves a vapor-lock mechanism. It is rather unclear how the simulation conditions utilized (protein restrained and simulation of only part of the protein) along with possible distortions in the structure due to crystallization conditions affected these results.
In this report, we present a molecular dynamics exploration of the entire MscS in a fully solvated lipid bilayer, with and without restraints. We describe permeation properties of the restrained structure of MscS and the spontaneous closure of the MscS transmembrane pore after release of restraints. In addition, we describe the effect of surface tension on the transmembrane channel as well as features of the membrane-protein interaction and the conduction properties of MscS.
|
METHODS |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTS AND DISCUSSIONCONCLUSIONSAPPENDIXSUPPLEMENTARY MATERIALACKNOWLEDGEMENTSREFERENCES |
|---|
Assembly of protein, lipids, and solvent
Using VMD (Humphrey et al., 1996) and its Orient plug-in, the MscS crystal structure (Bass et al., 2002; PDB code 1MXM) was oriented by aligning the protein's axis of symmetry along the pore parallel to the z axis. By using the VMD Membrane plug-in, the molecule was then immersed in a previously equilibrated rectangular patch of POPC membrane bilayer with solvated lipid headgroups. The initial membrane position along the z axis was determined by matching positions of headgroups with the seven Tyr27 side chains of the MscS transmembrane domain. The resulting placement of MscS in the membrane was in agreement with the suggestion in Bass et al. (2002) (Fig. 1 A). Lipids within 0.8 Å of the protein and those located in the pore were removed. The final number of lipids was 299, 159 on the periplasmic side and 140 on the cytoplasmic side.
The membrane-protein complex was then solvated using the VMD plug-in Solvate. Water molecules located in the hydrophobic part of the membrane and inside the protein's transmembrane pore were eliminated. The VMD plug-in Autoionize was used to place ions randomly, corresponding to a salt strength of 200 mM KCl. The resulting system had zero charge and 52,015 water molecules. Residues Asp, Glu, Lys, and Arg were assumed to be charged throughout the protein. A snapshot of the constructed system of a total of 224,340 atoms is shown in Fig. 2. The initial size of the system was 130 x 130 x 172 Å3 and shrank to 108 x 111 x 178 Å3 during the equilibration.
Molecular dynamics
All simulations were performed using the parallel molecular dynamics program NAMD 2.5 (Kalé et al., 1999) and the CHARMM27 force field for lipids and proteins (MacKerell et al., 1998a,b, 1992, Schlenkrich et al., 1996) along with the TIP3P model for water (Jorgensen et al., 1983).
A cutoff of 12 Å (switching function starting at 10 Å) for van der Waals interactions was assumed. Periodic boundary conditions were used with a flexible cell. Distances between periodic images of the protein were >17 Å in all cases. The particle-mesh Ewald method was used to compute long-range electrostatic forces without cutoff. The density of grid points for particle-mesh Ewald was at least 1/Å3.
An integration time step of 1 fs was used, permitting a multiple time-stepping algorithm (Grubmüller et al., 1991; Schlick et al., 1999) to be employed in which interactions involving covalent bonds were computed every time step, short-range nonbonded interactions were computed every two time steps, and long-range electrostatic forces were computed every four time steps. The use of a multiple time-stepping algorithm may affect the membrane pressure profile (Gullingsrud and Schulten, 2004). However, it allowed us to improve the performance of the simulations by a factor of 2, and thus reach the multinanosecond timescale required to observe structural relaxation of MscS.
Langevin dynamics was utilized to keep a constant temperature of T = 300 K in all simulations. The Langevin damping coefficient was set to 1 ps–1 unless otherwise stated. This value was chosen so as to speed up conformational changes of the channel by reducing the solvent viscosity. Constant pressure simulations at 1 atm were conducted using the hybrid Nosé-Hoover Langevin piston method with a decay period of 200 fs and a damping timescale of 50 fs. When surface tension was applied, components of the pressure tensor were set using the formula (Feller and Pastor, 1999; Gullingsrud et al., 2001)
 | (1) |
is the surface tension, and Lz is the size of the unit cell along the z axis, which is perpendicular to the plane of the membrane. The pressure along the z axis was kept at 1 atm and was set to 20 dyn/cm.
Analysis tools
To monitor the structural deformation of the protein, its root mean-square deviation (RMSD) was computed using VMD and coordinates saved every picosecond; the crystallographic structure served as the reference point, and only positions of protein backbone atoms were compared.
Pore radii were determined using HOLE (Smart et al., 1993). The HOLE program maximizes, using Monte Carlo simulated annealing, the radius of spheres located along and inside the pore defined by the van der Waals surface of the channel. The default all-atom set of van der Waals radii from AMBER (Weiner and Kollman, 1981) was used to allow comparisons with previously calculated radii profiles (Smart et al., 1997). Sample planes were spaced uniformly in the z direction with a 0.5 Å spacing. A note of caution should be added about the method utilized to measure the radius of the transmembrane pore as presented in Fig. 5. Our calculations give values that are smaller than stated in Bass et al. (2002) and subsequent publications, but are in agreement with values presented in Anishkin and Sukharev (2004). The disagreement could be due to the incorporation of hydrogen atoms and van der Waals radii of each atom in our calculations.
|
). For this purpose, anisotropies of the channel were taken into account (spherocylinders instead of spheres were used for optimization and estimation of radii) and a uniform separation between sample planes of 0.25 Å was used. In addition, 10 times more steps as normal in the HOLE Monte Carlo optimization were performed. HOLE estimates conductances using a formula adapted from Hille's equation (for details, see Smart et al., 1997),
 | (2) |
bulk is the resistivity of the permeant ions regarded as a function of the molar ion concentration c, s is the separation between the sample planes (s = 0.25 Å), A(z) is the area perpendicular to the channel direction vector, and the sum is performed over the sample planes along the z axis. HOLE also provides empirically corrected values of Gmacro using data from high resolution structures ("learning set") with known conductances (Gfirst), the net overall charge of the channel (Gsecond), and its length. Corrections to Gmacro based on the channel length were not taken into account because the length of the MscS transmembrane channel was not well represented in the "learning set" used in HOLE. The resistivity of 1 M KCl was assumed to be 0.08 
m (Smart et al., 1997) along with a linear dependence of on the ionic concentration. The final estimated conductance Gf was calculated as the average of Gfirst and Gsecond multiplied by 0.25 M, the molar concentration of KCl used in experiments (Martinac et al., 1987). The conductances summarized in Table 2 should be considered as estimates, correct only within a factor of 2.
|
|
Averaged ion density maps were generated using VMD and trajectories of simulations sim2a and sim2b (see below) saved every picosecond. The cell size was set to the maximum value reached during the respective simulation. The grid resolution was set to 1 Å in each direction of the cell.
System equilibration
We equilibrated the protein-lipid-solvent system for 3 ns. This simulation step will be referred to below as sim0. The equilibration resulted in a close packing of lipids against the protein, preventing a leakage of water molecules through the membrane-protein interface from the cytoplasmic to the extracellular side.
Simulation sim0 was divided into three substeps. In the first step, "melting" of lipid tails was performed. For this purpose, all atoms of the simulated system were frozen except those of the lipid tails, and 5000 energy minimization steps and 0.5 ns of dynamics under NVT ensemble conditions were carried out. At this point a Langevin damping coefficient of 5 ps–1 was used. In all subsequent calculations the system was simulated as an NpT ensemble.
Subsequently, we carried out step two of simulation sim0. During this step the protein backbone atoms were harmonically restrained with a force constant of k = 1 (Kcal/mol)/Å2. The step began with an energy minimization of 1000 steps and continued then with 1.42 ns of dynamics. Initially, a Langevin damping coefficient of 5 ps–1 was used, and during the last 0.3 ns it was reduced to a value of 1 ps–1, a value that was also assumed in all subsequent simulations. At 1.42 ns, the packing of lipids surrounding the protein was inspected and some lipids reordered by means of interactive molecular dynamics (Stone et al., 2001; Grayson et al., 2003). The system was then minimized for 1000 steps and molecular dynamics resumed for 0.373 ns. To prevent water molecules from reaching the protein-lipid interface, forces were applied that directed them back to the nearest water-lipid interface. These forces were switched on for intervals of 50 ps every 300 ps of dynamics in the second step of simulation sim0 only.
In the third and final step of simulation sim0, 0.7 ns of NpT ensemble dynamics were carried out without forces applied to the water molecules. A satisfactory packing of lipids resulted and no leakage of water molecules through the membrane-protein interface was detected, as stated already above.
After simulation sim0, four simulations (sim1, sim2a, sim2b, and sim3a) were carried out for a total of 10.35 ns. A summary of these simulations is shown in Fig. 3. In these simulations all components, i.e., all 224,340 atoms including protein, membrane, water molecules, and ions, were free to move, except in sim1, where the protein backbone was restrained. Indeed, the entire MscS was free to move in sim2a, sim2b, and sim3a, namely, both transmembrane and cytoplasmic domains that dynamically influenced each other. Such influence had been suggested by experiments and appears to be relevant for MscS stability and function (Miller et al., 2003a,b; Koprowski and Kubalski, 2003).
|
|
RESULTS AND DISCUSSION |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTS AND DISCUSSIONCONCLUSIONSAPPENDIXSUPPLEMENTARY MATERIALACKNOWLEDGEMENTSREFERENCES |
|---|
MscS with backbone restrained to the crystallographic structure
The first simulation (sim1) was designed to study MscS when restrained to the crystallographic structure. Indeed, in the case of this 1.45 ns simulation, the protein backbone was harmonically restrained with a force constant of k = 1 (Kcal/mol)/Å2 to the positions determined in the crystallographic structure. A snapshot of MscS at the end of this simulation is shown in Fig. 4, A and B. The pore was found to remain in the initially open state. This is concluded from the radius of the transmembrane pore measured as shown in Fig. 5 A. This figure presents the radius profile of the transmembrane pore of MscS. In the initial "crystal structure" the narrowest part of the transmembrane pore was found at the level of the seven Leu105 side chains (3.26 Å). During sim1, this part of the pore (neck) assumed a radius of 2.8 Å, a value that remained stable throughout sim1 (Supplementary Material, Fig. 20 A). In addition, a small displacement of the pore's neck toward the periplasm was observed (Fig. 5 A). Besides the channel opening (Fig. 5), we monitored in our simulations the pore constriction (Fig. 6), the RSMD values (Fig. 7), and water permeation (Fig. 8).
|
|
|
|
and studies that showed liquid-vapor oscillations of water in hydrophobic nanopores (Beckstein and Sansom, 2003, 2004). Formation of single files of water molecules was observed inside the bottom opening of the cytoplasmic domain of MscS during sim1. However, the size of the bottom opening was greatly reduced.
Spontaneous closure of the MscS pore
In a second simulation, sim2a (4 ns), we explored the behavior of MscS without restraints, starting from the final state of sim1. However, we did not release the restraints entirely at one instant, but rather released them gradually over 0.9 ns, choosing force constants of k = 0.5, 0.1, and 0.01 (Kcal/mol)/Å2, in steps of 300 ps. Subsequently, we carried out 3.1 ns of free dynamics.
A top view of the final state of this simulation is shown in Fig. 4 C. One can immediately recognize a dramatic motion of the protein that led to a complete closure of the transmembrane pore. The pore radius at the end of the simulation was 1.27 Å in its narrowest section, i.e., at Leu105. However, the value of the radius fluctuated between 0.5 and 2.4 Å (Supplementary Material, Fig. 20 B). The closing of the transmembrane pore was accompanied by a displacement of the pore's neck by 5 Å toward the periplasmic side (see Fig. 5 B).
RMSD values of MscS during sim2a increased to a level of 4 Å, as seen in Fig. 7 A. The RMSD values were found to differ for different zones of the protein as shown in Fig. 7 B. The smallest RMSD value was found for the pore zone, whereas the largest values were found for the cytoplasmic domain and for helices TM1–TM2. Fig. 7 A for the overall RMSD value and Fig. 7 B for the zone RMSD values show a different time-dependence, in that the total RMSD value increased throughout sim2a, whereas the zone RMSD values reached stable plateaus in the 2–3 Å range. This suggests that the increase in the overall RMSD value is due to relative motions of protein domains (Supplementary Material, Fig. 15), and not due to intrinsic local changes of protein structure.
The sudden change in the RMSD values of the pore region in sim2a (Fig. 7 B, between 6.5 and 7.5 ns) is due to a deformation observed in the pore-lining -helices and will be described in more detail below. The increase of the pore RMSD is trailing the increase of the RMSD values of TM1–TM2 and the cytoplasmic domain, suggesting a causal relationship.
A closer look at the dynamics of the pore in simulation sim2a revealed two interesting features. First, the bottom entrance of the transmembrane part of the pore, determined by the TM3B segment of TM3, moved up toward the periplasm (Supplementary Material, Fig. 19, sim2a), as already mentioned above. Second, the main feature of the pore closure was an asymmetric displacement of the seven Leu105 side chains along with the appearance of a kink in the -helical structure of TM3A in subunits A, C, and E (Fig. 6 and Supplementary Material, Fig. 19, sim2a).
A key deformation causing the closing of the pore is shown in Fig. 6 A. The TM3A helix of subunit C buckled shoving Leu105 toward the pore center and moving residues Gly101, Ala102, and Ala106 along with it (Fig. 6). The motion involved breaking and transferring of -helical hydrogen bonds. The helix TM3A is characterized by the sequence motif GxxGxxxGxxxxG, conserved in MscS homologs. This motif starts at residue Gly101 and ends at residue Gly113.
Two possible causes were identified for the deformation of pore lining structures: direct interaction with the lipid bilayer or tension transmitted through TM1 and TM2. As shown in Fig. 6 B, the deformation in subunit A (light blue) seems to be caused by the close proximity of lipids, whereas deformation in other subunits (C, red and E, orange) is not directly related to surrounding lipids (no contact is observed between lipids and TM3A in those cases). In both cases, the direct or indirect effect of the lipid bilayer tension on the overall structure of the protein could be overstated, since it is known that the force field parameters used in MD underestimate the correct area per lipid by 5% (Feller and Pastor, 1999; Gullingsrud and Schulten, 2004). This may explain the short timescale observed for the closure of the channel (Supplementary Material, Fig. 20 B).
The asymmetric movement of TM3A helices from different subunits seems to be compatible with the steric constraints experienced by seven helices that are in close proximity (Bezanilla and Perozo, 2002). However, distortions of helices resulted in a small lateral opening between TM3A helices of subunits A and B, compromising the integrity of the pore at the end of sim2a.
After release of harmonic restraints, the transmembrane zone of the pore remained empty and no conduction of water molecules could be discerned, as documented in Fig. 8 B. In addition, the distal opening of the cytoplasmic domain rapidly collapsed after release of restraints, and remained in a closed state during almost all of simulation sim2a (data not shown).
The closed state is favored under the chosen simulation conditions and should be the default state in a relaxed membrane (Bass et al., 2002). However, the spontaneous closure of the MscS transmembrane pore observed in simulation sim2a raised several questions that lead us to devise a third simulation, sim2b, in an attempt to keep the MscS pore open.
MscS assumes open state when surface tension is applied
The third simulation, sim2b (3.7 ns), applied surface tension to keep the MscS pore open, as was also done in experiments by Sukharev (2002). Coordinates of the system at 5.35 ns of sim2a were utilized as a starting point for sim2b (these coordinates correspond to the last step before complete release of restraints on the MscS backbone in sim2a, when the transmembrane pore is still in an open state; see Fig. 3). The simulation sim2b consisted of 0.4 ns of dynamics with the protein backbone harmonically restrained using a force constant of k = 0.01 (Kcal/mol)/Å2 and a surface tension of 20 dyn/cm. Then 3.3 ns of free dynamics were carried out with the same surface tension applied. A top view of the final state of sim2b is shown in Fig. 4 D. An increase in the transmembrane pore size can be easily recognized. Indeed, the pore radius at the end of sim2b is 3.75 Å at its narrowest part, at the level of residues Leu105 (Fig. 5 B). This value of the pore radius is even larger than the one observed in the original crystal structure (3.26 Å). Moreover, intermediate values of the pore radius at its narrowest section during sim2b were found to be 4 Å (Supplementary Material, Fig. 20 C). A minor displacement of the narrowest part of the pore toward the periplasm accompanied the increase in pore size, as shown in Fig. 5 B.
The surface tension applied in sim2b is approximately four times higher than what is needed to experimentally open MscS (Sukharev, 2002), and even higher than the rupture limit for a membrane. However, as mentioned before, the force-field parameters and simulation conditions (NpT) underestimate the experimental area per lipid (Feller and Pastor, 1999; Gullingsrud and Schulten, 2004), which may justify the large value used. Moreover, the high surface tension permitted us to induce the pore opening on the short timescale of molecular dynamics simulations. Simulations of MscL under surface tension have also utilized unphysiological tension values of 60 dyn/cm in vacuum (Gullingsrud et al., 2001) and 50–500 dyn/cm in a hydrated membrane environment (Colombo et al., 2003). In particular, opening of MscL in a hydrated membrane environment was observed only at 500 dyn/cm, which is more than an order-of-magnitude larger than the experimentally required magnitude. In none of our simulations did we observe membrane rupture.
The RMSD value for simulation sim2b is similar to that of sim2a, increasing to a level of 4 Å, as seen in Fig. 7 A. RMSD values for different zones of MscS reached rather stable plateaus at the end of sim2b, as shown in Fig. 7 C. However, RMSD values differ from each other for different zones of the protein. In particular, the RMSD value of the transmembrane helices TM1–TM2 is 1 Å larger than that of the pore, the cytoplasmic domain, and the transmembrane helix TM3A–TM3B. This can be explained as an effect of the applied surface tension and direct interaction of lipids with TM1 and TM2. As in sim2a, deformation of the TM3A–TM3B helix followed deformation of the TM1–TM2 transmembrane helices.
The in-plane protein areal expansion, A, has been estimated from experiments to be 8.4 nm2 (Sukharev, 2002). According to our simulations (sim2b and sim2a), the area in the x,y plane of the periodic cell varies from 11,629 Å2 to 12,491 Å2 (see Table 1), which represents a change in total area of 8.62 nm2, close to the experimental value; this change in area observed between the closed state of sim2a and the open state of MscS of sim2b includes the expansion of the membrane bilayer itself. Therefore, we have estimated the area occupied by the protein itself (see Table 1). The change in area occupied by the protein observed between the final states of simulations sim2a and sim2b was 2.81 nm2. The inferred area per lipid is in good agreement with results presented in Gullingsrud and Schulten (2004).
Pore dynamics in sim2b showed a rather different behavior than in sim2a. Although a movement of the pore's bottom edge toward the periplasmic side was observed, the pore structure tilted due to movements of the cytoplasmic domain. These movements are directly linked to the pore through TM3B. In addition, and due to the applied surface tension, the pore remained open at all times (Supplementary Material, Figs. 19, sim2b, and 20 C). A further distortion of the -helix in subunit G accompanied by the outward movement of its respective Leu105 side chain was observed (Supplementary Material, Fig. 19, sim2b.5 and sim2b.6) and is likely caused by direct interactions with the lipid bilayer (Fig. 6 B).
When surface tension was applied to keep the pore open, conduction of water molecules was enhanced and the transmembrane pore remained occupied by water molecules for almost the entire simulation sim2b (Fig. 8 C). On the other hand, the distal central opening of the cytoplasmic domain remained closed and conduction of water molecules through it was not observed.
MscS closed state under the influence of surface tension
The fourth and last simulation (sim3a) used as a starting point the last step of sim2a (see Fig. 3) and consisted of 1.2 ns of molecular dynamics with an applied surface tension of 20 dyn/cm. Our goal was to explore the effect of surface tension on the closed state of MscS obtained in sim2a.
The final state of sim3a is shown in Fig. 4 E. One can recognize that the surface tension led to a partial opening of the MscS transmembrane pore. As shown in the respective radius profile, the neck radius was 2.55 Å at the end of sim3a (Fig. 5); intermediate values for the neck radius of 2.8 Å were monitored (Supplementary Material, Fig. 20 D). The timescale of the MscS opening observed in sim3a is different and longer than the timescale of MscS closure observed in sim2a. Indeed, experiments suggest that the closure of MscS is faster than the opening (Sukharev, 2002).
Displacement of the pore's neck toward the cytoplasmic side along with a clear increase of radius in the periplasmic zone of the pore (Fig. 5 B) was also observed. The periplasmic opening seems to be caused by interaction of lipids with Tyr27 and with surrounding residues of some subunits. The effect of lipid-protein interaction on the periplasmic opening may be altered due to the presence of residues 1–26, missing in our model. Our finding is not in agreement (but not completely incompatible either) with the suggestion that the periplasmic side furnishes a pivotal point for the movement of TM1–TM2 in the gating mechanism (Bass et al., 2002; Bezanilla and Perozo, 2002). Since the behavior of different subunits is not conclusive, further simulations with smaller values of surface tension and longer simulation times are needed, along with a model for the periplasmic domain of MscS.
During sim3a, the RMSD value increased further as shown in Fig. 7 D. The RMSD values per zone also showed an increase, especially for the cytoplasmic domain and for the TM1–TM2 helices. It is important to note that the method used to apply surface tension induces a change in the area of the complete periodic cell. Of course, the water box was chosen big enough to absorb changes in MscS geometry; nevertheless, the generated flow of water molecules may exert a small perturbation on the cytoplasmic domain.
Finally, while trying to open the closed state in sim3a, conduction of water molecules was observed through the transmembrane segment of the MscS pore (Fig. 8 D). The shape of the pore permitted the formation of more than one file of water molecules. As in sim2a and sim2b, conduction of water molecules through the distal opening of the cytoplasmic domain was not observed.
Conduction of ions in MscS
Gating of MscS is monitored experimentally by ionic currents measured in patch-clamp experiments. Although presently molecular dynamics simulations cannot reach the timescale where natural conduction of ions arises, conductances can be estimated from the structure of the protein at different stages of the simulations using HOLE, as described in Methods. The resulting conductances of MscS are summarized in Table 2. The smallest value, 118.75 pS, is associated with the closed state at the end of sim2a (intermediates may have even smaller conductances); the largest value, 235 pS, is associated with the open state at the end of sim2b. The maximum value of the conductance predicted by HOLE is relatively small compared to the experimentally determined MscS conductance (1 nS). Thus, if our calculations are correct, it is likely that the open state of MscS corresponds actually to one of the intermediates arising in sim2b with a wider pore than that of the final state, or to an even wider state not observed on the timescale of our simulations. Clearly, further investigations are needed to address the conduction properties of MscS (see Roux et al., 2004, for a recent review on the subject of ion permeation). We note that overall ion conduction includes permeation of the cytoplasmic domain (Fig. 9).
|
|
 | (3) |
).
An exciting open problem is the role of the balloon-shaped cytoplasmic domain. In our simulations, spontaneous diffusion of ions through the cytoplasmic side openings of MscS was observed. For example, a Cl– ion transited from the bulk to the interior of the cytoplasmic domain in simulation sim2a. This event involved a series of interactions with conserved arginines (156, 184, 185, and 238) among other residues. The transit of a K+ ion from the inside of the cytoplasmic domain toward the bulk was also seen in sim2a. The positive ion interacted on its way with residues Glu134, Asp137, Asp159, and Glu187, among others. Of these, only Asp159 is conserved in more than 65% of the MscS homologs analyzed in Bass et al. (2002). The pathways of the ions are illustrated in Fig. 9.
Different degrees of openness were observed for the side openings of the cytoplasmic domain. In none of our simulations did the openings assume a completely closed state. Apparently, these side openings provide a more robust outlet and inlet than the distal opening of the cytoplasmic domain (Fig. 2).
Transmembrane helices of MscS and the gating mechanism
Upon insertion of the protein into the membrane bilayer, a strong interaction between the TM1–TM2 helices and the surrounding lipids became established. Indeed, inspection of the TM1–TM2 helices revealed an altered mobility after insertion. The tilt angle (with respect to the membrane normal) of helices TM1–TM2 was monitored during all simulations (data not shown) and revealed strong fluctuations. However, all subunits showed a decrease in their tilt angle of their TM1–TM2 domains (reductions of up to 8°) in sim2a, which seems to be in agreement with proposed movements associated with channel closure and the reorientation of -helical transmembrane domains (Bass et al., 2002; Bezanilla and Perozo, 2002). The tilt angle of helices TM3A and TM3B (residues 96–113 and 114–127) increased (up to 45° for TM3A and 90° for TM3B in some cases) in sim2a and sim2b.
Due to strong fluctuations and the flexibility of the protein, a clear pattern could not be established for the overall motion of the transmembrane helices on the timescale of our simulations. The TM1–TM2 domains of some subunits moved according to the suggested gating mechanism in Bass et al. (2002), as shown in Fig. 11. However, a concerted movement of all subunits in the same direction was not seen. In fact, some subunits moved in one direction whereas the others moved in the opposite direction (Fig. 11), likely due to the overall displacement of the protein in an irregular environment. The movement of the transmembrane helices was found to be caused by the interaction of lipid headgroups with charged and polar residues of the TM1–TM2 helices (Arg46, Arg54, Arg59, Arg74, Lys60, and Tyr75). The rearrangement of the membrane (discussed below) lead to a close interaction of the charged and polar residues of each transmembrane subunit of MscS with at least one lipid headgroup. Based on these results we infer that these interactions explain the movement of the transmembrane helices and are a key element in the gating mechanism. Voltage sensitivity may arise from a rearrangement of side chains of these residues as well as movement of the membrane headgroups and changes of the number of ions present in the protein-water-lipid interface. Further simulations need to be performed to address this point.
|
). Moreover, in four out of eight cases where the negatively charged residue Asp is not present in position 62, it is being replaced by residue Glu, also negatively charged. Similarly, in six out of nine cases where the positively charged residue Arg is not present in position 128, it is being replaced by Lys, also positively charged. These observations indicate a strong conservation of charge that will favor the formation of the respective salt bridge in similar structures. However, the specific role of these salt bridges in the gating mechanism is not clear. The movement of the TM1–TM2 domain could be limited by the presence of the salt bridge, reducing the sensitivity of the channel to membrane strain. On the other hand, the salt bridges provide direct links between the movement of the TM1–TM2 domains and the pore lining helices. Therefore, sensitivity of the channel to membrane strain could be increased due to their presence. Whether the protonation state for Asp62 and other residues changes throughout the opening or closing of the channel is a question that should be addressed in future studies of MscS gating. Preliminary simulations of MscS with all Asp and Glu residues protonated (low pH) have shown larger RMSD for the protein and inconclusive effects on gating of the pore.
|
|
2 Å of both the membrane and the protein occurred (Supplementary Material, Fig. 15 and Fig. 18 B). On the other hand, sim2b, in which a surface tension of 20 dyn/cm was applied, exhibited a fluctuating position of the protein (with respect to the membrane) that resulted in a small displacement of the protein toward the periplasm. In sim3a, the membrane and the protein moved together toward the cytoplasmic side. The relative motions described above suggest a faulty initial condition of the membrane, but may also be caused by the absence of the first 26 residues of MscS in our model.
The thickness of the membrane, computed as the average distance between phosphorus atoms of each layer, reached a value of 38 ± 1 Å in sim1, a value that remained without major variations in all subsequent simulations (in agreement with Gullingsrud and Schulten, 2004; data not shown). The lipid order parameter was computed and compared favorably to values presented in Gullingsrud and Schulten (2004), also showing good agreement with the experimentally determined values (Supplementary Material, Fig. 17).
Besides the relative movement of the protein and membrane under different external conditions, deformation of the membrane near the protein was observed in sim1 and subsequent simulations, as shown in Fig. 14. The arrangement of headgroups around the protein and its charged residues located in the protein-membrane interface revealed a dramatic thinning of the membrane bilayer close to the protein. As a result, side chains of charged residues reoriented and were not facing the hydrophobic zone of the membrane as suggested by the crystal structure, but rather interacted with membrane headgroups. This interaction is a key element in the gating mechanism and movement of transmembrane domains as discussed above.
|
|
CONCLUSIONS |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTS AND DISCUSSIONCONCLUSIONSAPPENDIXSUPPLEMENTARY MATERIALACKNOWLEDGEMENTSREFERENCES |
|---|
A solution to the problems described can be provided by molecular dynamics simulations, except that force fields may be error-prone and timescales covered are most often too brief. The shortcomings indicated can be overcome, e.g., when function may not depend on the accuracy of the force field and when external forces are applied that accelerate function without undue distortion. The latter path can be established only by trial and error. Unfortunately, MscS is a very large protein with over 2,000 amino acids. To study this protein in its lipid-membrane environment with a complex electrostatic force field requires one to simulate the dynamics of 200,000 atoms for 10 ns or longer, a considerable methodological challenge. The present study has met this challenge and described MscS accordingly in a series of scenarios that capture dynamics for MscS restrained to the crystal form, relaxation of the channel from the crystal form and its subsequent dynamics, as well as membrane-tension-induced channel behavior.
The simulations revealed key expected and unexpected characteristics of MscS that are likely essential for an eventual explanation of the mechanisms underlying the various functions of MscS. Being one of the earliest simulations of MscS, in fact the first that includes dynamics of the entire protein, the simulations mainly represent a fact-finding approach rather than one achieving an already reductionist picture of MscS. However, the behavior of MscS seen in the simulations proves to be thought-provoking and certainly will motivate further investigations.
The most dramatic feature observed in the simulations reported is the spontaneous closure of the MscS transmembrane pore when the crystallographic structure that exhibits an open pore was permitted to adjust itself to a membrane environment of 300 lipid and 52,000 water molecules. This is an extremely fortunate result since it is the first simulation of a spontaneous gating motion of a mechanosensitive channel. Undoubtedly, it will provide a sound basis for studying the mechanical coupling between protein and environment. The simulations above have taken the natural step of attempting to reopen the already closed channel and of attempting to prevent the closing by applying membrane tension, the latter proving highly successful. The atomic-level description permitted a detailed view of interactions and geometric transformations governing pore closing and opening. The role of charged and polar residues in the cytoplasmic side of the transmembrane domain was found to be extremely important. The membrane surrounding the protein adjusts and deforms, establishing close interactions with charged side chains, which were found to face lipid headgroups instead of the hydrophobic zone of the membrane as suggested by the crystallographic structure. Support for the suggested gating mechanism comes from the observation of TM1–TM2 helices moving in a hinge-like fashion, with the TM1–TM2 loop "swinging" and the TM2–TM3 loop acting as the hinge. Whether different units move in a concerted manner is not clear. However, the flexibility of these domains provides an efficient way for the protein to adapt to an irregular and changing membrane environment.
Beside its transmembrane pore, MscS sports a large, balloon-like chamber that acts most likely as a filter in front of the cytoplasmic entrance of the channel. This domain exhibited a dramatic role during the simulations. First, and most surprisingly, a central opening distal from the transmembrane channel closed quickly and remained closed during the simulations that kept the domain surrounded by its native aqueous environment. Second, the domain revealed, in keeping with its heptameric structure, seven lateral openings that experienced flowthrough of ions, the flow assisted through polar and charged amino acids. The openings remained large enough for ion passage throughout the simulations. Third and last, the arrangement of ions inside and around this cytoplasmic domain appeared to form a distinct pattern that should be relevant for its function.
The spatial rearrangements of MscS observed in our simulations suggest that this protein is very flexible, although most of its secondary structure elements remained intact (see Supplementary Material, Fig. 16). For example, the dynamics often broke the sevenfold symmetry of the crystallographic structure and results suggest that hypothetical mechanisms that assume symmetric reaction pathways for channel gating are ill-advised. What nature seems to hold in store is a protein that is very flexible and at the same time robust in regard to function. To achieve this with a multimeric system seems to be rather a miracle that definitely deserves much scrutiny.
Perhaps the most intriguing suggestion coming out of the simulations are interactions between the cytoplasmic domain and the transmembrane domain of MscS. Particularly interesting is the formation of salt bridges between the side chains located at the end of TM3B of one subunit and residues located in the loop between the TM1 and TM2 helices of the subsequent subunit. The connection between TM3B and the TM1–TM2 domains directly links the pore-lining helices to the movement of the transmembrane and cytoplasmic domains. It is unclear if this link immobilizes the TM1–TM2 domain or if it helps to open the pore, i.e., assists in gating.
We need to mention also key shortcomings of our simulations. First, missing residues in the crystal structure at the periplasmic zone may affect the movement of TM1–TM2 helices and the membrane position, whereas missing residues at the very end of the C-termini may stabilize the distal opening of the cytoplasmic domain. Second, the timescale of our simulations, the size of the membrane patch, and the symmetry of the periodic cell used limit the comparison of our results with experiments and add uncertainties. Third, changes in protonation states of Asp and Glu residues for different conformations of MscS were not yet considered.
Many questions still remain unanswered. Is the position of the membrane and its deformation near the protein an artifact of the initial conditions utilized? Is the small predicted conductance of the obtained open state a sign of an intermediate state? How is the pressure profile in the membrane modified by tension and deformation and how is it related to the interaction with charged and polar residues of the MscS transmembrane domain? Is this flexible protein affected directly by hydrostatic pressure and water flow? Why is there a concentration of negative ions in the periplasmic mouth of the channel and what is the relevance of the bimodal arrangement of ions inside the cytoplasmic domain of MscS? What is the effect of the membrane potential and how is it related to the protein-lipid interaction? What is the role of the observed salt bridges and how can the inactivation of the channel at low pH be explained?
A particular challenge to future work appears to be a re-engineering of MscS that eliminates much of the cytoplasmic domain. A small, re-engineered MscS may reveal much more readily its wide capacity of gating. Another challenge is to explain the role of the amazing cytoplasmic domain.
atom of residue Arg128 of one subunit and the C