Originally published as Biophys J. BioFAST on June 10, 2005.
doi:10.1529/biophysj.105.060368
Biophysical Journal 89:1669-1680 (2005)
© 2005 The Biophysical Society
Theoretical Studies of the M2 Transmembrane Segment of the Glycine Receptor: Models of the Open Pore Structure and Current-Voltage Characteristics
Mary Hongying Cheng *,
Michael Cascio 
and
Rob D. Coalson *
* Department of Chemistry, and Molecular Biophysics Program, University of Pittsburgh, Pittsburgh, Pennsylvania 15260; and Department of Molecular Genetics and Biochemistry, University of Pittsburgh School of Medicine, Pittsburgh, Pennsylvania 15261
Correspondence: Address reprint requests to Michael Cascio, E-mail: cascio{at}pitt.edu.
 |
ABSTRACT
|
|---|
The pentameric glycine receptor (GlyR), a member of the nicotinicoid superfamily of ligand-gated ion channels, is an inhibitory Cl– channel that is gated by glycine. Using recently published NMR data of the second transmembrane segment (M2) of the human 
1 GlyR, structural models of pentameric assemblies embedded in a lipid bilayer were constructed using a combination of experimentally determined constraints coupled with all-atom energy minimization. Based on this structure of the pentameric M2 "pore", Brownian dynamics simulations of ion permeation through this putative conducting open state of the channel were carried out. Simulated I-V curves were in good agreement with published experimental current-voltage curves and the anion/cation permeability ratio, suggesting that our open-state model may be representative of the conducting channel of the full-length receptor. These studies also predicted regions of chloride occupancy and suggested residues critical to anion permeation. Calculations of the conductance of the cation-selective mutant A251E channel are also consistent with experimental data. In addition, both rotation and untilting of the pore helices of our model were found to be broadly consistent with closing of the channel, albeit at distinct regions that may reflect alternate gates of the receptor.
|
INTRODUCTION
|
|---|
The glycine receptor (GlyR) is a member of the nicotinicoid receptor superfamily (also referred to as the Cys-loop receptor superfamily due to an absolutely conserved disulfide loop found in its extracellular domain) whose members are essential mediators in the propagation of electrical signals between cells at the nerve/nerve and nerve/muscle synapse. These receptors are membrane proteins that, in response to neurotransmitter binding, transiently open pores through the lipid bilayer in which they are embedded, allowing the passive movement of small ions down their electrochemical gradient. The GlyR, located predominantly in the spinal cord and lower brain, is an inhibitory neurotransmitter-gated channel which, as its name implies, is typically activated by glycine binding (1
–3). Binding of agonist gates a chloride-conducting channel that results in synaptic inhibition by hyperpolarization of the cell. Other inhibitory channels in this superfamily include the 
-aminobutyric acid type A and C receptors (GABAAR and GABACR, respectively). Excitatory members of the family include the nicotinic acetylcholine receptor (nAChR) and the serotonin type 3 receptor (5-HT3R).
Members of the nicotinicoid receptor superfamily have considerable structural and functional homology. These channels are formed by a pentameric arrangement of subunits, each sharing a common topology having a large extracellular domain (ECD) and four transmembrane domains (M1-M4). For the GlyRs, the site of agonist binding is assumed to be a water-filled cavity at the interface between the ECDs of neighboring subunits. The volume of the cavity expands and contracts during channel opening and closing (4,5). The recent determination of a crystal structure of the pentameric acetylcholine binding protein (AChBP) (6), a soluble homolog of the ECD of the nicotinicoid receptor superfamily, provides details about this domain and a template for modeling ligand binding at the subunit interfaces. Lysine scanning mutagenesis studies conducted on the ligand binding domain of the 
-subunit of nAChR confirmed that the structure of the AChBP is homologous to that of the ECD of nicotinicoid receptors (7). In fact, coupling of the AChBP with the pore domain of the serotonin type-3A receptor produced functional chimeric receptors (after mutagenesis of interacting contact loops for compatibility) (8).
The M2 domain (amino acids 252–270 in 1 GlyR are presumed to span the acyl chains) contributes to the lumen of the central ion-transport pore, which undergoes a conformational change to allow ion flow after the binding of an extracellular ligand (9). The nicotinicoid receptor superfamily generally has three rings of charge associated with their channel pore. The cytoplasmic ring of charge is negative in all members of the superfamily, whereas the charged rings near to the intracellular and extracellular channel mouth are positive for anion selective channels and negative for cation selective channels, respectively. Mutations to the GlyR that neutralize the extracellular ring (R271L and R271Q) produced a reduction in the sensitivities of glycine-activated currents and redistribution of the channel current to lower conductance levels (10,11). The charged ring near the intracellular mouth is critical to ion transport and thought to be the most constricted region in the channel (12–15). For the 1 GlyR, elimination of this charge (R252A) in TM2 gave rise to cells displaying no glycine-gated currents and precluded expression of functional GlyRs (16). A single point mutation in this constricted region of 1 (A251E) was sufficient to convert the channel from being anion selective to cation selective, whereas the deletion of P250 (SDM) from this mutant channel increased the cation current and allowed the permeation of divalent cations (13). Subsequent mutations of R271 (SDM + R271A, or SDM + R271E) changed the conductance and rectification ratio of the mutant cation-selective channel. All four cation-selective GlyR mutants have a minimum pore size comparatively larger than the corresponding anion-selective channels (12). The effect of the channel pore size on the ion permeation and selectivity is assumed to be associated with the ion dehydration process (12).
At present, there does not exist a complete experimental structure of any member of the nicotinicoid family except for the nAChR channel (17,18), for which cryoelectron microscopy structures are available (1OED and 2BG9 in the Protein Data Bank (PDB)). However, given the significant sequence and proposed topological similarities between family members, the structural architecture of any one of these neuroreceptors is believed to be archetypal for the family of ligand-gated channels. Recent NMR data on one TM2 transmembrane segment of the GlyR (19) as well as the cryoelectron microscopy structure of nAChR (17) can be used to construct a reasonable model of the GlyR channel. The goal of this study is to build a model of the pore-lining channel structure of GlyR using molecular dynamics as well as relevant experimental data, and to investigate the permeation mechanism via Brownian dynamics simulation. In addition, our molecular dynamics studies may provide some insight into the location of a potential gate(s) in the channel.
An outline of the article is as follows. In the Methods section we discuss the protocol for determining open as well as closed states of the GlyR pore. Furthermore, we outline a dynamic Monte Carlo method for carrying out Brownian dynamics (BD) simulations of ion permeation in our model system. In the Results section, we present our candidate for the putative open pore state and experimentally relevant current-voltage (I-V) curves calculated via BD simulation for this pore and for the pore of a closely related mutant channel, namely the cation-selective A251E GlyR. In the Discussion section, we address some of the limitations of our model, as well as subtleties in calculating current rectification ratios in the GlyR. We also describe preliminary studies of pore-closing transformations.
|
METHODS
|
|---|
Assembly of the pentameric GlyR M2 bundle
Putative open channel
The NMR-derived M2 transmembrane sequence (19) has 28 amino residues: APARVGLGITTVLTMTTQSSGSRASLPK (1MOT in the PDB bank). One residue at the N-terminus and four residues at the C-terminus were truncated to get a better helix representation. The resultant protein segment has the amino acid sequence: PARVGLGITTVLTMTTQSSGSRA, subsequently to be referred as TM2d. A standard N-terminus and C-terminus were patched onto the appropriate final residues in the TM2d segment. A pentameric GlyR M2 bundle for the putative open state was constructed using VMD (20) as follows:- TM2d was aligned with the channel axis z-direction to make it untilted.
- The untilted TM2d was rotated 15° x N (
) along its helix axis.
- In the open nAChR (21), the channel narrows fairly uniformly from the extracellular to the intracellular with the helix tilting 12°. Thus, the rotated GlyR TM2d sequence obtained in step 2 was tilted 12° to form a constricted section near the intracellular opening (12–15).
- Five copies of the TM2d segment obtained from step 3 were generated with 0°, 72°, 144°, 216°, and 288° rotations about the channel axis. The radial displacement was adjusted to have a minimum pore size
3 Å, which is close to the GlyR wild-type value 2.7 Å (22). Radii of the model channel were measured using the HOLE program (23).
- The open-channel candidates generated in steps 1–4 above were sorted based on the following experimental evidence or experimentally inspired suggestions: a), T258 (T6'), M263 (M11'), S267 (S15'), S270 (S18'), R271 (R19'), and A272 (A20') face the channel lumen (4,5,24); b), T264 (T12') does not face the channel lumen (5); c), R252 (R0') faces the channel lumen, as suggested by Keramidas et al. (12). Viable open-channel candidates were subjected to extensive all-atom energy minimization, as described below.
Putative closed channel
Putative closed forms of the channel were generated in two different ways. In one case an initial model was constructed in which the five TM2d sequences that form the putative open channel were rotated 15° (in a counterclockwise sense as viewed from the extracellular side of the channel). Alternately, the helices were untilted (i.e., aligned in the channel z-direction). In both cases the channel was then subjected to all-atom energy minimization.
All-atom energy minimization
Molecular dynamics (MD) simulations were performed using the NAMD2 program (25) with CHARMM-22 force field for proteins and lipids (26,27). The particle-mesh Ewald method was used to efficiently calculate full electrostatics with a periodic boundary condition of 70 x 70 x 70 Å. Smoothing functions were applied to both the electrostatics and van der Waals forces with a switching distance of 10 Å. The cutoff distance for the nonbonded interactions was 12 Å with the pair list distance extended to 14 Å. All covalent bonds between hydrogen and a parent atom (C, N, O, S) were fixed at the nominal bond length given in the CHARMM-22 parameter file. The time step was chosen to be 2 fs. Short-range nonbonded interactions were calculated every time step and the full electrostatics was evaluated every two time steps. The conjugate gradient and line search algorithm was used for energy minimization.
For computational efficiency, MD simulations were conducted in two steps. In the first step, a solvated lipid bilayer template was generated. A pentameric GlyR bundle was inserted into the center of a cylinder of lipid, formed by three concentric layers of a POPE lipid bilayer (cf. Fig. 1). POPE has been widely utilized as a lipid constituent to model biologically and experimentally relevant bilayers (28,29). The inner and outer radii of the POPE lipid cylinder are 20 and 34 Å, respectively. Fully equilibrated TIP3 waters were added to the system to form a hexagonal boundary condition of 70 x 70 x 70 Å. Any water molecules overlapping the lipid or protein were removed. To make the simulation system electrically neutral, 10 water molecules were selected randomly and replaced by Cl– ions. MD simulation to solvate the lipid bilayer template was performed with the GlyR bundle fixed. First the system was subjected to 2500 steps of energy minimization with a system pressure of 0 bar. Then the energy-minimized system was heated to 310 K for 10 ps with zero pressure applied. This procedure allows the system to relax and results in a change in the lipid packing. The relaxed system was further subjected to the Nose-Hoover constant pressure and temperature procedure (30,31) for 2 ns. A pressure of 1 bar was applied in all directions with constant surface area in the XY plane and a flexible cell in the channel axis z-direction. After 2-ns Nose-Hoover constant pressure and temperature simulation, the lipid bilayer had an area per lipid molecule of 53 Å2 (the estimated protein area is close to 1000 Å2), which falls within the range previously obtained for equilibrated POPE lipid patches (28). Finally, protein was removed and the PDB file for the resultant system, consisting of 10 Cl– ions, 100 lipids, and 4518 waters, was saved as a solvated bilayer template, to which a pentameric GlyR bundle can be docked.
View larger version (91K):
[in this window]
[in a new window]
|
FIGURE 1 The pentameric GlyR bundle inserted into the center of a cylinder of lipid formed by three concentric layers of a POPE lipid bilayer.
|
|
In the second step, the set of pentameric GlyR M2 bundles described above and prepared by VMD was inserted, one by one, into the center of the presolvated POPE bilayer template. In each case, the docked system was first energy minimized for 5000 steps and then subjected to a 200-ps equilibration simulation with protein fixed to improve the packing of the GlyR bundle inside the lipids. After that, the whole system was further energy minimized without any restraint until the relative energy change was <10–5 (32).
A dynamic Monte Carlo algorithm for performing Brownian dynamics simulations of ion permeation
Graf et al. (33,34) have developed a dynamic lattice Monte Carlo (DLMC) algorithm to investigate ion permeation through protein channels. DLMC is performed on a lattice using the rules of dynamic Monte Carlo, which generate kinetics equivalent to the high-friction limit of Brownian motion theory, i.e., Brownian dynamics (35). In the DLMC approach, configurations are generated by random changes of the ion positions. The total number of ions is characterized by N = NL + NR + NI + Nv. Here NL and NR are the fixed numbers of ions on the left and right boundaries, which are obtained by integrating the given boundary concentrations CL and CR over the volumes of the boundary layers. NI is the number of ions inside the system and Nv is the number of virtual ions. The total number of ions N is fixed and NI fluctuates. Nv is introduced only for counting purposes and is included to account for dynamic fluctuation of the number of ions in the interior of the system, and ensure the proportionality of Monte Carlo cycles to real time (characterizing time evolution under the corresponding multi-ion Smoluchowski equation). One Monte Carlo (MC) cycle consists of N steps. At each step, one ion k is randomly chosen to move ±1 lattice point in one direction (x, y, or z) if the chosen ion is not a virtual one. In this implementation, this new configuration is accepted if 
, where Rand is a uniform random number that lies in the interval [0, 1] and 
; 
W is the energy change between the two configurations based on the chosen particle k with charge qk (36):
 | (1) |
 | (2) |
The first term on the right-hand side of Eq. 2 is the electrostatic potential energy of ion k due to fixed charges on the protein/membrane and any externally applied electric potential: 
stat can be obtained by solving Poisson's equation:
 | (3) |
based on the complete (spatially dependent) dielectric constant profile and the specified boundary potential, with 
being the charge density due to the permanent charges and dipoles in the protein and membrane. The second term is the dielectric self-energy that accounts for the image potential due to the surface charge induced on dielectric boundaries by ion k. The third term represents Coulombic interactions between pairs of ions in an environment characterized by the (uniform) water dielectric constant 
w, and the fourth term corresponds to the image potential experienced by ion k due to the surface charge induced by ion j. The last two terms in Eq. 2 together represent ion-ion interactions in a dielectrically inhomogeneous environment. (Full computational details may be found (33,34,36)). For a simulation system with a uniform diffusion constant D for the ions, the total simulation time Ts is related to the number of the Monte Carlo cycles Nc as:
 | (4) |
where h is the grid size.
This dynamic Monte Carlo algorithm was adapted to calculate ion permeation in the GlyR channel. One feature of this GlyR model requires modification of the algorithm used in previous channel studies carried out via DLMC (33,34,36). In a fully realistic model of permeation through a narrow ion channel, allowance should be made for the fact that the diffusion constant characterizing the motion of an ion inside the channel is different from its bulk diffusion constant. When the diffusivity is nonuniform, its gradient must be incorporated into the drift term of the Brownian dynamic algorithm (37,38). To include this effect in the dynamic Monte Carlo algorithm, we add an additional single-particle effective potential to the overall configuration energy function (Eq. 2). That is, 
with
 | (5) |
where D0 and D(z) are the diffusivities in the bulk and at position z, respectively. Although this correction term can be applied to treat any given three-dimensional inhomogeneous diffusion constant profile, in the model considered in this article the diffusivity varies only along the channel direction. The diffusion constants and the associated ion displacement steps taken in one Monte Carlo step obey the relation:
 | (6) |
where h0 and hz are the associated ion displacement in the bulk and position z, respectively. Note that the procedure just outlined accounts exactly for the effect of inhomogeneous diffusivity on the BD algorithm (37,38).
Electrostatic calculations (solutions of the Poisson equation) were carried out on a 105 x 105 x 115 lattice with a lattice spacing of 1 Å. The actual dynamic Monte Carlo calculations were carried out via single ion moves of ±hz in the x-, y-, or z-direction. Due to movement of ions between different diffusion constant regions at different times, there is no "universal" lattice that characterizes the overall dynamics. This does not cause any inherent difficulties. In particular, the electrostatic energy change needed to decide whether to accept a trial move (cf. Eqs. 1 and 2) can be calculated by interpolation based on a multilinear interpolation algorithm (39). Furthermore, the condition of zero overlap (excluded volume) between an ion and the protein/membrane and between any two ions is easily enforced in continuous space, based on an assumed radius for each ion: the radii of Cl– and Na+ are taken to be 1.8 and 1.0 Å, respectively.
|
RESULTS
|
|---|
Putative open-channel structure and pore architecture
The structure of a multimeric assembly of the TM2 domain of human 1 GlyR in detergent micelles was determined by NMR (19). There are 20 different NMR-derived wild-type (WT)-GlyR TM2 structures deposited in the PDB (ascension No. 1MOT). All structures show the presence of a helical motif with restricted internal motions close to the N-terminus, and a random coil fragment at the C-terminus, which may be located outside the micelle during sample preparation (19). A well-defined -helix runs from V253 to T265. Of the 15 lowest energy NMR-derived TM2 structures considered (see Fig. 2 a for backbone traces), 75% of them show a coincident position of R252, which points in the same direction as M263 (Fig. 2 b). Our putative open channel was constructed using a truncated portion of the NMR structure (TM2d, residues 250–272) as a template for the transmembrane pore helix. It has been demonstrated that the charged ring at the intracellular mouth of the nAChR is exposed to permeant ions (14). Therefore, our putative open-channel structure was constructed based on the assumption that the similarly located R252 residues of the 1 subunits of the GlyR face the pore lumen and form the constricted region when the channel is open (as suggested by Keramidas et al. (12)). It should be noted that our pentameric GlyR bundle was constructed in a manner somewhat different from a previous model (32), in which the NMR-derived M2 sequences were docked into homopentameric channels with the T258 side-chain vector pointing toward the pore axis in both the open and closed structures.
View larger version (46K):
[in this window]
[in a new window]
|
FIGURE 2 (a) Backbone traces of the 15 lowest energy TM2 GlyR structures listed in the PDB bank (1MOT). These are characterized by a helical structure in the region close to the N-terminus (bottom) and a relatively large unstructured fragment at the C-terminal side (top). A solid -helix runs from V253 to T264. (b) Residue traces of the 15 lowest energy TM2 GlyR structures: 75% of these show a coincident position of R252 (shown inside the dotted circle).
|
|
As discussed in the Methods section above, to save computation time, a presolvated bilayer system was used as a template and a pentameric GlyR M2 bundle was inserted into its center. Waters were found to penetrate the lipid bilayer significantly, indicating effective solvation of the lipid bilayer. A typical MD snapshot of the docked simulation system is shown in Fig. 3. Viable open-channel candidates were subjected to intensive all-atom energy minimization. After energy minimization, the system with the smallest energy and the minimum pore radius not <2.5 Å was selected as our putative open channel. Fig. 4 a shows the putative open-channel architecture of the TM2d pentameric bundle as viewed from the intracellular surface, with lipid and water removed for display purposes. The TM2d structure in our energy-minimized pentameric bundle remained very similar to the NMR-derived structure. In the final energy-minimized structure of our putative open state, a Cl– ion was found to reside near the pore section formed by R252 (Fig. 4 b). This ion may stabilize the pore channel and affect the orientation of R252.
View larger version (136K):
[in this window]
[in a new window]
|
FIGURE 3 MD snapshot of the simulation system after the solvation stage. The pentameric GlyR bundle (in ribbon format) was embedded in a POPE lipid bilayer solvated with TIP3 water (red dots) on either side. The phosphorous atoms of the lipid headgroups are shown as gray spheres. N-terminus of M2 is on the intracellular surface, located at the bottom in this view.
|
|
View larger version (36K):
[in this window]
[in a new window]
|
FIGURE 4 Putative open structure after energy minimization. (a) N-terminus view. Lipids, ions, and water have been removed. Blue, cyan, brown, and yellow represent R252, A251, P250, and V253, respectively. (b) Possible binding site for Cl– near R252. Lipids and water have been removed. Cyan spheres represent Cl– ions. Arginines at positions 252 and 271 are indicated as blue stick figures.
|
|
In substituted cysteine accessibility method (SCAM) studies conducted on GlyR, Lobo et al. (5) observed that MTS reagents added extracellularly do not measurably react with residues below M263 in both the open and closed state (with the caveat that modified Cys residues alter channel function). M263C reacted with propylMTS, but not with decylMTS, suggesting that decylMTS reagents are too large to access M263C (5). Thus, it is inferred that M263C forms a relatively constricted region in both the open and closed states of the receptor. In our model channel, M263 was found to extrude into the pore lumen in the putative open state (cf. Fig. 5), with an average pore radius 4.0 Å, which would be comparable to the size of propylMTS (propylMTS has a volume of 102 Å3 (5)). Consistent with its presumed inaccessibility in SCAM studies, T264 in the model channel was found to be occluded behind M263 and did not face the pore lumen. Also, as observed in experimental studies (5,24), residues T258, S267, and S270 in the model channel were found to be accessible.
View larger version (27K):
[in this window]
[in a new window]
|
FIGURE 5 (a) Side view of M2 pore along channel axis. Critical residues in the pore-lining open state of the GlyR model channel are shown via space-filling models. (b) The model channel of the GlyR. The intracellular and extracellular protein are represented by funnel-like pores attached on either side of the transmembrane (membrane/protein) segment. The protein channel extends from z = –32 (Å) to z = 32 (Å). Two dipole rings ( ) account for the electrostatic potential generated by the hydroxyl groups T259 and T258 at z = –8 (Å), and G269 and S267 at z = 4 (Å), with the positive end 1 Å inside the channel and pointing toward the channel pore. The intermediate (INT.) and extracellular (EXT.) rings of charge R252 and R271 (+) are located 1 Å inside the channel at z = –18 (Å) and z = 7 (Å), respectively. The cytoplasmic (CYT.) ring of charge D247 (circles) is located 3 Å inside the channel at z = –24 Å.
|
|
Ion permeation
To test our model of the open channel, we compare our BD simulation results (implemented throughout using the dynamic Monte Carlo algorithm described in the Methods section) for current-voltage relationships with experimental studies conducted previously by one of the authors (40). In these experimental studies, functional recombinant homomeric human 1GlyR motifs were overexpressed in insect cell culture, solubilized, purified, and reconstituted into lipid vesicles via gel filtration. In subsequent black lipid membrane studies, the voltage dependence of single-channel current was measured.
For simplicity, we assume the homomeric 1 M2 GlyR channel to be cylindrically symmetric, thus significantly reducing the number of points that must be stored in the electrostatics look-up table (39), as well as the number of times that Poisson's equation has to be solved. Essentially the number of Poisson's equation solutions is reduced by "one dimension" (the azimuthal angle in cylindrical coordinates), which saves a factor of 10–100 depending on how fine a grid is required. Fig. 5 a shows the critical side chains of transmembrane pore lining amino acids, and Fig. 5 b illustrates the channel geometry and the position of residue charges and dipoles in our model channel. In Fig. 5 b, funnel-like protein density was attached on either side of the transmembrane pore to represent intracellular and extracellular protein. The membrane/protein/water system under study is dielectrically inhomogeneous. Whereas the bulk water is well characterized by a static dielectric constant of 80, the appropriate dielectric constant in the pore interior is not known experimentally and presumably it is lower than the bulk value. The value of the effective dielectric constant is estimated to range between 2 and 5 in nonpolar regions of protein interior (41) and various flexible polar groups in the protein suggest the dielectric constant of the protein is typically in the range of 2–10 (42,43). Despite all these complications, for simplicity, we assume here that the dielectric constant in the water bath w is the same as that inside the channel, namely w = 80 in all aqueous regions. Because the outer radius of the extracellular protein mass in nAChR is 35 Å (inner radius is 10 Å), the electrostatic effects of any water surrounding this extracellular protein are expected to be negligible (little effect on the self-energy of a single ion in the channel or the electrostatic interaction energy of two ions in our model channel was observed in simulations conducted under these conditions). Thus, we assume a uniform dielectric constant of m = 5 for the entire protein-membrane system, extending this region laterally outward from the protein pore (including in the extracellular and cytoplasmic domains) as indicated in Fig. 5 b.
The constructed model channel extends from z = –32 Å to z = 32 Å. Because the 1 GlyR-WT is pentameric, dipoles and residue charges are assigned five times, 72° apart. At present it is not possible to accurately compute pKa values for charged residues in GlyR due to structural and dielectric uncertainties. Thus, ring charges and dipole moments were empirically assigned based on our model-channel geometry, dielectric constant profile, and experimental IV curves. Two dipole rings were inserted into the channel to create an electrostatic potential, which mimics that generated by the hydroxyl groups T259 and T258, and S267 and G269 (see Fig. 5 caption for further details). Arginine rings at R252 and R271 were inserted at z = –18 Å and z = 7 Å, respectively. An aspartate ring (D247) was inserted at z = –24 Å. The effective dipole moment was assumed to be 0.3 D. The effective residue charges of D247, R252, and R271 were chosen to be –0.2, 0.3, and 0.2 e, respectively (e is the proton charge). The effective charge or dipole moment used in a continuum model is an approximate distillation of complicated microscopic phenomena. It is affected, e.g., by screening effects due to the neighboring amino acids or water. For this reason, the effective charges employed in this (or any other) study cannot be uniquely determined. The charge set specified above was motivated by the following considerations. Because R252 resides in a very constricted section of the pore, the occupancy of water molecules near this section is limited and thus the effective charge might be higher than that of R271, which resides in a water cavity (4,5). While keeping all the other charges and dipole moments the same, assigning the charge of R252 to be 0.2 vs. 0.3 e leads to variations in conductance of <20%. However, choosing the charge of R252 to be 0.3 e produced results consistent with MD investigation of the possible binding site for Cl– described in the previous subsection. Fig. 6 shows the electrostatic energy profiles associated with a single Na+ and Cl– ion as either of these passes through an otherwise unoccupied channel. Based on the choices of dielectric constant profile, charges and dipole moments specified above, the calculated Na+ energy barrier is found to be 15.6 kBT, somewhat smaller than the value of 27 kBT obtained by O'Mara et al. (44). We find an energy well for Cl– about –8.7 kBT near the intermediate ring, in contrast to the –16 kBT energy well obtained by O'Mara et al. (44).
View larger version (13K):
[in this window]
[in a new window]
|
FIGURE 6 Electrostatic energy profiles along the channel axis. The solid line and dashed line represent the electrostatic potential energy experienced by Na+ and Cl–, respectively, when they enter an empty channel. This energy is given by the first two terms of Eq. 2 at an applied voltage of 0 mV.
|
|
To our knowledge, the diffusion constant of an ion inside a protein channel has yet to be experimentally measured for any channel. Based on MD investigations of diffusion constants for K+ and Cl– inside the OmpF porin channel (38) and K+ in a Gramicidin A (45), the diffusion coefficient for both Na+ and Cl– was assumed to be D = 2.0 x 10–5 cm2/s in the bulk, linearly reduced from this bulk value at the channel entrance (zL = –32 Å and zR = 32 Å) to half the bulk value at the transmembrane channel mouth (zL = –27 Å and zR = 10 Å), and maintained at half the bulk value throughout the transmembrane domain –27 Å < z < 10 Å. Unless otherwise explicitly noted, for one MC step the displacement in the bulk is chosen to be 1 Å. Other sets of displacements, e.g., 1.4, 0.7, and 0.5 Å, were also tested; the difference in the calculated currents and concentrations were found to be within 15%. The total simulation time is 4 x 106 MC steps, which corresponds to 3.3 µs. Periodic boundary conditions were applied to the x- and y-directions, whereas for the z-direction fixed concentrations were imposed using a 1 Å layer thickness. To compare to the experimental data in Cascio et al. (40) the values CL = CR = 0.15 M were chosen. For each ion move that resulted in a change of the boundary layer concentration, one ion was randomly added to or removed from the boundary layer to keep the boundary concentration constant.
Brownian dynamics simulations indicated that chloride resides at discrete regions within the pore (Fig. 7): a peak in the Cl– concentration was found near R252, consistent with a possible "binding site" for chloride ions observed in our MD investigations (cf. Fig. 4 b), as well as in the BD simulations of O'Mara et al. (44). The computed IV curve obtained in the simulations with symmetrical 0.15 M NaCl solutions (the bath concentration in our simulation fluctuates as CL = 0.148 ± 0.007 M and CR = 0.145 ± 0.001 M when averaged over the range of applied voltages considered in our calculations) indicates a linear relationship and no clear current rectification (Fig. 8), as experimentally observed (40). The conductance obtained from simulation is 76 pS, which agrees well with the 79 pS obtained in experiments (40), thus supporting the choice of charges and dipole moments.
View larger version (17K):
[in this window]
[in a new window]
|
FIGURE 7 BD simulation result for the average number of Cl– ions inside the channel with an applied membrane potential of 60 mV. The two peaks correspond to the channel sections near the charged rings R252 and R271, as indicated.
|
|
View larger version (11K):
[in this window]
[in a new window]
|
FIGURE 8 Comparison of experimental data with numerical BD simulations. Open diamonds represent the experimental data obtained by Cascio et al. (40). Dots with error bars represent simulation data (based on eight runs). Line shows the least-square fit to the simulation data.
|
|
BD simulations were also performed with the boundary concentrations CL = 0.15 M (intracellular) and CR = 0.075 M (extracellular) to obtain the ion permeability ratio. The reversal potential and the permeability ratio are related through the Goldman-Hodgkins-Katz equation (12):
 | (7) |
where R is the ideal gas constant, F is Faraday's constant, T is the absolute temperature in Kelvin, and Vrev is the reversal potential (at which the total electric current is zero); 
and 
are the activities (12,46) of Na+, Cl– ions in the intracellular and extracellular solution, respectively. Two layers of boundary were used to accommodate the concentration setting and a displacement of 1.4 Å was used in the bath. The average concentrations in the intracellular and extracellular solutions were found to be 0.148 ± 0.006 M and 0.071 ± 0.001 M, respectively. Fig. 9 shows the calculated dilution I-V curve for our model channel. The reversal potential obtained from our BD simulations is 15.5 ± 0.8 mV. Using bath activity coefficients tablulated in standard tables (47), we calculated a permeability ratio of 
which falls within the range of values extracted directly from experiments, namely 24–28 (12,46). We note that other authors have defined the permeability ratio from Eq. 7 with activities replaced by the corresponding bulk ion concentrations (38,48,49). Using these bath concentrations in conjunction with our BD simulation-based reversal potential, gives 
. The factor of two variation for the permeability ratio (as compared to the value we calculated using tabulated activity coefficients) is not unreasonable given the various approximations and simplifications that underlie the BD model.
View larger version (11K):
[in this window]
[in a new window]
|
FIGURE 9 I-V curves obtained via BD simulation for asymmetric bathing solutions. Circles with error bars represent simulation data (based on 12 runs). Line shows the least-square fit of the simulation data. The average concentrations of the intracellular and extracellular solutions are 0.148 ± 0.006 (M) and 0.071 ± 0.001(M), respectively. The reversal potential corresponds to Vrev = 15.5 ± 0.8 mV.
|
|
Clearly our open-channel model of WT-GlyR is highly anion-selective, although rare events of spontaneous transport of Na+ were in fact observed. To further test this model, we carried out simulations of an A251E mutant GlyR. This single substitution introduces a negatively charged ring E251 adjacent to a ring of positive charge R252 and produces a preferentially cation-selective channel with an experimentally determined 
(12,13). Based on the suggestion by Keramidas et al. (12) and O'Mara et al. (44), the charge ring R252 was moved 4 Å deeper inside the channel and another ring of charge E251 was added 1 Å inside the channel. While keeping the pore radius and other charges and dipole moments the same, the charge of R252 and E251 were assigned to be 0.2 e and –0.2 e, respectively. In a symmetrical bath of 0.145 M NaCl, the experimentally determined single-channel currents at 60 and –60 mV are 0.66 ± 0.04 and 0.41 ± 0.01 pA (13) and the corresponding simulated values are 0.72 ± 0.15 and 0.29 ± 0.08 pA. In addition to the permeation of Na+, we observed counterion fluxes of Cl– (in the opposite direction to Na+, roughly in the proportion 
at a membrane potential of 60 mV). Numerical dilution experiments were performed with the boundary concentrations CL = 0.30 M (bath concentration fluctuates as 0.306 ± 0.006 M) and CR = 0.15 M (bath concentration fluctuates as 0.148 ± 0.001 M), and we obtained a reversal potential value of –14.0 ± 1.5 mV, from which we calculated a permeability ratio 
(based on the average bath activities), which is reasonably close to the experimental value and indeed displays preferential cation selectivity.
|
DISCUSSION
|
|---|
Modeling the pore from limited structural data
As described above, our model of the open state of the 1 GlyR pore region is consistent with a variety of experimental measurements on GlyR (both qualitatively, and, in direct comparison with I-V curves of the homomeric 1 GlyR, quantitatively). However, in the construction of our open-channel structure based on a pentameric arrangement of M2 transmembrane segments (cf. Fig. 10 a) there are several possible limitations that should be noted. First, the GlyR bundle considered in this work is only a simplified model of the channel. Even though it was assembled based on high-resolution NMR studies and the resultant channel structure agrees with many channel features deduced from scanning cysteine accessibility experiments, it is derived from only a small portion of the 421-amino-acid full-length subunit, and does not include any contribution from surrounding transmembrane segments, nor protein-protein contacts with the extracellular and intracellular loops and domains. Furthermore, it has been suggested that the N-terminal ends of M1 residues may also contribute to the channel lumen and gating in some states (50). Second, NMR-derived GlyR sequences show a freely random coil fragment at the C-terminus. Building a channel structure based on a freely random coil may impair the accuracy of the constructed model channel. However, the important amino acid sequences (e.g., R252-T264) that most strongly affect the permeation and selectivity of the channel in our simulations comprise a comparatively solid helix and stay relatively stable. Third, the static model channel employed to investigate the permeation neglects protein dynamics. The protein and membrane are in reality not static objects, but rather fluctuate on picoseconds and longer timescales, whereas an ion permeation event typically takes nanoseconds. Even if no major structural reorganization occurs during the course of ion permeation, these membrane/protein fluctuations influence the energetics of ion permeation in several ways. For example, they contribute to the polarizability of the membrane/protein medium, and pore-lining groups (e.g., carbonyl groups in gramicidin A) form directed electrostatic "bonds" with ions in the channel. Future studies will incorporate the dynamics of the protein into the model, i.e., feed up information extracted from microscopic or semimicroscopic models to refine the parameters that enter into coarse-grained continuum models such as BD (51,52).
View larger version (33K):
[in this window]
[in a new window]
|
FIGURE 10 (a) Putative open structure constructed based on the NMR-derived TM2d sequence; purple bonds represent M263, and pink bonds represent L255. (b) Putative closed structure constructed by rotating the five GlyR TM2d sequences 15° in the counterclockwise sense relative to the putative open structure. (c) Untilted structure. All figures are C-terminus view before energy minimization.
|
|
Current rectification
As shown in Fig. 8, our simulations provided no evidence of rectification, and were consistent with experimental determinations of reconstituted homomeric receptors in planar bilayers (40). A rather wide range of current rectification patterns has been reported in the experimental GlyR literature. Measurements in the native GlyR (spinal cord and cultured postnatal hippocampal neurons) exhibited whole-cell outward current rectification (53,54). In contrast, recombinant 1 GlyR in HEK293 cells elicited a single-channel inward rectification (3). Additionally, whole-cell currents obtained from the human GlyR 1 cDNA-transfected HEK-293 cell displayed no clear rectification and a strong inward rectification was obtained upon extracellular application of 1 µM cyanotriphenylborate (22). This inward rectification was attributed to blocking of the open channel by the negatively charged cyanotriphenylborate from the outside of the cell (22). In general, current rectification in ion channels is a complicated phenomenon and can depend on many parameters, such as receptor subunit composition, posttranslational modifications, and intrinsic channel structure as well as the surrounding solutions. For example, the outward rectification of NMDA receptors was attributed to channel block by external Mg2+ (55), whereas the inward rectification of the native nicotinic acetylcholine was attributed to the intermediate ring of residues at the channel inner mouth (56), and/or blocking the open channel by internal Mg2+ (57). Through extensive BD simulations performed on our model GlyR channel, we have found that a complicated interplay between the dielectric self-energy, electrostatic energy due to the protein structure, and the ion-ion interactions determines the I-V curve. Alterations in the channel structure, dipole moments, and the residue charges give rise to different degrees of current rectification. It is clear that current rectification requires an asymmetric Cl– ion free energy profile (53).
Closed states of the channel
In the nicotinicoid receptor superfamily, the physical location of the gate is fairly controversial. Cryomicroscopic studies of nAChR suggest that the narrowest constriction in the closed state of the receptor is a hydrophobic girdle in the vicinity of L251 (9') and V255 (13') (17,18). However, studies of the accessibility of cysteine mutants introduced into the TM2 domain of nAChR suggest that in the closed state the gate is close to the cytoplasmic face of the channel (58), and is approximately located near positions from –2' to 2' in the nAChR (15).
In addition, two different simple gating mechanisms have been proposed. Based on comparison of the extracellular domain of the nAChR with the AChBP, Unwin et al. (59) postulated that ligand binding induces a conformational change, resulting in an 15° rotational movement of the inner sheets in the ECD of the -subunits. The rotational motion in the inner sheet of -subunits is proposed to be propagated to yield a similar rotation in the TM2 helix, which forms the permeation pathway across the acyl chains of the bilayer. In an alternate model, Tang and co-workers (32) propose a gating mechanism wherein untilting the helix closes the channel. Based on energy minimization of differential packing of their NMR-derived structure of the pore-forming TM2 peptide, they postulate the energy barrier to Cl– flow is a constriction of the pore radius closer to Q266 (14') with channel gating effected by a tilting movement of the TM2 helices that increases the pore radius and forms a continuous hydrophilic surface.
We examined whether either simple helix rotation or helix untilting could limit anion permeation through our putative open model of the channel. In the first case, we rotated each transmembrane M2 monomer 15° counterclockwise (Fig. 10 b as viewed from the extracellular side). We observed a constricted region around L255, with the hydrophobic side chain extruding into the channel lumen, which effectively decreases the minimum pore radius close to 2.0 Å (Fig. 11). This location is broadly consistent with the SCAM studies in nAChR (18), which suggest a resting gate that is very close to the cytoplasmic surface.
View larger version (17K):
[in this window]
[in a new window]
|
FIGURE 11 The radius profiles of homopentameric TM2 channels after extensive energy minimization. Solid line represents the putative open channel. The line with solid diamonds represents the putative closed pore channel resulting from M2 rotation. The line with open diamonds represents pore channel upon untilting of the TM2 helices.
|
|
In our open-channel model, the helix tilt angle is 12°. We next tested whether aligning the helices with the bilayer normal (untilting) might also restrict ion flow. Because the volume of the water cavity region (e.g., M263–V277) has been observed to contract in the closed state (4,5), untilting was carried out with the constricted section in the open channel as the hinge to decrease the volume of the water cavity. In the untilted arrangement of the -helices (Fig. 10 c), M263 (11') forms the most constricted section, whose pore radius changes from 4.0 Å in the tilted configuration to 2.5 Å in the untilted one. The role of M263 as a possible gate may be supported by recent SCAM studies of the GlyR, which indicated that this residue borders a water cavity (5) and may form a constricted region in the channel and hence serve as a hydrophobic gate, where the permeating ion will be fully or partially deprived of its hydration shell (12). Furthermore, this region is consistent with the suggested leucine-valine gate observed near the midpoint of TM2 in cryomicroscopic studies of nAChR (17) and the hydrophobic gate at positions 5', 9', 12', and 16' observed by Kim and co-workers in their model of the closed form of the M2 pore of nAChR (60). These results are also broadly consistent with the location of a constricted region observed in SCAM studies of nAChR in its desensitized state (58).
Thus, both simple rotation and untilting act to constrict the channel in our simulations, albeit at different sites. It should be noted that our way to rotate or untilt helix differs slightly from Unwin et al. (59) or Tang et al. (32). Undoubtedly, these are rather simplistic models as gating is a complex phenomenon, with interconversions between multiple allosteric states as a function of the liganded state of the receptor (dependent on other factors such as lipid composition as well) that may not be so simply reduced to a single general description, i.e., rotation versus untilting. This is not a critique of either gating model, as the general models do distill phenomena observed at low resolution. The important consideration is whether our model is consistent with observed experimental data. Both gating models are consistent with select experimental data, so how can we accommodate both gates? One possibility is that these two broad mechanisms reflect alternate gates of the channel (61). The existence of nonidentical gates in the closed and desensitized states was originally proposed by Auerbach and Akk (62). As described above, initial SCAM studies conducted in the resting, closed state of this receptor showed that externally added reagents could penetrate and quickly access residues near the cytoplasmic surface at positions from –2' to 2' (15), similar to the constricted region observed in our rotation of the M2 helices. In subsequent studies conducted in the presence and absence of ligands, it was proposed that there may be a second, desensitization gate that forms a barrier to permeation around the highly conserved position 9' in the desensitized state (58), close to the constricted region observed in our (un)tilting studies of our open channel. This latter location is also consistent with the hydrophobic girdle gate proposed on the basis of cryomicroscopic studies (17). Given the absolute requirement of the nAChR channel for cholesterol and negatively charged phospholipids for activity (63) and the existence of multiple distinct allosteric states of the receptor as a function of lipid composition (64–66), an intriguing possibility is that the closed state resolved in cryomicroscopic studies of tubular ordered arrays of the nAChR in the absence of ligand may describe a desensitized-like closed state. Regardless, although accurate numerical simulation of the entire gating process of nicotinicoid receptors is not feasible at present due to the lack of structural resolution, the complexity of the multiple allosteric states of the receptor, the long (millisecond) timescales involved, etc., simple rotational and untilting actions are sufficient to constrict our model of the open GlyR channel and do suggest potential roles of these motions in closing the receptor.
|
CONCLUSION
|
|---|
We have constructed a pentameric bundle of the glycine receptor channel via molecular dynamics and related molecular modeling methodologies. The GlyR channel structure was based on the NMR-derived structure of the TM2 pore-lining segment in detergent micelles and was subjected to intensive energy minimization in an appropriate lipid bilayer membrane. The model channel is consistent with experimentally determined characteristics of the GlyR. Significantly, Brownian dynamics permeation studies yielded several current-voltage relationships in agreement with experimental determinations (12,13,40,46), providing some confidence that the constructed open-state structure is an appropriate model for functional simulations. In addition, the concentration profile obtained via BD simulation suggests that there may be a Cl– "binding site" near R252, which has implications not only for permeation mechanisms, but also for structural stability of the channel. Finally, we found that a 15° rotation of each M2 helix (counterclockwise as viewed from the extracellular side) generates a constricted region around L255 with the hydrophobic side chain extruding into the channel lumen, while untilting the open form of the helix causes constriction of the pore near the M263 (11') section. Thus, both types of conformational change generate final channel states that are essentially closed to ion flow. However, due to the extreme complexity of gating processes, the results presented herein constitute at best a hint of the protein transformations involved in gating of the GlyR.
|
ACKNOWLEDGEMENTS
|
|---|
The work of M. Cheng and R. D. Coalson was partially supported by National Institutes of Health grant P1938613 and U.S. Department of Defense Multidisciplinary University Research Initiative (ARO-MURI) grant DADD19-02-1-0227 and National Science Foundation grant CHR0092285. The work of M. Cascio was partially supported by National Institute on Drug Abuse grant DA016381.
Submitted on January 29, 2005;
accepted for publication May 25, 2005.
|
REFERENCES
|
|---|
1. Lynch, J. W. 2004. Molecular structure and function of the glycine receptor chloride channel. Physiol. Rev. 84:1051–1095.[Abstract/Free Full Text]
2. Cascio, M. 2004. Structure and function of the glycine receptor and related nicotinicoid receptors. J. Biol. Chem. 279:19383–19386.[Free Full Text]
3. Keramidas, A., A. J. Moorhouse, P. R. Peter, and P. H. Barry. 2004. Ligand-gated ion channels: mechanisms underlying ion selectivity. Prog. Biophys. Mol. Biol. 86:161–204.[CrossRef][Medline]
4. Lynch, J. W., N. L. Han, J. Haddrill, K. D. Pierce, and P. R. Schofield. 2001. The surface accessibility of the glycine receptor M2–M3 loop is increased in the channel open state. J. Neurosci. 21:2589–2599.[Abstract/Free Full Text]
5. Lobo, I. A., M. P. Mascia, J. R. Trudell, and R. A. Harris. 2004. Channel gating of the glycine receptor changes accessibility to residues implicated in receptor potentiation by alcohols and anesthetics. J. Biol. Chem. 279:33919–33927.[Abstract/Free Full Text]
6. Brejc, K., W. J. van Dijk, R. V. Klaasen, M. Schuurmans, J. van der Oost, A. B. Smit, and T. K. Sixma. 2001. Crystal structure of an ACh-binding protein reveals the ligand-binding domain of nicotinic receptors. Nature. 411:269–276.[CrossRef][Medline]
7. Sine, S. M., H. L. Wang, and N. Bren. 2002. Lysine scanning mutagenesis delineates structural model of the nicotinic receptor ligand binding domain. J. Biol. Chem. 277:29210–29223.[Abstract/Free Full Text]
8. Bouzat, C., F. Gumilar, G. Spitzmaul, H. L. Wang, D. Rayes, S. B. Hansen, P. Taylor, and S. M. Sine. 2004. Coupling of agonist binding to channel gating in an ACh-binding protein linked to an ion channel. Nature. 430:896–900.[CrossRef][Medline]
9. Barnard, E. 1992. Receptor classes and the transmitter-gated ion channels. Trends Biochem. Sci. 17:368–374.[CrossRef][Medline]
10. Rajendra, S., J. W. Lynch, K. D. Pierce, C. R. French, P. H. Barry, and P. R. Schofield. 1995. Mutation of an arginine residue in the human glycine receptor transforms beta-alanine and taurine from agonists into competitive antagonists. Neuron. 14:169–175.[CrossRef][Medline]
11. Langosch, D., B. Laube, N. Rundström, V. Schmieden, J. Bormann, and H. Betz. 1994. Decreased agonist affinity and chloride conductance of mutant glycine receptors associated with human hereditary hyperekplexia. EMBO J. 13:4223–4228.[Medline]
12. Keramidas, A., A. J. Moorhouse, K. D. Pierce, P. R. Schofield, and P. H. Barry. 2002. Cation-selective mutations in the M2 domain of the inhibitory glycine receptor channel reveal determinants of ion-charge selectivity. J. Gen. Physiol. 119:393–410.[Abstract/Free Full Text]
13. Moorhouse, A. J., A. Keramidas, A. Zaykin, P. R. Schofield, and P. H. Barry. 2002. Single channel analysis of conductance and rectification in cation-selective, mutant glycine receptor channels. J. Gen. Physiol. 119:411–425.[Abstract/Free Full Text]
14. Imoto, K., C. Busch, B. Sakmann, M. Mishina, T. Konno, J. Nakai, H. Bujo, Y. Mori, K. Fukuda, and S. Numa. 1988. Rings of negatively charged amino acids determine the acetylcholine receptor channel conductance. Nature. 335:645–648.[CrossRef][Medline]
15. Wilson, G. G., and A. Karlin. 1998. The location of the gate in the acetylcholine receptor channel. Neuron. 20:1269–1281.[CrossRef][Medline]
16. Lynch, J. W., S. Rajendra, K. D. Pierce, C. A. Handford, P. H. Barry, and P. R. Schofield. 1997. Identification of intracellular and extracellular domains mediating signal transduction in the inhibitory glycine receptor chloride channel. EMBO J. 16:110–120.[CrossRef][Medline]
17. Miyazawa, A., Y. Fujiyoshi, and N. Unwin. 2003. Structure and gating mechanism of the acetylcholine receptor pore. Nature. 423:949–955.[CrossRef][Medline]
18. Unwin, N. 2005. Refined structure of the nicotinic acetylcholine receptor at 4 Å resolution. J. Mol. Biol. 346:967–989.[CrossRef][Medline]
19. Yushmanov, V. E., P. K. Mandal, Z. Liu, P. Tang, and Y. Xu. 2003. NMR structure and backbone dynamics of the extended second transmembrane domain of the human neuronal glycine receptor 1 subunit. Biochemistry. 42:3989–3995.[CrossRef][Medline]
20. Humphrey, W., A. Dalke, and K. Schulten. 1996. VMD: visual molecular dynamics. J. Mol. Graph. 14:33–38.[CrossRef][Medline]
21. Unwin, N. 1995. Acetylcholine receptor channel imaged in the open state. Nature. 373:37–43.[CrossRef][Medline]
22. Rundstrom, N., V. Schmieden, H. Betz, J. Bromann, and D. Langosch. 1994. Cyanotriphenylborate: subtype-specific blocker of glycine receptor chloride channels. Proc. Natl. Acad. Sci. USA. 91:8950–8954.[Abstract/Free Full Text]
23. Smart, O. S., J. G. Neduvelil, X. Wang, B. A. Wallace, and M. S. P. Sansom. 1996. HOLE: a program for the analysis of the pore dimensions of ion channel structural models. J. Mol. Graph. 14:354–360.[CrossRef][Medline]
24. Shan, Q.,
TOP
ABSTRACT
INTRODUCTION
