A molecular dynamics study of an L-type calcium channel model

Gabriela Barreiro, Cristiano Ruch Werneck Guimarães and Ricardo Bicca de Alencastro,1

Physical Organic Chemistry Group, Departamento de Química Orgânica, Instituto de Química, UFRJ, Bloco A, sala 609, Cidade Universitária,Ilha do Fundão, CT, RJ 21949-900, Rio de Janeiro, Brazil


    Abstract
 Top
 Abstract
 Introduction
 Methods
 Results and discussion
 References
 
In this work, we propose a molecular model of the L-type calcium channel pore from the human cardiac {alpha}1 subunit. Four glutamic acid residues, the EEEE locus, located at highly conserved P loops (also called SS1–SS2 segments) of the {alpha}1 subunit, molecularly express the calcium channel selectivity. The proposed {alpha}-helix structure for the SS1 segment, analyzed through molecular dynamics simulations in aqueous-phase, was validated by the plotting of Ramachandran diagrams for the averaged structures and by the analysis of i and i + 4 helical hydrogen bonding between the amino acid residues. The results of the simulation of the calcium channel model with one and two Ca2+ ions at the binding site are in accordance with mutation studies which suggest that the EEEE locus in the L-type calcium channel must form a single high-affinity binding site. These results suggest that the Ca2+ permeation through the channel would be derived from competition between two ions for the only high-affinity binding site. Furthermore, the experimentally observed blocking of the Na+ flux at micromolar Ca2+ concentrations, probably due to the occupancy of the single high-affinity binding site for one Ca2+, was also reproduced by our model.

Keywords: {alpha}-helix structural analysis/Ca2+ binding site/calcium channel/molecular dynamics/molecular model


    Introduction
 Top
 Abstract
 Introduction
 Methods
 Results and discussion
 References
 
Voltage-gated calcium channels (VGCCs) are responsible for one of the communication mechanisms between cells. Ca2+ ions are major players in the intracellular signaling system that translates extracellular stimuli and regulates important phenomena such as muscle contraction, cell proliferation, gene expression, neurotransmitter release and other secretion processes, and cell death (Catterall, 1988Go; Imoto, 1993Go; Varadi et al., 1995Go, 1999Go).

One low (T type) and five high VGCC types (L, N, P, Q, R) have been studied through pharmacological and electrophysiological studies. Although only the pore-forming subunits of three members of T-type calcium channels have been cloned until now (Perez-Reyes, 1998Go; Perez-Reyes et al., 1998Go; Lee et al., 1999Go; Lacinova et al., 2000Go; Lory et al., 2000Go; McRory et al., 2001Go), the L-type subfamily has been characterized extensively by biochemical approaches. These studies revealed that the L-type calcium channel complex is a heteropentamer consisting of {alpha}1, ß, {alpha}2/{delta} and {gamma} subunits (Figure 1aGo). For calcium channels to be effective, Ca2+ ions must enter selectively through the pore of the {alpha}1 subunit, bypassing competition with other extracellular ions (Catterall, 1988Go; Imoto, 1993Go; Varadi et al., 1995Go, 1999Go; Randall and Benham, 1999Go). The predicted structure of the {alpha}1 subunit consists of four repeating motifs (MI–MIV), each motif comprising six hydrophobic segments (S1–S6). A highly conserved segment connecting the S5 and S6 transmembrane domains in each motif, termed the P loop or `SS1-SS2' region, is responsible for calcium selectivity in the pore region (Figure 1bGo) (Catterall, 1988Go; Varadi et al., 1999Go).



View larger version (30K):
[in this window]
[in a new window]
 
Fig. 1. (a) Molecular structure of the VGCC complex. Adapted from Varadi et al. (Varadi et al., 1999Go). (b) Proposed transmembrane arrangement of the VGCC {alpha}1-subunit. The dark and light gray cylinders are the SS1 and SS2 segments, respectively. Adapted from Catterall (Catterall, 1988Go). (c) Primary sequence of the pore-lining region. The EEEE locus is depicted in italic letters and the G homologous residues are also depicted in the dark gray vertical box.

 
Figure 1cGo shows the primary sequence of the P loop in each motif of the human cardiac {alpha}1 subunit. The SS1 segment in each motif is composed of amino acids predicted to be in {alpha}-helix arrangement (Varadi et al., 1995Go). The SS2 segment, on the other hand, is thought to form a different secondary structure, a random coil or a ß-strand structure (Varadi et al., 1999Go; Koch et al., 2000Go). The amino acids that connect the SS1 and SS2 segments are located in the SS2 region. They have been identified as MEGW, GEDW, FEGW and ATGE in motifs I–IV, respectively (Figure 1cGo) (Varadi et al., 1995Go, 1999Go). These residues are believed to form a ß-turn structure between the SS1 and SS2 segments (Varadi et al., 1999Go).

Four closely aligned glutamate residues, the EEEE locus, localized in the ß-turn region of the SS2 segment of each P loop segment molecularly express the calcium selectivity of the VGCCs. These negatively charged side chains form the Ca2+ selectivity filter (Yang et al., 1993Go; Varadi et al., 1995Go, 1999Go). Site-directed mutagenesis studies, in which the E residues were mutated to D, Q, A or K, have shown a considerable change into the monovalent and the divalent cation selectivity of the Ca2+ channel. These studies have also demonstrated that E residues in the SS2 regions collectively contribute to the high affinity binding of Ca2+ (Kim et al., 1993Go; Tang et al., 1993Go; Klöckner et al., 1996Go; Varadi et al., 1999Go). The E residues that form the ion selectivity filter are believed to be arranged in two close but non-equivalent planes, occupying trans positions (Varadi et al., 1999Go). On the other hand, in Na+ and K+ channels the selectivity filter residues are symmetrically arranged (Varadi et al., 1999Go).

Two different theories describe the calcium affinity-binding site. Whereas the first, called the multiple-binding site model, suggests that the EEEE locus is divided into two high-affinity calcium-binding sites (Almers and McCleskey, 1984Go; Hess and Tsien, 1984Go; McCleskey, 1999Go), the other suggests the existence of only one binding site of high-affinity (Yang et al., 1993Go; McCleskey, 1999Go). Yang et al. have shown that the neutralization of any of the E residues to a Q shifted the blocking of the monovalent current to a higher Ca2+ concentration (Yang et al., 1993Go). Since a single point mutation can shift the dissociation constant for Ca2+, the EEEE locus in L-type must form a single high-affinity binding site.

The three-dimensional structure of VGCCs is not available; thus, recent works aim to understand the highly selective mechanism of the calcium channel. Nonner et al. have recently examined the hypothesis that the selectivity of the L-type calcium channels is determined by mutual electrostatic screening and volume exclusion between the ions and the carboxylate oxygen atoms (Nonner et al., 2000Go). This study was able to predict Ba2+/Ca2+ and Na+/Li+ competition and Cl- exclusion. However, their selective filter model does not permit the understanding at the molecular level of the high calcium channel selectivity. Boda et al. have used an infinite cylinder model containing negatively charged glutamates, each modeled as a pair of half-charged oxygens, and Monte Carlo simulation to study the mechanism for calcium channels selectivity (Boda et al., 2000Go). Their results suggest that both the volume exclusion effects due to the narrow channel diameter, and the negative charge on the glutamates are important for the selectivity observed in the calcium channels. Nevertheless, an explanation of the calcium channel selectivity through the volume of the involved atoms still does not provide molecular information concerning the selectivity mechanism. In another work, Corry et al. have proposed a calcium channel model that can explain observed properties like the anomalous fraction effect and the mutation of the glutamate residues (Corry et al., 2000Go). They adopted a simple model, which resembles the structure of the KcsA bacterial potassium channel, and submitted it to electrostatic calculations and Brownian dynamics simulations to study the mechanism of ion permeation and selectivity in the channel. Their results have demonstrated that the behavior of the calcium channel can be explained by simple electrostatic interactions between ions, the channel boundary and the charges therein. However, again, explanations at the molecular level were not provided.

Doughty et al., on the other hand, have proposed two molecular models of the outer vestibule of the N-type VGCC pore (Doughty et al., 1995Go). Their models were constructed from ß-hairpin peptide segments of the S5–S6 loops of each of the four domains that form the channel. From these filter models, they derived the calcium channel dimensions necessary to allow the permeation of a hydrated calcium ion. In another work, Guy and Durrell have combined molecular modeling and mutagenesis experiments to determine the structure of Na+, Ca2+ and K+ channels (Guy and Durell, 1995Go). Among the approaches used, the most interesting, involving the VGCCs, was the use of a helix structure for the construction of the SS1 segment of the P loop.

Also, using the KcsA bacterial potassium channel crystal structure, Lipkind and Fozzard (Lipkind and Fozzard, 2001Go; the work of Lipkind and Fozzard has been published after the first submission of the present article) have proposed a molecular model of the outer vestibule and selectivity filter of the pore of the L-type calcium channel. Some structural and experimental evidences favor the idea that the KcsA channel is an evolutionary predecessor of sodium and calcium channels (Lipkind and Fozzard, 2000Go, 2001). Nevertheless, very important differences between potassium and calcium channels suggest that the selectivity filter of both channels are structurally dissimilar (Schetz and Anderson, 1993Go; Yang et al., 1993Go; Varadi et al., 1995Go, 1999Go) and, thus, that the KcsA channel may not be a good template to build a calcium channel model. Moreover, the L-type calcium channel model of Lipkind and Fozzard (Lipkind and Fozzard, 2001Go) neglects the experimental observation that the EEEE locus is arranged assymmetrically (Yang et al., 1993Go; Ellinor et al., 1995Go) and assumes a different amino acid sequence from that established in the literature (Varadi et al., 1995Go) to build the turn between the SS1 and SS2 segments.

The amino acids that connect the SS1 and SS2 segments are located in the SS2 region. They have been identified as MEGW, GEDW, FEGW and ATGE in motifs I, II, III and IV, respectively (Figure 1cGo) (Varadi et al., 1995Go, 1999Go). These residues are believed to form a ß-turn structure between the SS1 and SS2 segments (Varadi et al., 1999Go).

In this study, we propose a molecular model of the L-type calcium channel pore of the human cardiac {alpha}1 subunit. The model was constructed using the available experimental and theoretical results, and analyzed by molecular dynamics (MD) simulations. The MD simulations were used to study both the supposed {alpha}-helix structure for the SS1 segment in each motif, and the behavior of the channel model in the presence of one and two Ca2+ and also of one Ca2+ and one Na+ at the binding site. Our intention in doing so was to try to identify the presence of one or two high-affinity calcium-binding sites, and to investigate the mechanism of Ca2+ permeation and the blocking of Na+ current at higher Ca2+ concentrations.


    Methods
 Top
 Abstract
 Introduction
 Methods
 Results and discussion
 References
 
The main interest of this work was to study at the molecular level the selectivity filter of the L-type calcium channel, formed by four glutamates localized at the SS2 segment in each motif. The L-type calcium channel model was constructed from the primary sequence of the SS1 and SS2 segments (Heinemann et al., 1994Go; Varadi et al., 1995Go, 1999Go) shown in Figure 1cGo.

The first part of the model construction consisted of the election of the amino acids that would be important to reproduce the environment of the pore in the calcium channel. Thus, we considered all amino acid residues of the SS1 segment, and six residues of the SS2 segment: four which form the ß-turn between the SS1 and SS2 segments, plus two beyond that (MI: NFDNFAFAMLTVFQCITMEGWTD; MII: RSTFDNFPQSLLTVFQILTGEDWNS; MIII: DFDNVLAAMMALFTVSTFEGWPE; MIV: NNNFQTFPQAVLLLFRC- ATGEAW).

For the model construction, we initially evaluated the stability of the four SS1 segment {alpha}-helix structures in the aqueous-phase. The SS1 segment structures of all motifs were separately submitted to a MD simulation of 150 ps (50 ps for the equilibration stage and 100 ps for the data collection). Acetyl and N-methylamine blocking groups were used to cap the truncation points in these simulations. Each of the {alpha}-helix canonical structures ({phi}1 = –57° and {phi}1 = –47°) (Voet and Voet, 1995Go) were placed in a rectangle box with ~1525 water molecules, applying periodic boundary conditions. The stability of the four SS1 segment {alpha}-helix structures in the aqueous-phase was then analyzed in two parts. The first part was done by plotting Ramachandran diagrams for the averaged structures in solution, obtained from the PROCHECK program available on the Internet (Laskowski et al., 1993Go). The second part was performed by the analysis of helical hydrogen bonding (hHB) between the amino acid residues i and i + 4 (Vijayakumar et al., 1994Go), and pair-radial distribution functions between the oxygen of the water molecules and the carbonyl oxygen or the amide nitrogen in the (C=O)i···(H–N)i + 4 hydrogen bond. These analyses were done in order to evaluate the {alpha}-helix inherent structural weaknesses and the eventual water molecule insertions in the hHB of the canonical structure.

After the analysis of the {alpha}-helix stability, the pore of the channel was constructed. To accomplish that, the SS1 segment was built following the {alpha}-helix canonical dihedral angles. As mentioned before, the SS2 segment is thought to form a random coil or a ß-strand structure (Varadi et al., 1999Go; Koch et al., 2000Go). Since a random coil is a totally disordered and rapidly fluctuating set of conformations assumed by proteins, it is impossible to determine the dihedral values for the two residues of the SS2 segment beyond the turn region. Hence, these residues were assumed to be in a ß-strand conformation ({phi}1 = –119° and {phi}1 = 113°) (Voet and Voet, 1995Go). The ß-turn region, that connects the SS1 and SS2 segments was built following the dihedral angles {phi}2 = –60°, {phi}2 = 120°, {phi}3 = 90° and {phi}3 = 0° (Voet and Voet, 1995Go). Here, acetyl and N-methylamine blocking groups were also used to cap the truncation points.

Yang et al. (Yang et al., 1993Go) and Ellinor et al. (Ellinor et al., 1995Go) have verified in mutagenesis studies the distinct contributions of each of the E residues to high-affinity divalent cation binding. They have concluded that the pore would not be arranged in a symmetrical ring of four strictly equivalent E residues and proposed that the non-equivalence of these amino acids might arise from asymmetry in the positions of their {alpha}-carbon atoms. To explain this non-equivalence, Varadi et al. (Varadi et al., 1999Go), as mentioned above, have suggested that the E residues that form the ion selectivity filter should be arranged in two close but non-equivalent planes, occupying trans positions (Varadi et al., 1999Go). Whereas one of the planes would be formed by the E residues of motifs I and III, the other would be formed by the E residues of motifs II and IV. One alternative to accomplish that in the elaboration of a calcium channel model is to align the G residues of each motif (Figure 1cGo) (Varadi et al., 1999Go).

The four motif models were docked in a clockwise orientation viewed from the outside (Randall and Benham, 1999Go; Varadi et al., 1999Go) keeping the value of the narrowest diameter of the pore (~6 Å) (McCleskey, 1999Go) and analyzing the VDW and Coulomb energies. The calcium ion, which we used as a reference in the docking procedure, guided us in the orientation of the SS1 helices about their axis. To align the G residues (Figure 1cGo) (Varadi et al., 1999Go) of each motif model, we introduced pseudo-atoms at the center of mass of each G and defined a dihedral angle, the value of which was kept as close as possible to zero when the four motifs were docked.

To avoid undesirable interatomic contacts after the docking procedure and to allow the relaxation of the system, the side chains of all amino acids were energy minimized to a gradient norm of 0.01 kcal/mol/Å. We kept all backbone atoms fixed and submitted the system to a MD run of 100 ps (50 ps for the equilibration stage and 50 ps for the data collection stage), and energy minimization calculation to a gradient norm of 0.01 kcal/mol/Å.

The energy minimized structure was placed in a rectangle box with ~1550 water molecules, without periodic boundary conditions. To avoid evaporation, a 5 Å layer of water molecules surrounding the whole system was kept fixed during all simulations. In addition, to avoid anomalous behavior due to the absence of the discarded amino acids, the methyl carbon atoms of the acetyl and N-methylamine blocking groups were also kept fixed. The next step was the energy minimization of the free water molecules to a gradient norm of 0.01 kcal/mol/Å. In this procedure, the calcium ion and all amino acids of the system were kept fixed. After minimization, the free water molecules were submitted to a 20 ps MD simulation (10 ps for the equilibration stage and 10 ps for the collection data stage), a longer time than the relaxation time of water. Here, all other atoms of the system were kept fixed in order to let the water molecules adjust to the potential field of the solute. Finally, the complete system was energy minimized to a gradient norm of 0.01 kcal/mol/Å.

For the evaluation of the coordination of one calcium ion with the oxygen atoms of the EEEE locus, the system was submitted to a 200 ps MD simulation (100 ps for the equilibration stage and 100 ps for the data collection stage). In the last configuration of this MD run, a second ion (Ca2+ or Na+) was placed at the binding site, ~3.2 Å above the first ion. Keeping the ions fixed, these models were extensively relaxed through energy minimization and MD simulations. To analyze subsequently the state of coordination of the ions (two Ca2+, and one Ca2+ plus one Na+) in the channel model, the equilibrated configuration was first energy-minimized to a gradient norm of 0.01 kcal/mol/Å and then submitted to a 100 ps MD simulation (40 ps for the equilibration stage and 60 ps for the data collection stage), without constraints on the ions.

For the energy minimization and MD simulation steps we used the all-atom CVFF force field (Dauber-Osguthorpe et al., 1988Go) within the Discover program of Molecular Simulation Inc. To save computational time, spherical residue-based non-bonded interactions cut-offs of 10 and 15 Å were applied to the {alpha}-helix stability analysis and to the ion coordination studies, respectively. To smoothly `turn off' the interactions, we employed a quintic spline function from 8.5 to 10 Å and 13.5 to 15 Å, minimizing, in this manner, discontinuities in the potential energy surface. Whenever an atom moved more than half the length of the buffer region (between 10 and 11 Å and 15 and 16 Å, in each case), the neighbor list was updated. This ensured that no atoms outside the buffer region were able to move to positions close enough to interact. All MD simulations were performed at 310 K using the velocity scaling (Woodock, 1971Go) and the Berendsen (Berendsen et al., 1984Go) methods to control the temperature during the equilibration and data collection stages, respectively. The equations of motion were integrated every 1.0 fs using the Verlet Leapfrog algorithm (Hockney, 1970Go). The trajectory was sampled every 1 ps. All calculations were performed on a Silicon Graphics O2 R10000 workstation under the Irix® operational system.


    Results and discussion
 Top
 Abstract
 Introduction
 Methods
 Results and discussion
 References
 
Analysis of the dynamical stability of the proposed {alpha}-helix structure for the SS1 segments

The averaged {alpha}-helix structures in the aqueous phase for the SS1 segments, as obtained from the MD simulations, were first analyzed through their Ramachandran plots, using the PROCHECK program (Laskowski et al., 1993Go). The diagrams show that the great majority of the amino acids are localized in the most favored regions for an {alpha}-helix structure (motif I = 93.3%, motif II = 100%, motif III = 93.3%, motif IV = 100%) (Figure 2Go). Table IGo gives all averaged backbone dihedral angles for each motif. As shown by Table IGo, the averaged values are mainly close to the canonical {alpha}-helix dihedral angles ({phi} = –57° and {phi} = –47°).



View larger version (62K):
[in this window]
[in a new window]
 
Fig. 2. Ramachandran plots of the SS1 segments of motifs I–IV.

 

View this table:
[in this window]
[in a new window]
 
Table I. Averaged backbone dihedral angles for motifs I–IV
 
As suggested by Vijayakumar et al. (Vijayakumar et al., 1994Go), the analysis of the hHB between the carbonyl oxygen and the amide nitrogen of residues i···i + 4 is an important method to study the dynamical stability of secondary structures in proteins. As discussed by DiCapua et al. (DiCapua et al., 1990Go, 1991Go) hydrated {alpha}-helices may have an inherent structure weakness, which propitiates water insertion into the hHB. This structure weakness, which is probably manifested in the persistence length, the local geometry or the local electrostatic environment, may be seen through the hHB associated with these regions, which will be less strong and will oscillate with greater frequency and amplitude. The obtained hHB graphics for all four SS1 structures have shown that some amino acids have not kept this interaction during the time spanned by the simulation, probably due to an inherent structure weakness. Nevertheless, as shown below, the integrity of the {alpha}-helix structure was maintained for all SS1 segments because only a reduced number of hHB broke during the simulation.

In this work, the residue numbering does not follow the primary sequence numbering of the P loop (Koch et al., 2000Go). The motifs I, II and IV were numbered from 1 to 23, and the motif II from 1 to 25. Following the diagrams i···i + 4 of the hHB interaction (Brooks and Case, 1993Go) (Figure 3Go), it is possible to observe strong hHB (solid arrows), hHB that oscillates with greater frequency and amplitude (dashed arrows) and broken hHB (absence of an arrow).



View larger version (42K):
[in this window]
[in a new window]
 
Fig. 3. Backbone stability [(i···i + 4) hHB interaction (Brooks and Case, 1993Go)] of motifs I–IV. The solid arrows show the strong hHB interactions between the amino acids, the dashed arrows show the hHB interactions that oscillate with greater frequency and amplitude, and the absence of an arrow shows the broken hHB interactions during the MD simulations.

 
The analysis of the hydrogen bond stability reveals an oscillatory behavior or a hHB breaking at the terminal regions of all motifs. The weakness at these regions is probably derived from end effects (DiCapua et al., 1990Go, 1991Go). The hHB breaking at the terminal regions of motif I is graphically represented by the time evolution of the distance between the carbonyl oxygen and the amide nitrogen in the (C=O)i···(H–N)i + 4 hydrogen bond for the pair of residues N-1/F-5, F-13/T-17 and V-12/I-16 (Figure 4aGo, graphs 1, 3 and 5). The value of this distance in the canonical {alpha}-helix is ~3.0 Å (DiCapua et al., 1990Go, 1991Go). Thus, due to this inherent structure weakness, we have observed the insertion of water molecules between the (C=O)i···(H–N)i + 4 hydrogen bond for the pair N-1/F-5. This insertion is characterized by the first peak of the pair-radial distribution function for the carbonyl oxygen of N-1 and the oxygen atoms of the water molecules (Figure 4aGo, graph 2). Graphs 3 and 5 (Figure 4aGo) show that as soon as the simulation commences the (C=O)i···(H–N)i + 4 hydrogen bonds between F-13/T-17 and V12/I-16 break. As a result, a new i···i + 5 hydrogen bond between V-12/T-17 is formed (Figures 3 and 4aGoGo, graph 6). In this case, there is no water insertion between the (C=O)i···(H–N)i + 4 hydrogen bond for the pair F-13/T-17, as may be seen by the pair-radial distribution function for the carbonyl oxygen F-13 and the oxygen atoms of the water molecules (Figure 4aGo, graph 4).



View larger version (40K):
[in this window]
[in a new window]
 
Fig. 4. Time evolution of the broken hHB interactions [(C=O)i···(H–N)i + 4 interaction] for some pairs of residues; and of the pair-radial distribution functions between the oxygen atoms of the backbone carbonyl group and of the water molecules of motifs I (a) and II (b), respectively.

 
A similar behavior is observed for motif II. Graphs 7 and 8 (Figure 4bGo) show that the hHB between V-14/L-18 is replaced by a i···i + 5 hydrogen bond between V-14/T-19 (Figure 3Go). Figure 4bGo also shows that the hHB between N-6/S-10 breaks (graph 11) because of the influence of the side-chain hydroxyl group of S-10 (graph 12), which competes for the hydrogen bond with the amide nitrogen. We believe that the hHB breaking for the pair D-5/Q-9 (Figure 4bGo, graph 9) is probably derived from destabilization caused by the proximity of the proline residue (P-8). This is an example of a weakness caused by the presence of a strong {alpha}-helix breaker. This weakness propitiates the insertion of a water molecule in the hHB of D-5/Q-9, as characterized by the first peak of the pair-radial distribution function for the backbone carbonyl oxygen D-5 and the oxygen atoms of the water molecules (Figure 4bGo, graph 10).

The hHB breaking at the N-terminal region of motif III is graphically represented by the time evolution of the distance between the carbonyl oxygen and the amide nitrogen in the (C=O)i···(H–N)i + 4 hydrogen bond for the pair of residues D-1/V-5 (Figure 5aGo, graph 13). As discussed above, this disruption is derived from a weakness due to end effects, which propitiates the insertion of a water molecule in the hHB between D-1/V-5. It is demonstrated by the first peak of the pair-radial distribution function for the backbone carbonyl oxygen of D-1 and the oxygen atoms of the water molecules (Figure 5aGo, graph 14). As with the pair N-6/S-10 of motif II, the hHB between L-12/S-16 breaks because of the competition between the amide nitrogen and the side-chain hydroxyl group of S-16 for the hydrogen bond with the carbonyl oxygen of L-12 (Figure 5aGo, graphs 17 and 18). We believe that the HB interaction between the side-chain hydroxyl group of S-16 and the carbonyl oxygen of L-12 causes local destabilization of the {alpha}-helix structure, leading to the hHB breaking of the pair A-11/V-15 (Figure 5aGo, graph 15). As a result, a new i···i + 5 hydrogen bond between M-10/V-15 is formed (Figures 3 and 5aGoGo, graph 16).



View larger version (33K):
[in this window]
[in a new window]
 
Fig. 5. Time evolution of the broken hHB interactions [(C=O)i···(H–N)i + 4 interaction] for some pairs of residues; and of the pair-radial distribution functions between the oxygen atoms of the backbone carbonyl group and of the water molecules of motifs III (a) and IV (b), respectively.

 
The hHB breaking at the N-terminal region of motif IV is graphically represented by the time evolution of the distance between the carbonyl oxygen and the amide nitrogen in the (C=O)i···(H–N)i + 4 hydrogen bond for the pair of residues N-1/Q-5 (Figure 5bGo, graph 19). As this hHB breaks only at the end of the simulation, the first peak of the pair-radial distribution function for the carbonyl oxygen of N-1 and the oxygen atoms of the water molecules, which characterizes a water insertion, is not so significant (Figure 5bGo, graph 20). Graph 21 (Figure 5bGo) shows the hHB between the pair Q-5/Q-9. In this case, the proximity of the proline residue (P-8) does not lead to the breaking of this interaction, which only oscillates with greater amplitude and frequency (Figure 3Go). As a result of this, water molecules are not inserted in the hHB between the pair Q-5/Q-9. This is demonstrated by the reduced first peak of the pair-radial distribution function for the carbonyl oxygen of Q-5 and the oxygen atoms of the water molecules (Figure 5bGo, graph 22). As did DiCapua et al. (DiCapua et al., 1991Go), we have observed that, once inserted, the water molecules have a great tendency to bind at the carbonyl groups of the backbone (see graphs of pair-radial distribution functions, Figures 4 and 5GoGo). The results have shown that in the time span of the simulations, the {alpha}-helix secondary structure for the SS1 segment of each motif is stable in the aqueous phase. Now, if this structure is stable in this circumstance, the {alpha}-helix will be even more stable in a less hydrated environment, as in the typical environment of the L-type calcium channel.

Analysis of the behavior of the L-type calcium channel model in the presence of one and two Ca2+, and in the presence of one Ca2+ and one Na+ at the binding site

It has been observed experimentally (Yang et al., 1993Go) that at submicromolar concentrations of Ca2+, the Na+ ion can permeate freely through the calcium channel. At micromolar Ca2+ concentrations the Na+ flux is blocked, and at higher Ca2+ concentrations the Ca2+ flux is established. Two different theories describe the calcium movement through the channel. The multiple-binding site theory suggests the presence of two high-affinity calcium-binding sites (Almers and McCleskey, 1984Go; Hess and Tsien, 1984Go; McCleskey, 1999Go). In this theory, at micromolar Ca2+ concentrations, foreign ions, like the Na+ ion, are blocked, and a Ca2+ ion occupies one of the two high affinity sites. At higher Ca2+ concentrations, a second ion occupies the other site, and the electrostatic repulsion between them leads to the Ca2+ flux. On the other hand, the second theory suggests the existence of only one single high-affinity binding site (Yang et al., 1993Go; McCleskey, 1999Go) that is occupied by one Ca2+ ion at micromolar Ca2+ concentrations, blocking the Na+ flux. At higher concentrations (above millimolar), two Ca2+ ions compete with each other for the E residues of the single high-affinity binding site; the competition decreases the affinity of the pore for Ca2+, propitiating the Ca2+ permeation. As discussed above, Yang et al. (Yang et al., 1993Go) have shown that the neutralization of any of the E residues to Q shifted the blocking of the monovalent current to higher Ca2+ concentrations. Since a single point mutation can shift the dissociation constant for Ca2+, the EEEE locus in the L-type calcium channel must form a single high-affinity binding site.

In the present work, the model with one Ca2+ ion at the binding site (Figure 6a and bGo) was submitted to 200 ps of a MD simulation. As discussed by Ellinor et al. (Ellinor et al., 1995Go), during the simulation, we observed the coordination of the calcium ion with seven carboxylate oxygen atoms of the EEEE locus; i.e. the carboxylate oxygens of the E residues from motifs I, II and IV, and only one of the E residues from motif III. This is illustrated by the last configuration of the MD run (Figure 6cGo) and the time evolution of the distances between the calcium ion and the carboxylate oxygen atoms of the EEEE locus (the ideal coordination distance is ~2.4 Å) (Åqvist, 1992Go; Katz et al., 1996Go) (Figure 7Go). In accordance with the assumption of Varadi et al. (Varadi et al., 1999Go), we observed that even after the MD simulations, the EEEE locus is still arranged on two close, but non-equivalent planes (Figure 8aGo) occupying trans positions. Whereas one plane is formed by the E19 residues from motifs I and III, the other is formed by the E21 residues from motifs II and IV. Figure 8aGo also shows that all E residues are important to the Ca2+ binding in the pore. This result agrees with the mutation studies of Yang et al. (Yang et al., 1993Go) discussed above. In spite of the presence of several water molecules around Ca2+, the ion is coordinated only with the carboxylate oxygens during the time spanned by the simulation. The oxygen O2 of the E19 residue of motif III, which is the only one that does not interact with the ion, is stabilized by hydrogen bond interactions with water molecules and with the indol nitrogen of W23 of motif IV (Figure 6cGo). The W23 residue of motif II also seems to stabilize the binding site through hydrogen bond interactions between its indol nitrogen and the carboxylate oxygen of the E21 residue of the same motif (Figure 6cGo). Since the net charge of the binding site is –2, the hydrogen bond interactions between the glutamate and the tryptophan residues provide additional stabilization to it.



View larger version (131K):
[in this window]
[in a new window]
 
Fig. 6. The VGCC pore model with one Ca2+ ion at the binding site. (a) A vertical vision. (b) A horizontal vision from the bottom of the channel. (c) A close vision of the interaction between the Ca2+ ion and the carboxylate oxygen atoms of motifs I–IV and also of the presence of the W residues of motifs II and IV.

 


View larger version (28K):
[in this window]
[in a new window]
 
Fig. 7. Time evolution of the distances between the calcium ion and the carboxylate oxygen atoms of the EEEE locus (ideal coordination distance ~2.5 Å) (Åqvist, 1992; Katz et al., 1996).

 


View larger version (94K):
[in this window]
[in a new window]
 
Fig. 8. (a) The arrangement of the EEEE locus in two close but non-equivalent planes occupying trans positions. The E19 residues from motifs I and III form the upper plane (b) and the E21 residues from motifs II and IV form the lower plane (c).

 
As shown above, our results strongly suggest that the EEEE locus forms a single high-affinity binding site. If this is the case, the Ca2+ permeation at higher concentrations must take place when two Ca2+ ions compete with each other for the E residues. To investigate this assumption, we placed a second ion (CaB2+) in a position ~3.2 Å above the first ion (CaA2+), and, keeping the ions fixed, extensively relaxed the model through energy minimization and MD simulations. Figure 8bGo (left), shows the equilibrated configuration with the ions fixed to the starting positions. Whereas CaB2+ is coordinated to the carboxylate oxygens of the E19 residues from motifs I and III, CaA2+ is coordinated to the carboxylate oxygens of the E21 residues from motifs II and IV. Moreover, two water molecules are bonded to each calcium ion. This structure was then energy-minimized and submitted to a 100 ps MD simulation with no constraints in the Ca2+ ions. Figure 9aGo shows that after the energy minimization step (the initial configuration of the MD simulation) the distance between the Ca2+ ions increased to ~4 Å, as a result of the electrostatic repulsion. In this configuration, CaA2+ is now coordinated to the carboxylate oxygens of the E21 residue from motif IV, and to three other carboxylate oxygens: one from the E19 residue from motif I, one from the E21 residue from motif II, and finally, one from the E19 residue from motif III (Figure 9bGo). CaB2+, on the other hand, is coordinated to five carboxylate oxygens: four from E residues from motifs I and II, and one from the E19 residue from motif III (Figure 9cGo). As the simulation starts, CaA2+ is displaced from the binding site, and, at around 65 ps, all the coordination with the EEEE locus is almost lost, except for the coordination with one carboxylate oxygen of the E21 residue from motif IV (Figure 9bGo). As a result, the distance between the ions increases to ~5.5 Å (Figure 9aGo). Conversely, we observed that during the MD simulation CaB2+ keeps its interactions with the binding site (Figure 9cGo). The last configuration, as illustrated in Figure 8bGo (right), shows that whereas the displaced Ca2+ is heptacoordinated to one carboxylate oxygen and six water molecules, CaB2+ is octacoordinated to five carboxylate oxygen and three water molecules. Comparison of Figures 8a and bGo (right) seems to indicate that as CaA2+ is displaced from the EEEE locus, the coordination of one ion with the single high-affinity binding site starts to be re-established. Our results suggest that the Ca2+ permeation through the channel would be indeed derived from competition between two ions for the single high-affinity binding site.



View larger version (36K):
[in this window]
[in a new window]
 
Fig. 9. (a) Time evolution of the distance between the CaA2+ and CaB2+ ions after the energy minimization step (the initial configuration of the MD simulation). (b and c) Time evolution of the distances between the CaA2+ and CaB2+ ions and the carboxylate oxygen atoms of the EEEE locus, respectively. 120

 
To investigate if our model reproduces the experimentally observed blocking of Na+ flux when the single high-affinity binding site is occupied by one Ca2+ ion, we placed a Na+ ion at a position situated ~3.2 Å above the Ca2+ ion. As described, this model was extensively relaxed through energy minimization and MD simulations keeping both ions fixed. The equilibrated configuration with the Ca2+ and the Na+ ions fixed to the starting positions is shown in Figure 8cGo (left). This structure was energy-minimized and submitted to a 100 ps MD simulation with no constraints for the ions. Figure 10aGo shows that after the energy minimization step (the initial configuration of the MD simulation) the distance between the ions increased to ~3.7 Å. In this configuration, we observed that the Ca2+ ion is coordinated to five carboxylate oxygens: one from the E19 residue from motif I, one from the E21 residue from motif II, one from the E19 residue from motif III, and both the carboxylate oxygens from the E21 residue from motif IV (Figure 10bGo). The Na+ ion is also coordinated to five carboxylate oxygens (the ideal coordination distance is also ~2.4 Å) (Åqvist, 1992Go; Katz et al., 1996Go), two from the E19 residue from motif I, two from the E21 residue from motif II, and one carboxylate oxygen from the E19 residue from motif III (Figure 10cGo). After the simulation, the last configuration (Figure 8cGo, right) shows that the Ca2+ ion is octacoordinated to five carboxylate oxygens (Figure 10bGo) and three water molecules, whereas the Na+ ion is heptacoordinated to five carboxylate oxygens (Figure 10cGo) and two water molecules. The important point here is that, differently from what we observed when two Ca2+ ions occupy the binding site, Figure 10aGo shows that the distance between Ca2+ and Na+ does not change during the simulation, indicating that the monovalent ion is not able to displace the Ca2+ ion from its binding site. This last result is in agreement with the experimental observations of the blocking of the Na+ flux at micromolar Ca2+ concentrations.



View larger version (33K):
[in this window]
[in a new window]
 
Fig. 10. (a) Time evolution of the distance between Ca2+ and Na+ ions after the energy minimization step (the initial configuration of the MD simulation). (b and c) Time evolution of the distances between Ca2+ and Na+ ions and the carboxylate oxygen atoms of the EEEE locus, respectively.

 
Conclusions

We applied MD simulations to evaluate the three-dimensional structure and behavior of a molecular model of the L-type calcium channel pore.

In the first part of this work, we validated the averaged {alpha}-helix structures for the SS1 segments obtained from the MD simulations in the aqueous phase. We used Ramachandran plots, which showed that practically all amino acids were localized in the most favored regions for an {alpha}-helix structure. The analysis of the hHB between the carbonyl oxygen and the amide nitrogen of residues i···i + 4 allowed an evaluation of the degree of the {alpha}-helix stability. Because the obtained {alpha}-helix structures are very stable in the aqueous phase, it is reasonable to infer that these particular structures will be even more stable in a less hydrated environment, as is the case of the L-type calcium channel, a membrane protein.

The behavior of the calcium channel model with one Ca2+ ion at the binding site is in agreement with the results from mutation studies of Yang et al. (Yang et al., 1993Go), which suggests that the EEEE locus in the L-type calcium channel must form a single high-affinity binding site. We were able to confirm this result through the analysis of an MD simulation of the pore model in the presence of two Ca2+ ions at the binding site. We observed that at the end of the MD run the first Ca2+ ion is displaced by the second one from the EEEE locus and that the coordination of the second ion with the single high-affinity binding site starts to dominate the process.

Differently from what we observed with two calcium ions, Ca2+ was not displaced from its binding site by the Na+. This last result is in agreement with the experimentally observed blocking of Na+ flux at micromolar Ca2+ concentrations.

In conclusion, our results suggest that the glutamate residues form a single high-affinity binding site that coordinates with only one Ca2+ ion, and that the permeation through the channel is indeed derived from competition between two ions for the EEEE locus. On the other hand, Na+ is blocked because it cannot electrostatically displace the Ca2+ ion from the single high-affinity binding site.


    Notes
 
1 E-mail: bicca{at}iq.ufrj.br Back


    Acknowledgments
 
This research has received partial financial support from Brazilian agencies CNPq and FAPERJ (grants nos E26/151494/99 and E26/151898/2000). G.B. would like to acknowledge CAPES/Brazil and C.R.W.G. would like to acknowledge FAPERJ/Brazil for scholarships.


    References
 Top
 Abstract
 Introduction
 Methods
 Results and discussion
 References
 
Almers,W. and McCleskey,E.W. (1984) J. Physiol., 153, 585–608.

Åqvist,J. (1992) J. Mol. Struct. (THEOCHEM), 256, 135–152.

Berendsen,H.C.J., Postma,J.P.M., van Gunsteren,W.F., DiNola,A. and Haak,J.R. (1984) J. Chem. Phys., 81, 3684–3690.

Boda,D., Busath,D.D., Henderson,D. and Sokotowski,S. (2000) J. Phys. Chem. B, 104, 8903–8910.[CrossRef][ISI]

Brooks,C.L.,III and Case,D.A. (1993) Chem. Rev., 93, 2847–2502.

Catterall,W.A. (1988) Science, 242, 50–61.[ISI][Medline]

Corry,B., Allen,T.W., Kuyucak,S. and Chung,S.-Ho (2000) Biochim. Biophys. Acta, 1509, 1–6.[ISI][Medline]

Dauber-Osguthorpe,P., Roberts,V.A., Osguthorpe,D.J., Wolff,J., Genest,M., and Hagler,A.T. (1988) Proteins: Struct. Funct. Genet., 4, 31–47.[ISI][Medline]

DiCapua,F.M., Swaminathan,S. and Beveridge,D.L. (1990) J. Am. Chem. Soc., 112, 6768–6771.[ISI]

DiCapua,F.M., Swaminathan,S. and Beveridge,D.L. (1991) J. Am. Chem. Soc., 113, 6145–6155.[ISI]

Doughty,S.W., Blaney,F.E. and Richards,W.G. (1995) J. Mol. Graphics, 13, 342–348.[CrossRef][ISI][Medline]

Ellinor,P.T., Yang,J., Sather,W.A., Zhang,J.-F. and Tsien,R.W. (1995) Neuron, 15, 1121–1132.[ISI][Medline]

Guy,H.R. and Durell,S.R. (1995) J. Gen. Physiol. Ser., 48, 1–16.

Heinemann,S.H., Schlief,T., Mori,Y. and Imoto,K. (1994) Braz. J. Med. Biol. Res., 27, 2781–2802.[ISI][Medline]

Hess,P. and Tsien,R.W. (1984) Nature, 309, 453–456.[ISI][Medline]

Hockney,R.W. (1970) Methods Comp. Phys., 9, 136–211.

Imoto,K. (1993) FEBS Lett., 325, 100–103.[CrossRef][ISI][Medline]

Katz,K.A., Glusker,J.P., Beebe,S.A. and Bock,C.W. (1996) J. Am. Chem. Soc., 118, 5752–5763.[CrossRef][ISI]

Kim,M-S., Morii,T., Sun,L.-X., Imoto,K. and Mori,Y. (1993) FEBS Lett., 318, 145–148.[CrossRef][ISI][Medline]

Klöckner,U., Mikala,G., Schwartz,A. and Varadi,G. (1996) J. Biol. Chem., 271, 22293–22296.[Abstract/Free Full Text]

Koch,S.E., Bodi,I., Schwartz,A. and Varadi,G. (2000) J. Biol. Chem., 275, 34493–34500.[Abstract/Free Full Text]

Lacinova,L., Klugbauer,N. and Hofmann,F. (2000) Gen. Physiol. Biophys., 19, 121–136.[ISI][Medline]

Laskowski,R.A., McArthur,M.W., Moss,D.S. and Thorton,J.M. (1993) J. Appl. Crystallogr., 26, 283–291.[CrossRef][ISI]

Lee,J.H., Daud,A.N., Cribbs,L.L., Lacerda,A.E., Pereverzev,A., Klockner,U., Schneider,T. and Perez-Reyes,E. (1999) J. Neurosci., 19, 1912–1921.[Abstract/Free Full Text]

Lipkind,G.M. and Fozzard,H.A. (2000) Biochemistry, 39, 8161–8170.[CrossRef][ISI][Medline]

Lipkind,G.M. and Fozzard,H.A. (2001) Biochemistry, 40, 6786–6794.[CrossRef][ISI][Medline]

Lory,P., Monteil,A., Chemin,J., Leuranguer,V., Bourinet,E. and Nargeot,J. (2000) Therapie, 55, 249–254.[ISI][Medline]

McCleskey,E.W. (1999) J. Gen. Physiol., 113, 765–772.[Free Full Text]

McRory,J.E., Santi,C.M., Hamming,K.S.C., Mezeyova,J., Sutton,K.G., Baillie,D.L., Stea,A. and Snutch,T.P. (2001) J. Biol. Chem., 276, 3999–4011.[Abstract/Free Full Text]

Nonner,W., Catacuzzeno,L. and Eisenberg,B. (2000) Biophys. J., 79, 1976–1992.[Abstract/Free Full Text]

Perez-Reyes,E. (1998) J. Bioenerg. Biomembr., 30, 313–318.[CrossRef][ISI][Medline]

Perez-Reyes,E., Cribbs,L.L., Daud,A., Lacerda,A.E., Barclay,J., Williamson,M.P., Fox,M., Rees,M. and Lee,J.H. (1998) Nature, 391, 896–900.[CrossRef][ISI][Medline]

Randall,A. and Benham,C.D. (1999) Mol. Cell. Neurosci., 14, 255–272.[CrossRef][ISI][Medline]

Schetz,J.A. and Anderson,P.A.V. (1993) Biol. Bull., 185, 462–466.[Abstract/Free Full Text]

Tang,S., Mikala,G., Bahinski,A., Yatani,A., Varadi,G. and Schwartz,A. (1993) J. Biol. Chem., 268, 13026–13029.[Abstract/Free Full Text]

Varadi,G., Mori,Y., Mikala,G., and Schwartz,A. (1995) Trends Pharmacol. Sci., 16, 43–49.[CrossRef][ISI][Medline]

Varadi,G., Strobeck,M., Koch,S., Caglioti,L., Zucchi,C. and Palyi,G. (1999) Crit. Rev. Biochem. Mol. Biol., 34, 181–214.[Abstract/Free Full Text]

Vijayakumar,S., Vishveshwara,S., Ravishanker,G. and Beveridge,D.L. (1994) Analysis of hydrogen bonding and stability of protein secondary structures in molecular dynamics simulation. In Smith,D.A. (ed), Modeling the Hydrogen Bond. American Chemical Society, Washington, DC, pp. 175–193.

Voet,D. and Voet,J.D. (1995) Three-dimensional structures of proteins. In Biochemistry. John Wiley & Sons, New York, pp. 141–152.

Woodock,L.V. (1971) Chem. Phys. Lett., 10, 257–261.[CrossRef][ISI]

Yang,J., Ellinor,P.T., Sather,W.A., Zhang,J.-F. and Tsien, R. W. (1993) Nature, 366, 158–161.[CrossRef][ISI][Medline]

Received March 16, 2001; revised July 31, 2001; accepted September 15, 2001.