Molecular Dynamic Simulation of the Kv1.2 ... - ACS Publications

The head groups of lipid molecules have been pulled down and come near R2, which was found in Kv1.2 open state simulation work(16) and KvAP channel ...
0 downloads 0 Views 4MB Size
16966

J. Phys. Chem. B 2008, 112, 16966–16974

Molecular Dynamic Simulation of the Kv1.2 Voltage-Gated Potassium Channel in Open and Closed State Conformations Ming Han† and John Z. H. Zhang*,‡,§ Institute of Theoretical and Computational Chemistry, Key Laboratory of Mesoscopic Chemistry of MOE, School of Chemistry and Chemical Engineering, Nanjing UniVersity, Nanjing 210093, China, State Key Laboratory of Precision Spectroscopy, Department of Physics, East China Normal UniVersity, Shanghai 200062, China, and Department of Chemistry, New York UniVersity, New York, New York 10003 ReceiVed: September 5, 2008; ReVised Manuscript ReceiVed: October 20, 2008

We performed 45 ns atomistic force field based molecular dynamics (MD) simulations to explore the structure and dynamics of the Kv1.2 voltage-dependent potassium ion channel in both open and closed state conformations. The Kv1.2 crystal structure (PDB 2A79) based open state model and homologically derived closed state model are embedded in a hydrated lipid membrane. We present a detailed analysis of the protein stability, the environment of the gating-charge-carrying residues, the salt bridge within the voltage sensor, S4 segment helix movement, helix interaction between S4-S5 linker and S6, interaction between subunits, and interaction between protein and lipid membrane. The four subunits of the channel lost the symmetry from the starting structure during the simulation, especially the voltage sensors. According to our comparative analysis of the open and closed state conformations, S4 in our simulation behaves more near the screw helix model except for the tilt action. The S4 segment azimuthal rotation angle difference between two conformations shifts about 60° after being embedded in a membrane environment. The S4-S5 linker helix remains stable. The contact area between the S4 and S5 from the adjacent subunit increases after transition to the closed state from the open state. The current investigation provided valuable information to our understanding of the structure and dynamics of membrane-associated helices and the gating mechanism of the voltage-dependent potassium ion channel of the Kv channel. 1. Introduction Voltage-gated potassium (Kv) channels detect the changes in membrane potential and respond with conformational changes that lead to opening or closing of the K+ conductance pathway across the membrane.1 Kv channels are tetramers with four identical subunits, each containing six transmembrane segments, termed S1-S6. S5 and S6 form the central pore, and S1-S4 form the voltage sensor that is located around the pore. In particular, S4 has highly conserved arginine gating charges that have been identified as being responsible for coupling changes in membrane voltage to opening and closing of the pore.2,3 During depolarization, the channel moves the gating-chargecarrying residues nearer to the extracellular side of the cell membrane by electrostatic attraction. When the membrane is repolarized, gating-charge-carrying residues will go back. Kv channels play a key role in electrical signaling in the nervous system. There has been considerable progress in studies of Kv channels by both experimental1,4,5 and computational methods.6-10 However, many fundamental questions remain unanswered: How are the gating charges stabilized within the membrane? How far do the gating charges move? How is voltage sensor motion coupled to the motion of the channel gating? How are cooperative interactions between subunits mediated, and how does the protein channel interact with membrane lipids? The * To whom correspondence should be addressed. E-mail: john.zhang@ nyu.edu. † Nanjing University. ‡ East China Normal University. § New York University.

purpose of this study is to try to give some answers through atomistic molecular dynamics (MD) simulations. MD simulations have become a popular and powerful technique to study membrane proteins and lipids. MD simulations provide a computational approach to exploring structure/ function relationships of membrane proteins and the interactions of membrane proteins with their lipid bilayer environment.9-12 X-ray crystallographic studies of several Kv channels have revealed the molecular architectures,13,14 which provide computational basis for MD study. Kv1.2 is the first solved eukaryotic Kv channel, and there are already MD studies reported.15,16 However, the initial models in these works only include the open state conformation. From a physiological point of view, channel proteins carry out tasks that require significant conformational transition between the closed and open conformation. The transition is of the order of milliseconds, and thus cannot be covered in classical atomistic MD simulations. In our work, we used the available homology Kv1.2 closed state model17 and X-ray crystallographic structure for the open state conformation to further computational study of the sensing and gating mechanism of Kv1.2 by comparing the dynamical properties of the open and closed conformations. In the present study, we performed 45 ns atomistic MD simulations to explore the conformation dynamics of open and closed state conformations of the Kv1.2 channel in hydrated dipalmitoylphosphatidylcholine (DPPC) lipid membrane. By analysis of the protein stability, gating charge environment, salt bridges within the voltage sensor, S4 tilt and rotation, interaction between voltage sensor and pore domain, interaction between subunits, and interaction between protein and its environment, we hope to provide valuable clues to answer some of the

10.1021/jp807905p CCC: $40.75  2008 American Chemical Society Published on Web 12/04/2008

Kv1.2 Voltage-Gated Potassium Channel

J. Phys. Chem. B, Vol. 112, No. 51, 2008 16967

TABLE 1: System Composition and Box Vector of the Open and Closed State Systems conformation

no. of lipids

no. of lipids in U layer

no. of lipids in D layer

no. of water molecules

no. of Cl-

no. of K+

no. of atoms

box size (nm3)

open closed

522 528

245 255

277 273

28823 25758

40 40

40 40

123181 114286

14.7 × 14.7 × 7.7 14.4 × 14.3 × 7.6

questions mentioned above and to further the investigation of the nature of the voltage-gating mechanism. 2. Method The initial models come from the solved crystal structure of Kv1.2 (2A79)14 based open state structure as well as the homologically derived structure of the closed state conformation.17 We focus our study on the part from Ala162 to Thr421 of the channel provided in ref 17. All of the ionizable residues are in their ionization states expected at pH 7 based on the standard pKa values. Both the open and closed state models were inserted into the pre-equilibrated DPPC bilayer by aligning the protein’s axis of symmetry with the bilayer normal according to a protocol mentioned in ref 18. The solvation and MD simulations were run using GROMACS package v.3.3.1.19,20 K+ ions were initially placed at sites S0, S2, and S4 in both the open and closed state models, where they were suggested to be the most stable positions by the molecular dynamics free energy perturbation technique.21 Other potassium and chloride ions were added to the system by random replacement of water molecules, until a final concentration of about 0.1 M was achieved. The number of lipid molecules, ions, water, and total atoms in each simulated system and the size of the box are listed in Table 1. We select the Gromos force field with Berger’s lipid parameters.22 The simple point charge (SPC)23 water model was used. The calculation of electrostatic forces utilized the PME implementation of the Ewald summation method.24 A cutoff of 10 Å was used for the van der Waals interactions. Energy minimizations were performed employing a steepest descent algorithm following the setup to relax any steric conflicts generated during system setup. The system was then subjected to a 1 ns restrained MD run, during which the nonhydrogen protein atoms and K+ ions in the filter were restrained with a force constant of 1000 kJ mol-1 nm-2, to enable the water and lipid bilayer to pack properly around the protein. The simulations were run in the NPT ensemble. In the first 1 ns unrestrained run, both temperature and pressure were controlled by the

Figure 1. RMSD of CR of the whole protein, voltage sensor, and pore (black, red, and green lines, respectively) of open and closed state conformations from the starting structures.

Berendsen method.25 In the final 45 ns production run, the temperatures of the protein, lipid, and solvent (waters and ions) were separately coupled using the Nose´-Hoover thermostat26,27 at 323 K, with a coupling constant of Γt ) 0.1 ps. System pressures were semi-isotropically coupled using the ParrinelloRahman barostat28 at 1 bar, via a coupling constant of Γp ) 1 ps and compressibility of 4.5 × 10-5 bar-1. The bonds of nonwater molecules and water molecules are constrained by using the fast LINCS29 and SETTLE30 algorithms, respectively. The integration time step was 2 fs, and the coordinates and velocities were saved every 4 ps for subsequent analysis. No electric field was added in this system during this simulation. Simulations were performed on a 64-node PC cluster with dual Intel Xeon5150 processors. Secondary structure is analyzed by DSSP.31 Helix tilt angle and azimuthal rotation angle are analyzed by using the definition and calculation method in ref 32. Solvent-excluded surface (SES) is calculated by MSMS software.33 Other data are analyzed by GROMACS or locally written codes. Structural diagrams were prepared using VMD software34 developed at UIUC. 3. Result and Discussion The potassium ion channel needs one or several milliseconds to transit between closed state and open state. This time span is too long for atomistic MD simulation. In this paper, we concentrate on the comparison of the structure and conformational dynamics between open and closed states. The recently solved crystal structure of Kv1.2 in the open state conformation14 and the availability of the closed state conformation by the homology method17 provide us a chance to carry out MD simulations to investigate the gating mechanism of the voltagegated channel at atomic levels of detail. 3.1. Structure Stability. The equilibration of the systems was monitored through the potential energy, CR root-meansquare deviation (RMSD) of the atom coordinate, and proteinlipid interaction. The potential energy of both open and closed conformations was stabilized after 15 ns. The RMSD from the starting structure as a function of time was analyzed to assess the degree of conformational drift, as shown in Figure 1. The

Figure 2. Number of H-bonds between the Kv1.2 protein and lipid molecules as a function of simulation time.

16968 J. Phys. Chem. B, Vol. 112, No. 51, 2008 RMSD of the whole protein of the open state conformation jumps at about 7.4 ns, slowly increases to 0.32 nm within the first 22 ns, and then fluctuates within 0.28 and 0.34 nm. The initial jump in RMSD is presumed to reflect the relaxation of the protein upon transfer from a crystal to a bilayer environment. For the closed state, the RMSD fluctuates between 0.26 and 0.31 nm after 21 ns. The RMSDs of the pore are relatively small, indicating little conformational drift, while the voltage sensor domain has a significantly higher fluctuation than the whole protein and the pore. The voltage sensor is essentially an independent domain inside the membrane. It is a highly mobile domain to make it possible to carry charged amino acids easily through the membrane electric field. The RMSD of every subunit differs between each other, and all final RMSDs fluctuate within 0.22-0.39 nm in both open and closed conformations after reaching equilibrium. In Jogini and Roux’s work,16 the open conformation state with shorter loops was studied; the RMSDs of the channel backbone and pore domain relative to the initial structure are 4.5 Å and less than 2.0 Å, respectively, and lateral whole-body motions of the sensors were observed. The bigger RMSDs compared to our work here perhaps result from the X-ray solved structure, while our models were optimized in ref 17. However, in our simulation, we did not find the obvious lateral motions of the transmembrane part of the sensors. The RMSD of the closed state has a decreasing trend near the end of the simulation. We extended our simulation for another 5 ns and found it still fluctuates around the plateau value. Thus, proteins in two simulation systems are stable after 22 and 21 ns, respectively. We check the interaction between the protein and the lipid bilayer by counting the number of H-bonds. It needs more than 30 ns to stabilize the protein-lipid interaction, as shown in Figure 2. The whole simulation system of channel protein embedded in lipid membrane needs a long time to equilibrate, which agrees with Deol et al.’s computational analysis that simulations of about 20 ns duration do not fully sample protein and lipid interactions.12 We selected the final 10 ns configurations for our following static property analysis. The CR root-mean-square fluctuations (RMSFs) were calculated to check the fluctuations in the structure as a function of the residue number. Four subunits are plotted separately in Figure 3. The peaks in the RMSF values are observed at the N-terminal, at the C-terminal, and in the extracellular and intracellular surface loops. The voltage sensor fluctuates more than the pore domain, which agrees with the RSMD results. The mean fluctuation of the core transmembrane region of the pore shows a low mobility with a value of about 0.05 nm. The large loop between the S1 helix and the S2 helix on the extracellular side fluctuates most markedly because a short R helix within this loop in the solvent environment separates away from the rest of the protein. It was confirmed that the removal of the short R helix within the S1-S2 loop has little effect on voltage-dependent opening.35 Each subunit has a similar overall pattern, with the fluctuations ranging from 0.05 to 0.23 nm for the open state and from 0.04 to 0.25 nm for the closed state. However, the S2-S3 loop, S3, S4, and S3-S4 loop in the open state conformation have greater fluctuations than in the closed state, especially the extracellular half-part of S3 and S4. This is because S3 and S4 have closer interaction with the lipid headgroup in the open state while they are nearer to the core of the membrane and the pore domain in the closed state. In the closed state, the transmembrane core of the S4 helix fluctuates in a similar range as the core transmembrane region of the pore

Han and Zhang

Figure 3. RMS fluctuation as a function of residue number. Four subunits are plotted separately, and each curve is shifted from the previous one by 0.2 nm in the vertical axis for clarity. The transmembrane region is indicated by shadows.

except the subunit B. In particular, the S4-S5 loop shows very stable behavior in both open and closed conformations. The fluctuations on the secondary structure can be determined and visualized according to DSSP.31 The overall secondary structure pattern of Kv1.2 is maintained, as shown in Figure 4. The transmembrane region except for S3 and S4 is relatively stable; especially the secondary structure of the transmembrane region of the pore domain remains almost completely R-helical except for the PVP turn in both open and closed state conformations. Although the S4-S5 linker is solvent-exposed, it keeps a remarkably stable helical structure. The differences between the four subunits and between the closed and open state conformations are mostly located in the solvent-exposed loops on the extracellular and intracellular side, which do not have a definite secondary structure. S3 and S4 show big fluctuations. Only in subunit B of the open state, the Pro265 residue demarcates S3 into two helices as in KvAP.13 In other subunits and in the closed state, S3 keeps one helix. In particular, in the open state conformation, the intracellular half of the S4 helix adopts a 310-helix pattern, and remains stable in subunits A and D. In the closed conformation, it keeps the initial R-helix in two of the subunits and loses partial helicity in other subunits. The transition between the R-helix and 310-helix is presumed to have a “concertina effect” to stretch S4 without necessitating rigid body rotation of the S4.35 3.2. Voltage Sensor Structure and S4 Movement. The four subunits lost the symmetry of the starting structure during simulation, especially the voltage sensors, as shown in Figure 5. Four antiparallel transmembrane helices, S1-S4, pack

Kv1.2 Voltage-Gated Potassium Channel

J. Phys. Chem. B, Vol. 112, No. 51, 2008 16969

Figure 5. Top view of the Kv1.2 channel. Loops are not shown except the S3-S4 loop and S4-S5 linker. Blue, red, green, and purple are for subunits A, B, C, and D respectively. The left is the open state conformation, and the right is the closed state conformation.

Figure 4. Secondary structure as a function of simulation time calculated by DSSP. The subunits A, B, C, and D are listed from the bottom to top, respectively. The residues in each subunit run from Ala162 to Thr421.

counterclockwise to form the voltage sensor. It is in charge of detecting the transmembrane potential change to make gating voltage dependent by charge-carrying residues.36,37 In Kv1.2, six gating-charge-carrying residues including five arginine residues, R294, R297, R300, R303, and R309, and one lysine residue, K306, termed R1, R2, R3, R4, K5, and R6, are located at every third position along S4, as shown in Figure 6. The environment of the gating charges and the movement of S4 are the two major controversial issues in the literature.5 By atomistic

level simulation, we can see the relative distribution of the residues, solvent, and lipids clearly. By comparison of the two end conformations in the dynamic hydrated lipid environment, we can further study the movement of S4 and the gating mechanism. In the open state conformation of this simulation, R1 is located in a lipid environment, as shown in Figure 7a and b. It is near the membrane surface and extends to the phosphorus atom layer of the lipid bilayer. R2 is located in a relatively complicated environment. It is deeper and near the domain surrounded by S1, S2, S3, and S4. In subunit A, R2 forms a stable H-bond with the DPPC headgroup and glutamate E183 (S1 helix), as shown in Figure 6; in subunit B, R2 forms a transient H-bond with DPPC; in subunit C, R2 forms stable H-bonds with DPPC after 25 ns; in subunit D, R2 forms transient H-bonds with DPPC and E183. The head groups of lipid molecules have been pulled down and come near R2, which was found in Kv1.2 open state simulation work16 and KvAP channel simulation work too.38 R3, R4, and K5 are located within the domain surrounded by S1, S2, S3, and S4 and face helices S1 and S2, where they have salt bridge interactions with acidic amino acids. A salt bridge between R3 and E183 remains stable in four subunits. Another salt bridge R3-E226(S2 helix) keeps stable in subunits A and B during the simulation. In subunit C, R3 has stable H-bonds with DPPC and the lipid headgroup is pulled down. R4 and K5 are located in a protein environment. R4 forms stable H-bonds with E226 in four subunits. The H-bond between R4 and E183 is stable only in subunit C. K5 forms stable salt bridges with E236 (S2 helix), and with aspartate D259 (S3 helix) in four subunits. R6 is located in a lipid membrane environment, as shown in Figure 7a and b. In addition to the ionized hydrogen bonds between basic and acidic residues, positively charged amino acids form H-bonds with polarized amino acids. R4 forms stable H-bonds with Ser176 (S1 helix) in four subunits. The difference in salt bridge formation in four subunits results from the easy mobility of the voltage sensor, the long side chains of the gating-charge-carrying residues, and the complicated environment. In the closed state conformation, R2 is totally buried by the domain surrounded by S1, S2, S3, and S4. R1 is inside this domain but very near the pore domain. R1 and R2 are far from lipids that can be seen in Figure 7d and e. R1 and E226 form stable H-bonds in four subunits. In subunit D, R1 and E183 form a salt bridge. R2 and D259 form stable H-bonds in four subunits. In subunit B, R2 forms a stable H-bond with E226. In subunits C and D, R2 forms a H-bond with E236. R3 is exposed to intracellular aqueous solution and has stable interac-

16970 J. Phys. Chem. B, Vol. 112, No. 51, 2008

Han and Zhang

Figure 6. Distribution of gating-charge-carrying residues and the negative charged residues within the voltage sensor in subunit A at 45 ns. The left is the transmembrane part of the voltage sensor in the open state conformation, and the right is in the closed state conformation.

Figure 7. Cumulative number of water, lipid head groups, and methylene groups around the gating-charge-carrying residues as a function of distance from the residue mass center.

tion with lipid molecules. The salt bridge between R3 and D259 remains stable in subunit A and transient in subunits B and C. R4, K5, and R6 are located near the intracellular membrane surface. There are two water filled crevices on both sides of the membrane. After simulation begins, the water from both sides of the membrane penetrates into the two crevices. The number of water molecules can be seen in Figure 7c and f. In the open state conformation, R1, R2, R3, and R4 interact with water in an external aqueous crevice and K5 interacts with water in an internal aqueous crevice. In the closed state conformation, R2 separates the internal and external water crevices, which is different from the initial homology model17 of R1 at the bottom of the extracellular water-accessible vestibule. There are three major models of voltage sensing to describe the movement of the voltage sensor: transporter model, helical screw model, and paddle model according to the different interpretations about the environment of charge-carrying residues, the extent of S4 movement, and the nature of S4

movement.5 We calculated the main properties of S4 movement from the closed to open state by comparison of the closed and open state conformations. The vertical distance of the CR of the gating-charge-carrying residues, the tilt angle change, and the azimuthal rotation angle change of the backbone of S4 from the closed to open state along the simulation time are shown in Figures 8-10, respectively. The vertical movements of the CR of R1, R2, R3, R4, K5, and R6 after 35 ns are 0.93, 1.14, 1.15, 0.78, 0.33, and -0.08 nm, respectively. The result of the protein models in ref 17 is 0.8, 1.05, 1.04, 0.74, 0.19, and -0.39 nm. The big change after the model in ref 17 is embedded into a hydrated membrane environment is R6 because it is attracted by the lipid membrane and dragged near the membrane water interface in the open state system, as shown in Figure 7. In our simulation, the center of the S4 segment moves about 0.672 nm, which is close to the averaged residue CR movement of 0.65 nm in ref 17. The tilt angle and azimuthal rotation are calculated by the method provided in ref 32. The

Kv1.2 Voltage-Gated Potassium Channel

Figure 8. Vertical movement of the gating-charge-carrying residues from closed to open state conformation as a function of simulation time.

Figure 9. Change of the tilt angle change from closed to open state conformation as a function of simulation time.

Figure 10. Azimuthal rotation angle change from closed to open state as a function of simulation time.

tilt angle (the angle between the axis of S4 and the normal of the membrane) difference between the closed and open state in the initial models is 20°. After a big fluctuation, the tilt angle difference keeps stable around 25° after 24 ns. Although the averaged tilt angle 62° in the closed state is larger than 37° in the open state, S4 has less part to extend to the lipid environment because the S4 in the closed state translates closer to the pore domain. When we describe the azimuthal rotation, S4 tilt is assumed to take place in the same plane comprised of the normal of the membrane and the axis of the S4. The azimuthal rotation angle (the angle

J. Phys. Chem. B, Vol. 112, No. 51, 2008 16971

Figure 11. Change in contact area between the S4-S5 linker and S6 as a function of simulation time.

Figure 12. Fluctuation in the contact area of neighbor subunits as a function of simulation time.

given by the direction of the tilt with respect to a vector orthogonal to the helix axis that passes through the CR of R1) difference between the closed and open state decreases from ∼175 to ∼120°, as shown in Figure 10. We selected R2 as a reference residue, and the final stable value is 120° too. The S4 azimuthal rotation angle shifts 20 and 40° in the opposite direction in the open and closed states, respectively. Thus, hydrated membrane in the simulation system does affect the protein models. From the analysis above, the Kv1.2 channel in our simulation is more near the helical screw model. It agrees with the main points of the helical screw model;39-41 that is, the S4 helix does translation and rotation together to respond to depolarization. The difference is that tilt of the S4 helix does not happen in the helix screw model, while in our simulations the axis of the S4 helix tilts about 25° between the two end state conformations. In the paddle model,13,42,43 the S4 helix does tilt a big angle. However, the rigid paddle unit and as big as 15-20 Å vertical movement of the S4 helix are not found in our simulation. Does every Kv channel have the exact same voltage sensing mechanism, or is it different according to the different chemical structure? This needs more real and computational experiments to explore. 3.3. Interaction of Helices during Gating. Conformational changes in the voltage sensor domain are transferred to the ion conduction pore via the helical linker, and result in opening or closing of the intracellular gate of the ion conduction pathway. The S4-S5 linker and intracellular end of S6 in the same subunit are important for coupling voltage

16972 J. Phys. Chem. B, Vol. 112, No. 51, 2008

Han and Zhang

Figure 13. Top view of the Kv1.2 channel and the nearby lipids (lipids that have atoms within the circle of radius 4 nm from the center of the pore in the X-Y plane: (a) in an open system; (b) in a closed system. The hydrogen atoms are depicted in white, carbon atoms cyan in protein and yellow in lipids, nitrogen atoms blue, and oxygen atoms red. (c and d) Two-dimensional number-density maps of DPPC lipid in the open and closed state systems, respectively.

Figure 14. Order parameter of DPPC tails. The two tails of the lipid molecule are calculated separately.

sensing to gating.5 The S4-S5 linker is an amphipathic R-helix and keeps parallel to the membrane plane in the interface of the membrane and aqueous solution. The

hydrophobic residues are located in the side facing the lipid membrane, while hydrophilic residues are located in the side facing water solution,44 which makes this helix stable, as shown in DSSP analysis in Figure 4. There is a Pro405Val406-Pro407 motif in the S6 helix, which forms a curve in this helix to make the intracellular half-part of S6 reach close to the S4-S5 linker. To quantitatively measure the helix contact, we calculated contact areas, which is an effective way to study the helix packing in membrane protein. The contact area measured by solvent-accessible surface (SAS) is shown in Figure 11. The solvent-accessible surface is traced out by the probe sphere center as it rolls over the helices. To calculate solvent-accessible surface, we have used the double cubic lattice method (DCLM)45 which is implemented in GROMACS. The contact area between the S4-S5 linker and S6 is ∼6.0 nm2 in both the open and closed state conformations, although contact area largely fluctuates in the closed state conformation. We compared the contact area of the final frame at 45 ns by the solvent-excluded surface (SES) by MSMS software.33 The solvent-exclusion surface is

Kv1.2 Voltage-Gated Potassium Channel composed of points at which probe spheres approach closest to atoms of the protein.46 The calculated contact area by solvent-excluded surface is about 3.69 nm2 in the open state conformation and about 3.89 nm2 in the closed state conformation. Thus, the side chains form extensive contact with each other in both the open and closed state conformations. The nearly same contact area in the two conformations indicates that the two parts stick to each other closely in the whole process of gating. From this close contact, the S4-S5 linker couples the movement of the voltage sensor to the channel gating action. The voltage sensors are located at the four corners of the pore domain independently. How do four voltage sensors control the gating cooperatively? The inter-subunit interaction between the sensor part and pore domain may play an important role.5 S1 and S4 of one subunit make contact with S5 of the adjacent subunit. The contact area of the two neighbor subunits calculated by solvent-accessible surface after 35 ns keeps ∼8.10 nm2 in the open state conformation and ∼10.27 nm2 in the closed conformation, as shown in Figure 12. The intersection area between the voltage sensor and pore of the 45 ns frame calculated by solvent-excluded surface is 4.95 nm2 in the open state and 5.87 nm2 in the closed conformation. The difference of the contact area between the two conformations mainly comes from decreased S1-S5 contact area and increased S4-S5 contact area. By solvent-accessible-surface calculation, S1-S5 is decreased by ∼1.41 nm2 and S4-S5 is increased by ∼2.70 nm2 from the open to closed conformation. These measurements agree with our observation of the snapshots of open and closed conformations. From the comparative analysis, we guess, during repolarization, S4 tilts toward the membrane, the intracellular end of S4 goes across the middle between S1 and S5 from the adjacent subunit, and the intracellular ends of S1 and S5 separate further, while the extracellular end of S1 keeps close to the end of S5; thus, the intersection between S1 and S5 is decreased, the intersection between S4 and S5 is increased, and S4 makes some contact with S6. The increased contact between S4 and the pore domain perhaps increases the stability of the protein in the membrane environment. Although the contact area between the S1 and S5 decreased in the closed state, the Thr184 residue on S1 keeps contact with Phe348 on S5 from the neighbor subunit in both open and closed conformations, which is thought to make the voltage sensor paddle transmit force effectively to the S4-S5 linker helix.35 3.4. Lipid Bilayer Environment. The membrane is a disordered environment. Lipid molecules have strong electrostatic interaction with protein, as shown in Figure 2. The pulled down lipid head discussed before in this paper makes the membrane distorted. Lipid molecules fill in the space between the voltage sensor and the pore domain to provide stability to the channel by steric interaction. We calculated the twodimensional number density of the final frame of our simulation at 45 ns, as shown in Figure 13. In the X-Y plane, the lipids are distributed evenly except the lipids that enter the protein zone, and we did not find that lipid molecules in the concave area between the voltage sensors are the most dense, which is mentioned in ref 35. The ordering of the lipid acyl tails is measured by order parameter according to S ) 0.5〈3 cos2 θ 1〉, where θ is the angle of the molecular axis and the bilayer normal and the molecular axis is defined as the vector connecting the two acyl chain carbons Cn-1 and Cn+1, and shown in Figure 14. The order parameter of the open state conformation

J. Phys. Chem. B, Vol. 112, No. 51, 2008 16973 is smaller than the closed state perhaps because the protein’s lipid-facing interface of the open state conformation is more irregular. 4. Conclusion The Kv1.2 voltage-dependent potassium ion channels of both the open and closed conformations embedded in fully hydrated lipid bilayers have been studied by atomistic force field based MD simulations. The potential energy, RSMD, and H-bond interaction between protein and lipid molecules are measured to confirm the simulation time is long enough for sampling. We analyzed the protein stability during simulation, the gating-charge-carrying residue environment, and salt bridge within the voltage sensor. The four subunits lost the symmetry from simulation starting structure, especially the voltage sensors. By comparison of the open and closed state conformations, we calculated the property of S4 segment movement. During depolarization, the mass center of S4 has a vertical distance of 0.672 nm, the tile angle change is ∼25°, and the azimuthal rotation angle change is ∼120°. Thus, the model in our simulation is more near the helix screw model except for the small tilt angle. The S4-S5 linker keeps stable in the whole process of simulation and sticks to the intracellular end of the S6 segment. The contact area between the S4-S5 linker is nearly the same value in the open and closed state conformations. The contact area between neighbor subunits shows that the S1-S5 intersection decreases while S4-S5 intersection increases from the open to closed state. Interaction between protein and lipids takes a long time to reach equilibrium. After equilibrium, lipids distribute evenly around the protein. The smaller order parameter of the lipid tail in the open state conformation perhaps results from the greater irregularity of the protein’s lipid-facing interface of the open state conformation. Until the closed conformation is solved, how the closed state conformation differs from the open state conformation and the details of gating remain unknown. The gating process perhaps includes several intermediate conformations and is not so simple like what we describe by only two end states. We have provided some valuable information on the study of the Kv channel. Certainly more real experimental and computational work will be needed to further the research in this field. Acknowledgment. This work is partially supported by the National Science Foundation of China (Grant No. 20773060) and National Basic Research Program of China (Grant No. 2004CB719901). References and Notes (1) Hille, B. Ion Channels of Excitable Membranes, 3rd ed.; Sinauer: Sunderland, MA, 2001. (2) Bezanilla, F. Physiol. ReV. 2000, 80, 555. (3) Sigworth, F. Q. ReV. Biophys. 1994, 27, 1. (4) Yellen, G. Nature 2002, 419, 35. (5) Tombola, F.; Pathak, M. M.; Isacoff, E. Y. Annu. ReV. Cell DeV. Biol. 2006, 22, 23. (6) Sansom, M. S. P.; Shrivastava, I. H.; Bright, J. N.; Tate, J.; Capener, C. E.; Biggin, P. C. Biochim. Biophys. Acta 2002, 1565, 294. (7) Roux, B. Curr. Opin. Struct. Biol. 2002, 12, 182. (8) Ash, W. L.; Zlomislic, M. R.; Oloo, E. O.; Tieleman, D. P. Biochim. Biophys. Acta 2004, 1666, 158. (9) Gumbart, J.; Wang, Y.; Aksimentiev, A.; Tajkhorshid, E.; Schulten, K. Curr. Opin. Struct. Biol. 2005, 15, 423. (10) Roux, B.; Schulten, K. Structure 2004, 12, 1343. (11) Tieleman, D. P. Computer simulations of transport through membranes: passive diffusion, pores, channels and transporters. Proc. Aust. Physiol. Soc. 2006, 37, 15.

16974 J. Phys. Chem. B, Vol. 112, No. 51, 2008 (12) Deol, S. S.; Domene, C.; Bond, P. J.; Sansom, M. S. P. Biophys. J. 2006, 90, 822. (13) Jiang, Y.; Lee, A.; Chen, J.; Ruta, V.; Cadene, M.; Chait, B. T.; MacKinnon, R. Nature 2003, 423, 33. (14) Long, S. B.; Campbell, E. B.; MacKinnon, R. Science 2005, 309, 897. (15) Treptow, W.; Tarek, M. Biophys. J. 2006, 90, L64. (16) Jogini, V.; Roux, B. Biophys. J. 2007, 93, 3070. (17) Pathak, M. M.; Yarov-Yarovoy, V.; Agarwal, G.; Roux, B.; Barth, P.; Kohout, S.; Tombola, F.; Isacoff, E. Y. Neuron 2007, 56, 124. (18) Kandt, C.; Ash, W. L.; Tieleman, D. P. Methods 2007, 41, 475. (19) Lindahl, E.; Hess, B.; van der Spoel, D. J. Mol. Model. 2001, 7, 306. (20) van der Spoel, D.; Lindahl, E.; Hess, B.; van Buuren, A. R.; Apol, E.; Meulenhoff, P. J.; Tieleman, D. P.; Sijbers, A. L. T. M.; Feenstra, K. A.; van Drunen, R.; Berendsen, H. J. C. Gromacs User Manual, version 3.3; 2005 (www.gromacs.org). (21) Aqvist, J.; Luzhkov, V. Nature 2000, 404, 881. (22) Berger, O.; Edholm, O.; Jahnig, F. Biophys. J. 1997, 72, 2002. (23) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Hermans, J. Interaction models for water in relation to protein hydration. In Intermolecular forces; Pullman, B., Ed.; D. Reidel Publishing Company: Dordrecht, The Netherlands, 1981; p 331. (24) Darden, T.; York, D.; Pedersen, L. J. Chem. Phys. 1993, 98, 10089. (25) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684. (26) Nose´, S. Mol. Phys. 1984, 52, 255. (27) Hoover, W. G. Phys. ReV. A 1985, 31, 1695. (28) Parrinello, M.; Rahman, A. J. Appl. Physiol. 1981, 52, 7182.

Han and Zhang (29) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. J. Comput. Chem. 1997, 18, 1463. (30) Miyamoto, S.; Kollman, P. A. J. Comput. Chem. 1992, 13, 952. (31) Kabsch, W.; Sander, C. Biopolymers 1983, 22, 2577. (32) Ozdirekcan, S.; Etchebest, C.; Killian, J. A.; Fuchs, P. F. J. J. Am. Chem. Soc. 2007, 129, 15174. (33) Sanner, M. F.; Olson, A. J.; Spehner, J. C. Biopolymers 1996, 38, 305. (34) Humphrey, W.; Dalke, A.; Schulten, K. J. Mol. Graphics 1996, 14, 33. (35) Long, S. B.; Tao, X.; Campbell, E. B.; MacKinnon, R. Nature 2007, 450, 376. (36) Seoh, S. A.; Sigg, D.; Papazian, D. M.; Bezanilla, F. Neuron 1996, 16, 1159. (37) Aggarwal, S. K.; MacKinnon, R. Neuron 1996, 16, 1169. (38) Sands, Z. A.; Sansom, M. S. P. Structure 2007, 15, 235. (39) Ahern, C. A.; Horn, R. Trends Neurosci. 2004, 27, 303. (40) Durell, S. R.; Shrivastava, I. H.; Guy, H. R. Biophys. J. 2004, 87, 2116. (41) Lecar, H.; Larsson, H. P.; Grabe, M. Biophys. J. 2003, 85, 2854. (42) Jiang, Y.; Ruta, V.; Chen, J.; Lee, A.; MacKinnon, R. Nature 2003, 423, 42. (43) Ruta, V.; Chen, J.; MacKinnon, R. Cell 2005, 123, 463. (44) Long, S. B.; Campbell, E. B.; MacKinnon, R. Science 2005, 309, 903. (45) Eisenhaber, F.; Lijnzaad, P.; Argos, P.; Sander, C.; Scharf, M. J. Comput. Chem. 1995, 16, 273. (46) Greer, J.; Bush, B. L. Proc. Natl. Acad. Sci. U.S.A. 1978, 75, 303.

JP807905P