Ca2+ Selectivity of the Sarcoplasmic Reticulum Ca2+−ATPase at the

12282

J. Phys. Chem. B 2007, 111, 12282-12293

Ca2+ Selectivity of the Sarcoplasmic Reticulum Ca2+-ATPase at the Enzyme-Water Interface and in the Ca2+ Entrance Channel Feng Xiang,† Robert I. Cukier,‡ and Yuxiang Bu*,†,‡ Key Laboratory for Colloid and Interface Chemistry of Ministry of Education, The Modeling & Simulation Chemistry DiVision, Shandong UniVersity, Jinan 250100, P. R. China and Department of Chemistry, Michigan State UniVersity, East Lansing, Michigan 48823 ReceiVed: May 19, 2007; In Final Form: August 1, 2007

The sarcoplasmic reticulum (SR) Ca2+-ATPase, a P-type transmembrane protein, can transport Ca2+ from the cytoplasmic to the luminal side over other cations specifically. The proposed Ca2+ entrance channel, composed of the main-chain carbonyl oxygen and side-chain carboxyl oxygen atoms of the amino acids, opens on the enzyme surface, just above the biphospholipid layer membrane-water interface, where Trp residues are frequently found. In this work, the physicochemical nature of Ca2+ selectivity over Mg2+ on the surface of the SR Ca2+-ATPase has been investigated using the density functional theory (DFT) method. The selection process can be regarded as the first step of the specificity of the enzyme to transport Ca2+. Subsequently, the specificity of the entrance channel to conduct Ca2+ over other cations has also been explored. As revealed by thermodynamic analyses, either the aromatic or the aliphatic amino acid residues distributed on the surface of Ca2+-ATPase have a bigger affinity to Mg2+ than to Ca2+, resulting in a concentration decrease of free Mg2+ in the local region. Thus, Ca2+ can transport into the Ca2+-entrance channel more easily. Whereafter, for a small quantity of Mg2+ entering this channel accompanying the Ca2+ current, the strong electrostatic interactions between Mg2+ and the ligands will limit the activity of this metal ion, which facilitates the weakly bonded Ca2+ passing through the channel at a relatively high rate, as suggested by the “sticky-pore” hypothesis. Furthermore, the corresponding theoretical investigations have demonstrated that the increase of the ligand electronegativity can enhance their discrimination between these two cations effectively.

1. Introduction P-type ATPase is an important group of enzymes that transports ions across the plasma membrane using the energy released from ATP hydrolysis.1 For example, the sarcoplasmic reticulum (SR) Ca2+-ATPase membrane protein, which plays a significant role in controlling the muscle contraction and signaling the metabolism of various cells, is a type II P-type ATPase transporting Ca2+ over other cations across the biological membrane. Experimentally, the Ca2+-ATPase selectively transports Ca2+ in the presence of 103-105-fold higher concentrations of Mg2+.2 However, the selective mechanisms of the Ca2+-ATPase to the calcium ions against the magnesium ions remain ambiguous, and much fewer studies on the basis of the physicochemical characteristics have been reported. Therefore, it is very challenging for us to investigate how Ca2+ could bind to its target protein against such a high background concentration of Mg2+ and what are the factors to manipulate this process. One of the ubiquitous features of membrane proteins is the preference of tryptophan residue for the regions corresponding to membrane surfaces that presumably arises from enhanced stability due to distinct interfacial interactions.3,4 The physical basis for this preference is believed to be dominated by tryptophan’s flat rigid shape that limits access to the hydrocarbon * To whom correspondence should be addressed. E-mail: byx@ sdu.edu.cn. † Shandong University. ‡ Michigan State University.

core. Furthermore, the side chain of tryptophan is an indole, whose representative π electronic structure and associated quadrupolar moment (aromaticity) favor residing in the electrostatically complex interface environment.3 The X-ray structure of Ca2+-ATPase (PDB code: 1SU4) has found some tryptophan residues located at the enzyme-water interface, which also corresponds to the membrane-water interface, and more importantly, Trp50, Trp107, and Trp932 are in the vicinity of the Ca2+ entrance channel proposed by Toyoshima’s5 and Lee’s4 groups, respectively (Figure 1a). As reported previously, the cation-π interaction is a strong noncovalent binding force,6 and the side chains of the amino acids with the electron-rich π ring can stabilize a positive charge as efficiently as solvation by water.7a Therefore, it can be expected that these tryptophan residues may be the effective selector to distinguish Ca2+ from Mg2+ before they enter into the entrance channel. Notably, as shown in Figure 1a, there are other amino acid residues, such as Glu/Asp, Arg/Lys, and Gln/Asn, etc., located on the enzyme surface which could also play important roles in the specificity of the Ca2+-ATPase to transport Ca2+. To the best of our knowledge, investigation on the discrimination of the protein cavities on the enzyme surface between Ca2+ and Mg2+ has hardly been reported previously, and physicochemical information about the ion selectivity is also very scarce. However, the available literature concerning the cation with acidic, neutral, and aromatic side chain of the amino acid residues is abundant from experimental and theoretical aspects, which can supply us the appropriate methods to deal

10.1021/jp073883q CCC: $37.00 © 2007 American Chemical Society Published on Web 10/03/2007

Ca2+ Selectivity of the Sarcoplasmic Reticulum

J. Phys. Chem. B, Vol. 111, No. 42, 2007 12283

Figure 1. (a) The region of the Ca2+-ATPase-H2O interface close to the Ca2+ entry pathway and (b) The two Ca2+ high-affinity binding sites in the 2.6 Å X-ray structure of the SR Ca2+-ATPase (PDB code: 1SU4). The two Ca2+ ions are shown in brown, while the water molecules are represented by the red asterisks. In(a), the route depicted by the arrow, that is, the protein cavity filled with the water molecules, is the possible Ca2+ entrance channel (see refs 4 and 5). The green horizontal line shows the likely position of the membrane-H2O interface, suggested by ref 4. In addition, the region above the horizontal line is the cytoplasmic side, and that below the horizontal line corresponds to the biphospholipid layer.

with the present system and help us understand the interactive characteristics of the cations and ligands.7-15 For example, Dudev et al. addressed that the switch in the carboxylate-binding mode from bidentate to monodentate can affect the selectivity of a protein cavity for a given metal cofactor.8a Previous theoretical studies have also analyzed the interaction of alkali9,10a,11 and alkaline earth metal cations12 and the transition-metal cations copper mono- and dication10b,13 with glycine. On the other hand, Dougherty and co-workers reviewed the important nature of the cation-π interactions in biological structures and artificial receptors.14 Recently, an ab initio calculation on π-cation-π sandwich complexes was performed to investigate the additivity of cation-π interactions.7b Experimentally, the binding energies of a number of metal cations with phenol and indole in the gas phase were explored by radiative association kinetics analysis, supplemented by density functional calculations.15 Obviously, the above investigations exploring the interaction between biologically active metal ions and various amino acids are good guidance for us in the use of computational methods and basal theories. In this work we elucidate the key factors that influence ion selectivity near and in the calcium entrance channel of Ca2+ATPase using the DFT method. Specially, we evaluate if and how the Ca2+ selectivity over Mg2+ in Ca2+-ATPase depends on (i) the nature of the metal ion, (ii) the size and charge of the ligands, (iii) the different coordination forms (the first- vs second-shell coordination) of the ligands to the cations, and (iv) the binding modes (mono- vs bidentate binding) of the negatively charged carboxyl groups. In addition, the valence selectivity of Ca2+ over Na+ and K+ has also been mentioned. 2. Methods 2.1. Models. Although tryptophan, glutamate acid/aspartate acid, and arginine/lysine could exist in more than one ionization state, our focus is on the most likely charge state at physiological pH, so the corresponding models of these amino acid residues are in the neutral, monoanionic, and monocationic forms, respectively. Accordingly, the side chain of Trp is modeled by the indole molecule (C8NH7), the side chains of Asp/Glu and Arg/Lys are modeled by the acetate (CH3COO-) and methylamine cations (CH3NH3+), while that of Asn/Gln is modeled by formamide (HCONH2). In proteins, Mg2+ is predominantly six coordinate,16 as in aqueous solution, while Ca2+ is usually seven coordinate.17 Thus,

complexes containing these two natural metal cations were modeled as MgL6 and CaL7 (L denotes the coordinated molecules), respectively. A general reaction scheme of the binding of a certain metal cation, Ca2+ for example, to the same amino acids in the different forms can be described as follows ∆Ga

A + Ca‚W7 98 A‚Ca‚W7 ∆Ga*

A + Ca‚W7 98 A*‚Ca‚W7*

(1) (2)

where the reactant A is the modeled amino acid residue while the products are the metal-occupied proteins representing the different binding forms. In detail, when A denotes the sole ligand, the products correspond to the two different binding forms, i.e., the ligand binds in the first or second shell of the metal cations. Correspondingly, when two kinds of amino acid residues bind to the hydrated metal cation in a concerted manner, the products correspond to the four-residue binding forms: one treating both of the residues in the first shell, one treating both of them in the second shell, and two each treating only one of them in the first shell and another in the second shell. Obviously, the different products corresponding to the same reactants are the isomeric compounds. In this work, the difference in the binding free energies of these two reaction equations is used to detect the optimal binding geometries for the certain metal ion. On the other hand, a general reaction scheme of the binding of Ca2+ and Mg2+ to the same ligands in the same binding modes can be described as follows ∆Ga

A + Ca‚W7 98 A‚Ca‚W7 ∆Gb

A + Mg‚W6 98 A‚Mg‚W6

(3) (4)

Here, the difference in the binding free energies, which can be used to explore the different affinities of the ligands to Ca2+ and Mg2+, of these two reactions can be achieved by the following equation

12284 J. Phys. Chem. B, Vol. 111, No. 42, 2007

Xiang et al.

∆∆G ) ∆Ga - ∆Gb ) [G(A‚Ca‚W7) - G(A) - G(Ca‚W7)] [G(A‚Mg‚W6) - G(A) - G(Mg‚W6)] ) [G(A‚Ca‚W7) - G(A‚Mg‚W6)] [G(Ca‚W7) - G(Mg‚W6)] (5) Notably, this equation is used to detect the discrimination of the protein cavities between the different metal cations, Ca2+ and Mg2+. 2.2. DFT Calculations To elucidate the geometry and energy characteristics of the complexes, three calculation steps were adopted. First, all complexes were optimized at the HF/3-21G level using the Gaussian 03 code;18 then, the full geometry reoptimizations for these complexes were performed employing the B3LYP19,20 method with the 6-31+G(d,p) basis set. This function was chosen as it reproduces the experimentally observed metal-aliphatic amino acid parameters8b and provides precise results comparable to those generated by MP2 method for complexes with π-electron characteristics.7c Finally, the relevant energies and other quantities were refined using the B3LYP/6-311++G(2d,p) single-point calculations at the B3LYP/ 6-31+G(d,p) geometries. For each fully optimized structure, B3LYP/6-31+G(d,p) vibrational frequency was computed to verify that the molecule was at the minimum of its potential-energy surface (PES). No imaginary frequencies were found in any of the metal complexes. To test the degree of dependence of the results with respect to the functional chosen, the correlative calculations for parts of the complexes were performed using the hybrid BHandHLYP functional21 in combination with the same basis sets as those used in B3LYP calculations. For reactions 1-4 the interaction energy, ∆E, was evaluated as follows

∆E ) ΣE(product) -ΣE(reactant)

(6)

where ΣE(reactant) and ΣE(product) are the total energies of the reactants and products corrected by ZPVE. To obtain the true interaction energies, the basis set superposition error (BSSE) corrections were considered using the counterpoise procedure.22 Additionally, the Gibbs free energy of formation, depicted as ∆G for the above reactions, was determined according to the following formula

∆G ) ΣG°(product) - ΣG°(reactant)

(7)

where the Gibbs free energies were estimated by the total electronic energies and the corresponding thermal corrections. To describe the cation-ligand interactions accurately, the natural bond orbital (NBO) analysis23 used for monitoring the spn hybrid character of the σ-bonding, σ*-antibonding, and lone pair (LP) orbitals was also calculated for the discussed complexes at the B3LYP/6-311++G(2d,p) and BHandHLYP/ 6-311++G(2d,p) levels. Simultaneously, the CHELPG (electrostatic potential charges from electrostatic potentials generalized)24 atomic charges were calculated at the same level. Here, the VDW radii of Ca2+ and Mg2+ required in the computation were cited from the CRC Handbook of Chemistry and Physics.25 All of the calculations were performed using the Gaussian 03 suite of programs.18 3. Results 3.1. Discrimination of the Aromatic Indole Ligand Between Ca2+ and Mg2+. Before the metal ions reach a protein

cavity they could exist in the hydrated forms. As mentioned above, the tryptophan residue, which is the ubiquitous molecule at the interface of Ca2+-ATPase-water, plays very important roles structurally and functionally. Then, can indole, the side chain of Trp, be the effective selector for Ca2+ and Mg2+? If it is the case, what coordination mode would indole take; i.e., could indole bind to the cations by displacing a metal-bound water molecule or just bind in the second-coordination shell of the metal cations? To answer these questions we compared the two different coordination modes: (i) the indole molecule binds to the metal cations in the first coordination shell, namely

Ind + (M‚Wn)2+ ) [(Ind‚M‚Wn-1)‚W]2+

(1a)

and (ii) the indole molecule exists in the second coordination shell of the metal ions

Ind + (M‚Wn)2+ ) [Ind‚(M‚Wn)]2+

(1b)

where M denotes the metal ions and n ) 7 when M is Ca while n ) 6 for Mg; Ind denotes the indole molecule, and W denotes the water molecule. What should also be emphasized is the use of brackets: the molecules in parentheses are the ones binding in the first coordination shell of the metal ions, and the molecules in the square brackets but out of the parentheses are those binding in the second shell. The same may apply to other reaction equations mentioned below unless indicated elsewhere. To simplify the representation, the coordinations of the sole indole molecule to Ca2+ and Mg2+ are denoted as the Trp case. The structures of the fully optimized [(Ind‚Ca‚W6)‚W]2+ and [Ind‚(Ca‚W7)]2+ complexes are shown in Figure 2a, and the corresponding [(Ind‚Mg‚W5)‚W]2+ and [Ind‚(Mg‚W6)]2+ complexes are depicted in Figure S2 of the Supporting Information. The differences in the binding free energies (∆∆G) of the two metal cations corresponding to the same binding mode are listed in Table 1. Simultaneously, the detailed enthalpies, entropies, free energies, and binding energies (after the corrections of BSSE and ZPVE) are also listed. As shown for the electrostatic potential surface in Figure 3, the indole molecule has three electron-rich zones (color-coded in red), namely, the delocalized π-electrons ring of benzene and pyrrole and the electronegative N atom, representing the possible binding sites of the indole molecule to hydrated Ca2+ and Mg2+. In contrast, the most electropositive regions (color-coded in blue) occur near the plane of the molecule around the hydrogen atoms. Notably, the electron densities above the phenyl and pyrrolyl rings are much larger than that above the N atom. On the basis of the electrostatic potential map, when indole binds in the first shell of the metal cation, three binding modes are designed, that is, the metal-N, metal-pyrrole ring, and metal-benzene ring forms. The corresponding reaction numbers are represented as 1a-1, 1a-2, and 1a-3, respectively. Subsequently, plenty of calculations have been performed to find the minima on the PES corresponding to the different binding modes. However, not all the desired geometries have been optimized successfully after many endeavors. The metal-benzene ring binding modes are not found for both of the metal cations, i.e., only reactions 1a-1 and 1a-2 can take place. In fact, when binding to the first shell, all the electron-abundant zones of indole participate in the reactions, some of which bind to the positively charged metal cations, and others stabilize the complexes with the O-H‚‚‚N lone pair or O-H‚‚‚π interactions. Obviously, when indole participates in the reaction as the sole protein ligand it prefers to bind Mg2+ rather than Ca2+ either in the first-shell mode or the second-shell one (∆∆G )

Ca2+ Selectivity of the Sarcoplasmic Reticulum

J. Phys. Chem. B, Vol. 111, No. 42, 2007 12285

Figure 2. Parts of B3LYP-optimized geometries for the Ca2+-products of (a) Trp case, (b) Trp-Gln case, (c) Trp-Glu case, and (d) Trp-Trp case. Several unmentioned geometries in the text are shown in Figure S1 of the Supporting Information. Ca2+ is denoted as brown, O atom is represented as red, N atom is denoted as blue, C atom is denoted as gray, and H atom is denoted as white-gray, respectively. Distances of the H-bonds and metal-O bonds are in Å. In the nomenclatures of the products, the former two terms denote the corresponding reaction numbers, and the last one denotes the isomeric products for the same reaction number. In addition, the italic letter m represents the monodentate coordination of acetate in the first shell, while b denotes the bidentate coordination of acetate.

12286 J. Phys. Chem. B, Vol. 111, No. 42, 2007

Xiang et al.

Figure 3. Electrostatic potential surfaces of indole calculated by DFT/ B3LYP method (left) and the corresponding geometry of the optimized indole (right). The color code spans from -1.000e-2 (red) to 1.000e-2 (blue).

TABLE 1: B3LYP Calculated Free Energies (∆G), Enthalpies (∆H), Entropies Terms (T∆S), Binding Energies (∆E), and Difference in the Binding Free Energies (∆∆G) for the Trp Case in the Gas Phase (in kcal/mol) reaction no.

coord. metal

1a-1b

Ca2+

1a-2c 1b

Mg2+ Ca2+ Mg2+ Ca2+ Mg2+

∆G

∆H

T∆S

∆E

∆∆Ga

-1.78 -5.58 -4.89 -7.67 -8.82 -13.61

-14.79 -17.44 -18.14 -18.58 -20.32 -25.38

-13.00 -11.85 -13.24 -10.91 -11.50 -11.76

-14.84 -17.40 -18.49 -18.66 -20.35 -25.00

-3.80 -2.78 -4.79

a ∆∆G ) ∆G (Mg2+ complexes) - ∆G (Ca2+ complexes). b The hydrated metal cations bind the indole ligand in the cation-N and O-H‚‚‚π forms. c The hydrated metal cations bind the indole ligand in the cation-π (pyrrole ring) and O-H‚‚‚π forms. Detailed geometrical information is shown in Figures 2 and S2.

-3.80/-2.78 and -4.79 kcal/mol for reactions 1a and 1b, respectively). This may arise from the differences of the electrostatic interaction of these complexes because in the gas phase this factor plays a prominent role in cation-π, cationN, O-H‚‚‚N lone pair, O-H‚‚‚O lone pair, and O-H‚‚‚π binding modes existing in the complexes studied in this work, although the charge-transfer term may also be important.14 Among these five interactions the former two terms appear in the first-shell binding mode, while the other three exist in all the binding modes. Relative to Mg2+, Ca2+ is much larger (the ionic radii of Ca2+ and Mg2+ are 1.00 and 0.72, respectively26), having a lower charge density and thus weaker electron acceptability. Consequently, for reaction 1a, the Ca2+-dipole (for cation-N) and Ca2+-quadrupole (for cation-π) interactions are smaller than those in Mg2+ complexes. At the same time, the stronger Mg2+-O couplings have weakened the O-H bonds consisting of the O-H‚‚‚π interactions, resulting in the more favorable formation of these π bonds, which are mainly the dipole-quadrupole and quadrupole-quadrupole interactions. For example, for 1aCa1 and 1aMg1 both have two obvious O-H‚ ‚‚π interactions, and the corresponding H-bond lengths in 1aCa1 are 2.64 and 2.40 Å, while those in 1aMg1 are only 2.51 and 2.24 Å, indicating the stronger binding of the O-H‚‚‚π bonds in Mg2+ complex. Correspondingly, the O-H stretching vibrations in 1aMg1 relative to those in 1aCa1 have blue shifted 58 and 20 cm-1 (Table 2). In addition, the variation of the O-H‚ ‚‚π interactions can be reflected by the spn hybridization of the O-H bonds;27 that is, the bond is weakened and elongated along with the increase of the p character of the σ orbital. The (OH) σ orbital of these complexes is the linear combination of the two orbitals localized at O and H. Moreover, the hybrid character of the (O-H) σ orbital is determined by the change of the orbital localized at O because the orbital localized at the hydrogen atom is of 100% s character in all complexes. As

shown in Table 2, the calculated spn hybrid characters 2.75/ 2.97, which correspond to the two O-H‚‚‚π bonds, of the orbital localized at O in 1aCa1 decrease to 2.57/2.69 in 1aMg1. Similarly, for reaction 1b, the preferred coordination of indole to Mg2+ can also be understood from the above explanations. To elucidate the ligand discrimination more elaborately, two different types of charge distributions, namely, NBO23 and CHelpG24 charges, have also been determined. As expected, all the calculated charge distributions indicate that the electron transfers take place from the indole rings to the hydrated metal cations. This transfer is defined as ∆q and listed in Figure 4a. The detailed charge distributions are shown in Table S1 of the Supporting Information. Concretely to say, both the NBO and CHelpG charge analyses show more electrons transferring from the indole ring to Mg2+ than to Ca2+ (Figure 4a), in accordance with the stronger electron acceptability of Mg2+. For example, the NBO charge analyses (listed in Table S1) have indicated that in complex 1aMg1 0.06 electrons transfer from the indole molecule to Mg2+ while in 1aCa1 the transfer amount is 0.04. The CHelpG charge analyses have exhibited much more electrons transferring from indole to Mg2+ than to Ca2+ (0.35 vs 0.27). These facts demonstrate a stronger interaction between Mg2+ and the electron-rich indole ring. Notably, when the indole molecule binds on the second coordination shell, the ligand discrimination between Ca2+ and Mg2+ is better than the first-shell mode (∆∆G ) -4.79 vs -3.80 and -2.78 kcal/mol). This may be understood from the big steric effect of indole. Simultaneously, the second-shell mode supplies the bigger binding strength and Gibbs free energy. For example, ∆E and ∆G of 1bCa1 are bigger than those of 1aCa1 and 1aCa2 (see Table 1). In addition, although exhibiting the different bonding strength to the ligands, the Mg2+ and Ca2+ complexes have several similar geometrical and thermodynamic characteristics for reactions 1a and 1b: (i) comparison of ∆G and ∆E shows that the indole molecule prefers to bind in the second coordination shell of the metal ions, indicating that the gas-phase reaction is driven rather by favorable dipole-quadrupole and quadrupolequadrupole interactions (for (O)-H-π) than the chargequadrupole interaction (for cation-π). On the other hand, desolvation of the metal ions is difficult in the gas phase; (ii) although all reactions have the unfavorable entropic characteristics due to the ordered arrangement of the water molecules around the metal cations, the negative and big enthalpy terms determine whether the reactions are enthalpy driven; (iii) when the indole ligand binds in the first coordination shell, the cation-π binding mode is more favorable than the cation-N form, indicating that the electron density above the π ring is much larger than that above the N atom, consistent with the electrostatic surface shown in Figure 3. To test the degree of dependence of the results with respect to the functional chosen, all complexes in the Trp case have been recalculated using the BHandHLYP functional and the same basis sets as those used in B3LYP calculations. Detailed data are shown in Figures S3 and S4 and Tables S2 and S3 of the Supporting Information. Obviously, the coordination modes and main bond lengths calculated at B3LYP and BHandHLYP functionals are similar. More importantly, the energetic calculations, NBO analyses, and charge transfers determined at the B3LYP level have exhibited the same results as those from BHandHLYP calculations, that is, the indole molecule has stronger affinity to Mg2+ than to Ca2+ in all binding modes. Considering the consistency of the results calculated by these

Ca2+ Selectivity of the Sarcoplasmic Reticulum

J. Phys. Chem. B, Vol. 111, No. 42, 2007 12287

Figure 4. The charge transfer from the protein ligands to the hydrated metal cations (∆q) of all products calculated at the B3LYP/6-311++g(2d,p) level. The horizon axis denotes the products corresponding to the different reactions. (a) The charge transfer from the single indole ligand to the hydrated metal cations; (b) The charge transfer from the indole and formamide ligands to the hydrated metal cations; (c) The charge transfer from the indole and acetate ligands to the hydrated metal cations; and (d) The charge transfer from the two indole ligands to the hydrated metal cations.

two functionals, the B3LYP results are reliable and used to analyze the systems below unless indicated elsewhere. In conclusion, the aromatic indole molecule prefers to bind Mg2+ rather than Ca2+ in the second-shell mode in the gas phase. Additionally, when the indole ligand binds in the first coordination shell, the metal-π mode is more favorable than the metal-N one. 3.2. Concerted Discrimination of the Aromatic Indole and Other Aliphatic Amino Acid Ligands between Hydrated Ca2+ and Mg2+. The investigations above have shown that the indole molecule prefers to bind Mg2+ in the second coordination shell. However, on the surface of the Ca2+-ATPase, other

residues (Glu/Asp, Gln/Asn, Arg/Lys, etc.) may also have important roles on the ligand discrimination between Ca2+ and Mg2+. For example, as shown in Figure 1a, Asn111 interacts with Trp107 directly and Arg110 exists in close proximity to Trp107. Simultaneously, there are several acidic amino acid residues (Glu51, Glu55, Glu109, etc.) located near the possible Ca2+ entrance channel that help to define the surface of the Ca2+-ATPase. In addition, basic residues exist, such as Arg63, to stabilize the region composed of the acidic residues mainly. Thus, these residues and Trp may work in concert to select the appropriate metal cation. According to these facts we modeled the following reaction processes

TABLE 2: Stretching Vibrational Frequencies (in cm-1) and the Corresponding NBO Analyses for the O-H‚‚‚π Bonds in Trp Casea 1aCa1

1aMg1

1aCa2

1aMg2

Str(O-H‚‚‚π)

3688 3823

3630 3803

3649 3856

3566 3833

NBO[O-H(O)]

73.25%/2.75 74.74%/2.97

71.96%/2.57 72.88%/2.69

75.05%/3.02 74.68%/2.96

71.25%/2.48 73.21%/2.74

1bCa1

1bMg1

3541 3557 3683 71.95%/2.57 72.21%/2.61 73.36%/2.76

3444 3532 3584 70.57%/2.40 71.22%/2.48 71.67%/2.54

a In NBO analyses the percentages before the slash refer to the calculated weight of the orbital localized at the atom in the bracket and the data behind the slash refer to the spn hybridization of the orbital localized at the atom in the bracket. Otherwise, calculations of the vibrational frequencies are at the B3LYP/6-31+G(d,p) level and those of NBO data are at the B3LYP/6-311++G(2d,p) level.

12288 J. Phys. Chem. B, Vol. 111, No. 42, 2007

Xiang et al.

Ind + (M‚Wn)2+ + A ) [(Ind‚M‚A‚Wn-2)‚W2]2+

(2a)

Ind + (M‚Wn)2+ + A ) [Ind‚(A‚M‚Wn-1)‚W]2+

(2b)

TABLE 3: B3LYP Calculated Free Energies (∆G), Enthalpies (∆H), Entropies Terms (T∆S), Binding Energies (∆E), and Difference in the Binding Free Energies (∆∆G) for the Trp-Gln Case in the Gas Phase (in kcal/mol)

Ind + (M‚Wn)2+ + A ) [A‚(Ind‚M‚Wn-1)‚W]2+

(2c)

reaction no.

coord. metal

2a-1b

Ca2+

Ind + (M‚Wn)2+ + A ) [Ind‚(M‚Wn)‚A]2+ HCONH2, CH3COO-,

(2d) CH3NH3+

where A denotes the or group, which is representative of the side chain of Gln/Asn, Glu/Asp, or Arg/Lys, respectively. As a versatile π-donor, the indole side chain of tryptophan can also form a π-cation-π model, which is ubiquitous in vivo.28,29 At the interface of Ca2+-ATPase-water it can be noted that the two indole side chains, belonging to respective Trp107 and Trp932, are in close proximity to each other. In fact, analysis of the crystal structure of the Ca2+-ATPase in the Ca2+-bound state has revealed that the distance between the nearest carbon atoms of these two side chains is only 4.80 Å. Notably, there is another aromatic side chain of Phe809 whose centroid of the benzene ring is only 5.09 Å away from that of Trp932. Thus, investigations of the π-cation-π interactions are also important for the ion selectivity at the enzyme-H2O interface. Similar to the designs above, the following reaction equations are considered

Ind + (M‚Wn)2+ + Ind ) [(Ind‚M‚Wn-2‚Ind)‚2W]2+

(2e)

Ind + (M‚Wn)2+ + Ind ) [Ind‚(Ind‚M‚Wn-1)‚W]2+

(2f)

Ind + (M‚Wn)2+ + Ind ) [Ind‚(M‚Wn)‚Ind]2+

(2g)

To simplify the representations, the different ligand binding forms, that is, the concerted indole-formamide, indole-acetate, indole-methylamine cation, and indole-indole coordinations, are denoted as Trp-Gln, Trp-Glu, Trp-Arg, and Trp-Trp cases, respectively. 3.2.1. Trp-Gln Case. Compared to the Trp case, the TrpGln case exhibits a stronger affinity to Mg2+ than to Ca2+ (-5 to -14 kcal/mol for ∆∆G of reactions 2a-2d vs -3 to -5 kcal/mol for those of reactions 1a and 1b shown in Tables 1 and 3), that is, participation of the formamide molecule to the coordination enhances ligand discrimination between Ca2+ and Mg2+. To verify this point we calculated the sole binding of formamide to these two cations. As shown in Tables 1 and 4, ∆∆G of reaction 1c is more negative than that of reaction 1a (-5.36 vs -3.80 and -2.78 kcal/mol). This arises from, on one hand, the bigger electronegativity of formamide relative to indole as well as, on the other hand, the stronger electron acceptability of Mg2+ than Ca2+, which results in the bigger sensitivity of Mg2+ to the increase of the electron donation of the ligand. For example, in Tables 1 and 4 the Mg2+-formamide interaction is 20.97 and 18.88 kcal/mol more approved than the two kinds of Mg2+-indole interactions, while the corresponding Ca2+-formamide interaction is 19.41 and 16.30 kcal/mol more approved than the Ca2+-indole ones. As expected, these factors have resulted in the considerable geometrical changes of the complexes. For products with a similar binding mode, for example, 1aCa1/1aMg1, 2aCa1/2aMg1, and 2cCa1/2cMg1, which exhibit the Ca2+-N and O-H‚‚‚π coordination forms, introduction of the formamide ligand has lengthened the distances of the metal-π systems to different degrees. As shown in Figures 2 and S2, for the Ca2+ complexes the Ca2+-centroid bond lengths of 1aCa1, 2aCa1, and 2cCa1 are 2.81, 2.83, and 2.86 Å,

2a-2c 2b-1d 2b-2e 2c-1b 2c-2f 2c-3c 2d

Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+

∆G

∆H

T∆S

∆E

∆∆Ga

-17.98 -22.61 -23.78 -30.77 -21.17 -34.78 -26.84 -34.56 -14.40 -25.91 -19.77 -26.05 -19.37 -27.44 -24.40 -32.92

-42.75 -46.25 -47.33 -52.80 -46.27 -57.02 -49.63 -56.58 -40.72 -48.37 -43.59 -47.95 -43.29 -49.24 -46.48 -55.35

-24.76 -23.63 -23.53 -22.02 -25.09 -22.23 -22.78 -22.01 -26.31 -22.45 -23.82 -21.90 -23.90 -21.79 -22.06 -22.42

-42.53 -44.94 -47.29 -52.16 -46.02 -57.01 -49.95 -56.51 -40.05 -47.49 -43.37 -47.34 -43.38 -48.54 -45.43 -53.57

-4.63 -6.99 -13.61 -7.72 -11.51 -6.28 -8.07 -8.52

a ∆∆G ) ∆G (Mg2+ complexes) - ∆G (Ca2+ complexes). b The hydrated metal cations bind the indole ligand in the cation-N and O-H‚‚‚π forms. c The hydrated metal cations bind the indole ligand in the cation-π (benzene ring) and O-H‚‚‚π forms. d The formamide molecule is located above the pyrrole ring. e The formamide molecule is located above the benzene ring. f The hydrated metal cations bind the indole ligand in the cation-π (pyrrole ring) and O-H‚‚‚π forms. Setailed geometrical information is shown in Figures 2 and S2.

TABLE 4: B3LYP Calculated Free Energies (∆G), Enthalpies (∆H), Entropies Terms (T∆S), Binding Energies (∆E), and Difference in the Binding Free Energies (∆∆G) for Reactions of the Respective Formamide and Acetate Ligands Binding to Hydrated Ca2+ and Mg2+ in the Gas Phase (in kcal/mol) reaction coord. no.a metal 1c 1d 1e

Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+

∆G

∆H

T∆S

-21.19 -26.55 -187.06 -206.87 -193.05 -207.62

-31.86 -36.63 -203.47 -214.67 -209.68 -219.70

-10.67 -10.08 -16.40 -7.79 -16.62 -12.07

∆E

∆∆Gb

-31.50 -5.36 -36.29 -203.68 -12.25 -215.93 -210.06 -9.32 -219.38

a Reaction numbers 1c, 1d, and 1e represent the following reactions, respectively: HCONH2 + (M‚Wn)2+ ) [(HCO ˜ NH2‚M‚Wn-1)‚W]2+, CH3COO- + (M‚Wn)2+ ) [(CH3COO ˜ -‚M‚Wn-1)‚W]2+, and CH3COO+ (M‚Wn)2+ ) [(CH3CO ˜O ˜ -‚M‚Wn-2)‚2W]2+ where O ˜ denotes the bounded atom and M is the coordination metal cation. b ∆∆G ) ∆G (Mg2+ complexes) - ∆G (Ca2+ complexes).

respectively. Thus, the attraction of the formamide O atom to Ca2+ has lengthened the Ca2+-centroid bonds only 0.02 and 0.05 Å. However, for the Mg2+ complexes the corresponding bond lengths of 1aMg1, 2aMg1, and 2cMg1 are 2.20, 2.40, and 2.37 Å, i.e., the attraction of the formamide O to the positively charged magnesium has lengthened the Mg-centroid bonds 0.20 and 0.17 Å, respectively, in accordance with the strong Mg2+formamide interaction. Additionally, for Ca2+ and Mg2+ complexes with a similar binding mode analyses of the charge distributions have also revealed considerable electron transfers from the ligands to Mg2+ than to Ca2+ (see Figure 4b) due to the larger ability of Mg2+ to accept electrons. It is interesting to find that the different binding mode of the ligands (the first- or second-shell mode) has supplied the different reaction preference of the protein ligands to the metal cations. In detail, the optimal one corresponds to reaction 2b1, whose product exhibits the geometrical character of the firstshell formamide ligand located above the second-shell pyrrole

Ca2+ Selectivity of the Sarcoplasmic Reticulum ring (Figures 2 and S2). This reactive preference could also be explained by the electrostatic interaction. With the big electron donation, the formamide prefers to bind to the metal cations directly in the first-shell mode. Simultaneously, the big steric effect of indole hampers its binding in the primary coordination sphere. Although exhibiting a different reactive preference, the Ca2+ and Mg2+ complexes have several similarities because Ca and Mg are in the same group of the Periodic Table: (i) all reaction processes are thermodynamically favorable and participation of the formamide molecule increases the metal-ligand affinities because of the enhanced electrostatic interactions; (ii) when indole binds in the first coordination shell, the energetic calculations shown in Table 3 also exhibited a bigger electronegativity of the electron-rich π ring than the N atom with ∆G [metal-π] > ∆G [metal-N]. Introduction of the formamide ligand plays dual roles: on one hand, the big electronegativity of its oxygen atom makes it a good discriminator to Ca2+ and Mg2+; and on the other hand, its different binding modes (the first- vs second-shell mode) have supplied the optimal coordination forms of the ligands to the metal cations. 3.2.2. Trp-Glu Case. As the most common amino acid ligands coordinated to metal cations in metalloproteins, the carboxylatecontaining residues, viz., deprotonated aspartates and glutamates, are able to bind the metal cation either monodentately (via one of the carboxylate oxygens) or bidentately (with both carboxylate oxygens). Furthermore, these two different coordination modes could affect the geometrical and energetic preferences of the metal-binding sites. 8a Thus, the different binding modes (mono- and bidentate forms) of the acidic side chains of Asp/ Glu are both considered for eqs 2a′-2d′. To differentiate from reactions 2a-2d, the corresponding reactions describing the Trp-Glu case are represented as 2a′-2d′. Simultaneously, when acetate, modeling the side chain of Asp/Glu, binds in the first shell, its mono- and bidentate forms are represented as Am′ and Ab′, respectively. As pointed out, when acetate coordinates in the second shell, only bidentate forms are listed because compared to the monodentate binding their increased hydrogen bonds can stabilize the products better.8a Detailed geometrical and thermodynamic information is shown in Figures 2, S1, and S2 and Table 5. In analogy to the Trp case, the Trp-Glu case has enhanced the ligand discrimination between Ca2+ and Mg2+ with ∆∆G of reactions 2a′-2d′ being in the range from -13 to -19 kcal/ mol (Table 5). Obviously, this enhancement comes from the strong electrostatic interaction between negatively charged acetate ligand and positively charged metal cations (see ∆∆G of reactions 1d and 1e shown in Table 4). In fact, the negatively charged acetate is a such a strong electron donor that the metalcarboxylate O interaction can be regarded as an ionic bond, while the metal-formamide O, meta--π, and O-H‚‚‚π coupling are weakly electrostatic interactions. Simultaneously, analyses of ∆q have shown that, relative to Ca2+, most of the reactions exhibit the increased charge transfer from the ligands to Mg2+, consistent with the stronger electrostatic interaction between Mg2+ and the ligands (Figure 4c). Analyses of the thermodynamic data reveal that the different carboxylate-binding modes (mono- vs bidentate forms; first- vs second-shell forms) have obvious influences on the free-energy gain. As shown in Table 5, of all the reaction, 2bm′-1 has the most negative ∆∆G value (-18.81 kcal/mol), indicating that this reaction is the optimal one to discriminate between Ca2+ and Mg2+. Corresponding to this reaction, its products 2b′Ca1-m

J. Phys. Chem. B, Vol. 111, No. 42, 2007 12289 TABLE 5: B3LYP Calculated Free Energies (∆G), Enthalpies (∆H), Entropies Terms (T∆S), Binding Energies (∆E), and Difference in the Binding Free Energies (∆∆G) for the Trp-Glu Case in the Gas Phase (in kcal/mol) reaction coord. no.a metal 2am-1c 2am-2d 2am-3e 2ab 2bm-1f 2bm-2g 2bb 2c-1c 2c-2d 2d

Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+

∆E

∆∆Gb

∆G

∆H

T∆S

-183.07 -196.41 -183.03 -196.50 -182.45 -198.78

-211.43 -221.61 -211.17 -220.20 -210.42 -221.89

-28.34 -25.19 -28.13 -23.69 -27.96 -23.10

-211.49 -13.34 -221.05 -212.49 -13.47 -220.11 -211.61 -16.33 -222.33

-185.19 -204.00 -187.19 -202.45 -189.66 -202.58 -173.22 -189.54 -174.96 -192.39 -180.50 -198.32

-211.93 -229.82 -212.56 -228.03 -217.80 -227.12 -205.59 -217.54 -202.17 -217.44 -207.62 -224.14

-26.73 -25.81 -25.36 -25.57 -28.13 -24.53 -32.35 -28.00 -27.20 -25.03 -27.11 -25.81

-214.00 -229.49 -214.46 -227.60 -218.73 -227.03 -205.67 -215.47 -202.92 -216.70 -207.78 -221.50

-18.81 -15.26 -12.92 -16.32 -17.43 -17.82

a The subscript m denotes the monodentate binding mode of the acetate ligand, while b denotes the bidentate one. b ∆∆G ) ∆G (Mg2+ complexes) - ∆G (Ca2+ complexes). c The hydrated metal cations bind the indole ligand in the cation-N and O-H‚‚‚π forms. d The hydrated metal cations bind the indole ligand in the cation-π (pyrrole ring) and O-H‚‚‚π forms. e The hydrated metal cations bind the indole ligand in the cation-π (benzene ring) and O-H‚‚‚π forms. f The acetate molecule is located above the pyrrole ring, and the unbound O atom is cis to the pyrrolyl N. g The acetate molecule is located above the benzene ring, and the unbound O atom is trans to the pyrrolyl N. Detailed geometrical information is shown in Figures 2, S1, and S2.

and 2b′Mg1-m have monodentately bound acetate, which is located above the pyrrole ring in the first shell and indole in the second shell. Simultaneously, the unbound acetate O atom is cis to the pyrrolyl N (depicted in Figures 2 and S2). As for the bidentately coordinated acetate, the corresponding ∆∆G of reaction 2bb′ is -12.92 kcal/mol, much smaller than the monodentately bound form due to the following reasons. First, when binding in the first shell, relative to the monodentate form, the bidentate binding of acetate may have a bigger steric effect, which destabilizes the metal-ligand binding. Second, in the monodentate form the uncoupled carboxylate oxygen has formed two O-H‚‚‚O bonds, which can stabilize the metal-ligand interaction further (Figures 2 and S2). Calculations on the sole acetate ligating in the mono- and bidentate forms have also shown similar thermodynamic characteristics to those shown in reactions 1d and 1e of Table 4. In addition, the BHandHLYP calculations have shown results similar to those of the B3LYP calculations (Figures S3 and S4 and Table S4 of the Supporting Information). The other reactive characteristics of the Trp-Glu case are similar to those of the Trp-Gln one and will not be discussed in detail. 3.2.3. Trp-Arg Case. An attempt to optimize the concerted binding of methylamine cation, the representation of the side chains of Arg and Lys residues, and hydrated metal ions to indole resulted in the rapid detaching of CH3NH3+ from the complexes composed of indole and hydrated metal cations, indicating that a single π system like indole cannot afford sufficient electronegativity to bind these two kinds of cations. On the other hand, this fact has suggested that the indole-

12290 J. Phys. Chem. B, Vol. 111, No. 42, 2007

Xiang et al.

TABLE 6: B3LYP Calculated Free Energies (∆G), Enthalpies (∆H), Entropies Terms (T∆S), Binding Energies (∆E), and Difference in the Binding Free Energies (∆∆G) for the Trp-Trp Case in the Gas Phase (in kcal/mol) reaction no. 2e 2f-1b 2f-2c 2g-1d 2g-2e 2g-3f 2g-4g

coord. metal

∆G

∆H

T∆S

∆E

∆∆Ga

-6.06 -11.31 -9.08 -15.98 -12.14 -21.97 -10.88 -20.45 -12.38 -22.22 -12.26 -20.33

-32.34 -35.55 -33.97 -38.56 -36.39 -46.09 -35.73 -44.15 -37.16 -46.05 -36.34 -44.34

-26.27 -24.23 -24.88 -22.56 -24.24 -24.11 -24.84 -23.69 -24.77 -23.82 -24.07 -24.01

-32.82 -35.30 -34.64 -38.89 -37.48 -46.57 -36.63 -44.45 -38.24 -46.59 -37.50 -44.54

-5.25

Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+ Ca2+ Mg2+

-6.90 -9.83 -9.57 -9.84 -8.07

a ∆∆G ) ∆G (Mg2+ complexes) - ∆G (Ca2+ complexes). b The hydrated metal cations bind the first-shell indole ligand in the cation-N and O-H‚‚‚π forms. c The hydrated metal cations bind the first-shell indole ligand in the cation-π (benzene ring) and O-H‚‚‚π forms.d-g The 2 second-shell indole ligands bind the hydrated metal cations in the four different forms, and detailed geometrical information is shown in Figures 2, S1, and S2.

methylamine cation interaction is weaker than that between indole and the hydrated metal cations. Similarly, in the calculations of Ma and Dougherty they found that the binding energy of CH3NH3+-indole is -25.7 kcal/mol, smaller than that of Na+-indole (-32.6 kcal/mol).14 3.2.4. Trp-Trp Case. Due to the importance of the π-cation-π interaction widely existing in the metalloprotein, reactions 2e-2g are calculated to explore whether concerted coordination of two aromatic indole molecules to Ca2+ and Mg2+ can discriminate these two cations effectively. The products of the corresponding reactions are shown in Figures 2, S1, and S2, and the thermodynamic data are displayed in Table 6. Notably, due to the big steric effect of the indole ligand, optimizations have not found the products of reaction 2e, in which both indoles bind in the first shell. Compared with the Trp case whose ∆∆G varies from -3.00 to -5.00 kcal/mol, the Trp-Trp case has shown better discrimination between Ca2+ and Mg2+ with ∆∆G varying between -5.00 and -10.00 kcal/mol, i.e., the ligands prefer to bind the hydrated Mg2+ rather than hydrated Ca2+. This is understandable because, relative to the sole indole ligating, the two indole molecules exhibit the bigger affinity to Mg2+. Notably, in this case the reaction with the biggest discrimination between Ca2+ and Mg2+ corresponds to 2g-3 in which the two indole molecules bind in the second shell. Again, analyses of the charge transfer reveal results consistent with the thermodynamics data (Figure 4d). Obviously, the bigger electrostatic interaction and charge transfer between indole and Mg2+ should be the main resource of the bigger affinity of the ligands to this cation. 4. Discussion The energetic analyses show that despite the size, charge, and binding modes, the coordinated protein ligands located on the surface of the Ca2+-ATPase prefer to bind smaller size Mg2+ rather than Ca2+. Furthermore, the nature of both the metal and the protein ligands can affect ion selectivity.

4.1. Dependence on the Nature of the Metal Ions. Being in the same group of the Periodic Table, Ca2+ and Mg2+ both have strong acidities to bind the negatively charged ligands. However, relative to Ca2+, Mg2+ has a smaller ionic radius,26 having higher charge density and thus stronger electron acceptability. As a result, Mg2+ with the smaller size prefers to exhibit a coordination number of six,16 while the bigger effective ionic radius of Ca2+ determines its coordination number to be seven in Ca2+ proteins.17 Consequently, the electrostatic interactions between Mg2+ and the electron-rich ligands, including Mg2+dipole (for cation-O and cation-N forms) and Mg2+quadrupole (for cation-π) interactions, are bigger than those in the Ca2+ complexes. As reported previously,14 the electrostatic interaction is the primary factor to determine the metal-ligand interaction. Thus, for the same coordination mode, the ligand has shown a bigger affinity to Mg2+ than Ca2+. On the other hand, Mg2+ is a better electron acceptor than Ca2+ and has a more sensitive response to the changes of the ligand electronegativity. Thus, when the electronegativity, i.e., the electron donation, of the protein ligands increases, the ligand discrimination between Ca2+ and Mg2+ (∆∆G) is also enhanced correspondingly. 4.2. Dependence on the Coordinated Protein Ligands. The thermodynamics calculations in this work have shown that the physicochemical characteristics of the different ligands have essential influences on the ion selectivity. 4.2.1. Effect of Ligand ElectronegatiVity. As revealed by Tables 1-6, the electronegativity of the ligands can affect the ligand discrimination between Ca2+ and Mg2+. For a certain cation the larger the electronegativity the ligand has, the stronger the electrostatic interaction the ligand-cation complex will exhibit. For the investigated ligands the order of the electronegativity is as follows: CH3COO- > HCONH2 > C8NH7 (Tables 1 and 4), which determines the order of the reaction trends (∆G) and binding strengths (∆E) of the corresponding reactions: Trp-Glu case > Trp-Gln case > Trp-Trp case > Trp case. In fact, for the negatively charged acetate its electronegativity is so big that in complex 2c′Mg2, a proton transfers from water to acetate O atom (Figure S2). It is essential to point out that the enhanced electronegativity results in the bigger affinity of the ligands to Mg2+ rather than to Ca2+ because of the metallic nature, as shown in Table 4. Consistent with this, calculations in this work have revealed that the order of the discrimination of the studied cases between Ca2+ and Mg2+ is Trp-Glu > Trp-Gln > Trp-Trp > Trp. 4.2.2. Effect of Ligand-Binding Forms. The calculations in the present work have indicated that the binding forms of the ligands, that is, the first- and second-shell binding, also have influences on the protein discrimination between Ca2+ and Mg2+. For example, for Trp and Trp-Trp cases, whose ligand is an indole molecule, the coordination reactions with the biggest ion selectivity correspond to those in which the protein reactants bind to the hydrated metal cations in the second shell. In detail, for the Trp case, the ∆∆G values of the first- and second-shell binding forms are, respectively, -3.80 (-2.78) and -4.79 kcal/ mol. Similarly, in the Trp-Trp case the reaction with both protein ligands binding in the second shell is 2.94 kcal/mol more favored than that with the first- and second-shell protein binding form. In contrast, for the Trp-Gln and Trp-Glu cases, in which the ligands are the aromatic indole and aliphatic amino acid residues, the coordination reactions with the biggest ion selectivity correspond to those whose products have the monodentate aliphatic residue in the first shell and aromatic indole in the second shell. Clearly, the strong electrostatic interactions

Ca2+ Selectivity of the Sarcoplasmic Reticulum

J. Phys. Chem. B, Vol. 111, No. 42, 2007 12291

Figure 5. The sketch map of the ion selectivity in the Ca2+-channel. (a) The stronger affinity of the negatively charged residues located at the enzyme-H2O interface to Mg2+ facilitates the Ca2+ permeation through the Ca2+-ATPase. (b) The stronger affinity of the negatively charged residues located in the entrance channel to Mg2+ facilitates the Ca2+ permeation through the Ca2+-ATPase. (c) When one Ca2+ ion is electrostatically bound to the negatively charged residues located in the selectivity channel, it can only be removed with the aid of Coulomb repulsion from another Ca2+ ion.37(d) (d) The lesser repulsion from Na+ is unable to displace a resident Ca2+ ion, thus Ca2+ blocks Na+ currents.37(d)

between the electron-rich O atoms of the aliphatic ligands and the metal cations have made these ligands prefer to bind the metal cations directly, while for the π-electron-abundant indole, its big steric effect hampers its binding in the first shell. In addition, for the negatively charged acetate ligand different binding modes (mono- or bidentate binding) exist, which is also an important factor to determine ion selectivity. When the monodentate acetate switches to the bidentate one (this phenomenon is known as a carboxylate shift30 and is very important in the metalloproteinases), the ligand discrimination between Ca2+ and Mg2+ will decrease with the ∆∆G values of reaction 2bm′ changing from -18.81 and -15.26 kcal/mol to -12.92 kcal/mol of reaction 2bb′. This is mainly due to, on one hand, the smaller steric effect of the monodentate acetate in the first coordination shell and, on the other hand, stabilization of the H bonds composed of the free acetate O atom in the monodentate form with the water H atoms. 4.3. Implications to the Actual Ca2+ Selectivity of Ca2+ATPase. The calculated results have revealed that the amino acid residues on the surface of the Ca2+-ATPase, in spite of the different size and charges, prefer to bind Mg2+ rather than Ca2+ because of the bigger electrostatic interaction between the magnesium cation and the coordinated ligands (all calculated ∆∆G values are negative), which is consistent with the experimental conclusion that at the physiological temperature the Mg2+-ligand interactions are stronger than Ca2+-ligand ones.31 Thus, these residues can discriminate Ca2+ from Mg2+ effectively. Clearly, when hydrated Ca2+ and Mg2+ move to the vicinity of the Ca2+-ATPase, the protein pockets at the enzyme-water interface bind Mg2+ more tightly than Ca2+, resulting in the concentration decrease of Mg2+ in local regions. Then Ca2+ transports into the entrance channel of the Ca2+ATPase favorably, and selective transportation of Ca2+ against the higher concentration of Mg2+ is achieved. This process is illustrated in Figure 5a. Obviously, in this process H2O-ligand exchange will take place partially or fully and is energetically favored, as shown by the thermodynamic data. For example, ∆G of the reactions to describe the water-indole, waterformamide, water-monodentate acetate, and water-bidentate acetate exchanges are -1.78 (-4.89), -21.19, -187.06, and -193.05 kcal/mol, respectively. According to the above analyses, the amino acid residues located on the surface of the Ca2+-

ATPase may play dual roles. On one hand, their existence on the hydrophilic surface may stabilize the enzyme through interaction with the environmental molecules. On the other hand, they can select the proper Ca2+ ions to perform further physiological functions. Accompanying the Ca2+ current there are probably few Mg2+ ions to go into the entrance channel, which is composed of the main-chain carbonyl oxygen and sidechain oxygen of Gln/Asn, Glu/Asp, and Thr mainly.4,5 As revealed by the data in Table 4, the formamide, modeling the side-chain oxygen of Gln/Asn and the main-chain carbonyl oxygen, and acetate (mono- or bidentate), modeling the sidechain oxygen of Glu/Asp, molecules exhibited a larger affinity to Mg2+ relative to Ca2+. As depicted by Figure 5b, when few magnesium cations enter into the entrance channel with the Ca2+ current they will be bound tightly by the amino acid ligands. Inevitably this will hamper Mg2+ transportation, but Ca2+ can transport through the transmembrane region with a high transportation rate more easily. These properties were first explained with the so-called “sticky-pore” hypothesis32 in which ions that are bound with higher affinity pass through the channel more slowly and so have a lower transportation. Other very important factors to distinguish Ca2+ and Mg2+ are the size of the cavities formed by the binding sites and the preferred coordination number of these two cations. In the X-ray structure of the SR Ca2+-ATPase two Ca2+ high-affinity sites I and II exist (Figure 1b). Both of the sites contain seven coordination oxygen atoms, some of which bear only a partially negative charge (the main-chain carbonyl; the side-chain carbonyl or hydroxyl; the water molecules) and others bear a fully negative charge on the carboxyl groups. Analyses of the crystal structure revealed that the mean Ca2+-O distances of the two sites are 2.36 and 2.42 Å, consistent with the average heptacoordinated Ca2+-O bond length (2.41 Å) from the Cambridge Structure Database (CSD) analysis,33 while the average hexacoordinated Mg2+-O bond length is 2.07 Å.33 Thus, the high-affinity sites are good matches to the Ca2+ ion. The above factors elucidated the origins of size selectivity of the cations. In the cytoplasm the phenomenon of Ca2+ selectivity against much higher background concentrations of monovalent Na+ and K+ is also found,2,34 which arouses much interest to investigate the valence selectivity of the cations.35-37 Most of these studies focused on three aspects (for simplifica-

12292 J. Phys. Chem. B, Vol. 111, No. 42, 2007 tion, only the comparison between Ca2+ and Na+ is performed). First, it is proposed that the selectivity of Ca2+ ions is a result of ions competing to achieve charge neutrality in a selectivity filter with finite space. In this model Ca2+ ions have the same charge-neutralizing effect as two Na+ ions while occupying less of the limited volume of the filter. Essentially, not enough Na+ ions can be squeezed into the filter to achieve charge neutrality.35,36 Second, the Brownian dynamics simulations performed by Corry et al.37 revealed that the electrostatic attraction of the protein is all that accounts for ion permeation and selectivity. In the transportation process repulsion between two resident Ca2+ ions is found to speed up their exit, as illustrated in Figure 5c. Because the divalent calcium ions are more strongly attracted by the channel they can displace the sodium ions to occupy this region. Once there, Ca2+ can only be moved by the repulsion from another divalent ion and not by the lower repulsion from Na+ (Figure 5d). Third, the preferred coordination number of Na+ (6)38 and the average Na+-O bond length (2.40 Å)39 are unmatched to the Ca2+ channel. 5. Conclusions We investigated the nature of the different amino acid residues interacting with the biologically active Ca2+ and Mg2+ and the selective transportation of Ca2+ over Mg2+ at the enzyme-water interface and in the entrance channel of the SR Ca2+-ATPase at the DFT level of theory. The calculated results demonstrate that the electronegative protein cavities, composed of either aliphatic or aromatic amino acids, at the enzyme-water interface and in the entrance channel prefer to bind Mg2+ rather than Ca2+ due to the stronger electrostatic interactions between Mg2+ and ligands. Further investigations show that the affinities of amino acids to the hydrated cations depend on the nature of the metal cations and the electronegativity and binding mode (the first- vs second-shell binding; mono- vs bidentate binding) of the amino acid ligands. In detail, relative to Ca2+, Mg2+ has a stronger electrostatic interaction to the same ligands because of its smaller radius and bigger charge density. On the other hand, when the electronegativity of the ligand increases, the ligand exhibits a bigger affinity to Mg2+. Simultaneously, discrimination of the ligands between Ca2+ and Mg2+ is also enhanced. For the studied cases, the order of the ligand discrimination between Ca2+ and Mg2+ are as follows: TrpGlu > Trp-Gln > Trp-Trp > Trp. In addition, the ligand binding modes have essential effects on ligand discrimination. In the Trp and Trp-Trp cases reactions with the biggest ion selectivity are those which exhibit the second-shell indole binding, while for the Trp-Gln and Trp-Glu cases reactions with the biggest ion selectivity are those whose products have the second-shell indole and the first-shell monodentate aliphatic amino acids. More importantly, the calculations supplied feasible mechanisms to explain the size selectivity of Ca2+ over Mg2+. Concretely to say, when these two kinds of hydrated cations move to the vicinity of the Ca2+-ATPase, the protein pockets at the enzyme-water interface bind Mg2+ more tightly than Ca2+, resulting in the concentration decrease of Mg2+ in the local region. Then Ca2+ transports into the entrance channel of the Ca2+-ATPase favorably. As for the few Mg2+ ions, which may go into the entrance channel accompanying the Ca2+ current, the big affinity of the ligands to Mg2+ hampered the transportation of this cation freely, so the incompactly bound Ca2+ can transport through the transmembrane region of the Ca2+-ATPase at a relatively high rate. This phenomenon is consistent with the “sticky-pore” hypothesis reported previously. However, in the SR Ca2+-ATPase the mechanisms of Ca2+

Xiang et al. selectivity over other cations may be more complicated than the results revealed by our calculations. Furthermore, more environmental factors, such as solvent effects and long-distance interactions, may also have influences on the ion specificity of the Ca2+-ATPase. The correlative work is in progress. Acknowledgment. This work was supported by the NSFC (20633060 and 20573063 to Y.B.), NIH (grant no. GM62790 to R.I.C.), NCET (to Y.B.), and Virtual Lab for Computational Chemistry & SCC of CNIC-CAS, MCBILIN at MSU, HPCC at SDU. Supporting Information Available: B3LYP-optimized geometries for the Mg2+ complexes and NBO and ChelpG atomic charges on the different fragments of all products. This material is available free of charge via the Internet at http:// pubs.acs.org. References and Notes (1) Glynn, I. M. J. Physiol. 1993, 462, 1. (2) Drake, S. K.; Lee, K. L.; Falke, J. J. Biochemistry 1996, 35, 6697. (3) Yau, W. M.; Wimley, W. C.; Gawrisch, K.; White, S. H. Biochemistry 1998, 37, 14713. (4) Lee, A. G.; East, J. M. Biochem. J. 2001, 356, 665. (5) Toyoshima, C.; Nakasako, M.; Nomura, H.; Ogawa, H. Nature 2000, 405, 647. (6) Frontera, A.; Quin˜onero, D.; Garau, C.; Costa, A.; Ballester, P.; Deya`, P. M. J. Phys. Chem. A 2006, 110, 9307. (7) (a) Felder, C.; Jiang, H.; Zhu, W.; Chen, K.; Silman, I.; Botti, S. A.; Sussman, J. L. J. Phys. Chem. A 2001, 105, 1326. (b) Liu, T.; Zhu, W.; Gu, J.; Shen, J.; Luo, X.; Chen, G.; Puah, C. M.; Silman, I.; Chen, K.; Sussman, J. L.; Jiang, H. J. Phys. Chem. A 2004, 108, 9400. (c) Tan. X.; Zhu, W.; Cui, M.; Luo, X.; Gu, J.; Silman, I.; Sussman, J. L.; Jiang, H.; Ji, R.; Chen, K. Chem. Phys. Lett. 2001, 349, 113. (8) (a) Dudev, T.; Lim, C. Acc. Chem. Res. 2007, 40, 85. (b) Dudev, T.; Lim, C. J. Am. Chem. Soc. 2006, 128, 10541. (9) Jensen, F. J. Am. Chem. Soc. 1992, 114, 9533. (10) (a) Hoyau, S.; Ohanessian, G. Chem. Eur. J. 1998, 4, 1561. (b) Hoyau, S.; Ohanessian, G. J. Am. Chem. Soc. 1997, 119, 2016. (11) Wyttenbach, T.; Witt, M.; Bowers, M. T. Int. J. Mass Spectrom. 2001, 243, 182. (12) Pulkkinen, S.; Noguera, M.; Rodrı´guez-Santiago, L.; Sodupe, M.; Bertran, J. Chem. Eur. J. 2000, 6, 4393. (13) Bertran, J.; Rodriguez-Santiago, L.; Sodupe, M. J. Phys. Chem. B 1999, 103, 2310. (14) Ma, J. C.; Dougherty, D. Chem. ReV. 1997, 97, 1303. (15) Ryahov, V.; Dunbar, R. C. J. Am. Chem. Soc. 1999, 121, 2259. (16) Jernigan, R.; Raghunathan, G.; Bahar, I. Curr. Opin. Struct. Biol. 1994, 4, 256. (17) Gouaux, E.; MacKinnon, R. Science 2005, 310, 1461. (18) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A., Jr.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian 03, revision B.05; Gaussian, Inc.: Wallingford, CT, 2004. (19) (a) Becke, A. D. J. Chem. Phys. 1993, 98, 5648. (b) Becke, A. D. Phys. ReV. 1988, A38, 3098. (20) Lee, C.; Yang, W.; Parr, R. G. Phys. ReV. 1988, B37, 785. (21) Becke, A. D. J. Chem. Phys. 1993, 98, 1372. The BHandHLYP keyword in Gaussian employs a functional according to the formula, BHandHLYP: 0.5Ex(HF) + 0.5Ex(LSDA) + 0.5∆Ex(Becke88) + Ec(LYP), which is not precisely the formulation proposed by A. D. Becke in his paper. (22) Frisch, M. J.; Del Bene, J. E.; Binkley, J. S.; Schaefer, H. F., III. J. Chem. Phys. 1986, 84, 2279.

Ca2+ Selectivity of the Sarcoplasmic Reticulum (23) (a) Glendening, E. D.; Reed, A. E.; Carpenter, J. E.; Weinhold, F. NBO Version 3.1. (b) Reed, A. E.; Curtiss, L. A.; Weinhold, F. Chem. ReV. 1988, 88, 899. (24) Breneman, C. M.; Wiberg, K. B. J. Comput. Chem. 1990, 11, 361. (25) CRC Handbook of Chemistry and Physics, 82nd ed.; Lide, D. R., Ed.; CRC Press: Cleveland, 2001-2002. (26) Marcus, Y. Chem. ReV. 1988, 88, 1475. (27) Alabugin, I. V.; Manoharan, M.; Peabody, S.; Weinhold, F. J. Am. Chem. Soc. 2003, 125, 5973. (28) de Vos, A. M.; Ultshc, M.; Kossiakoff, A. A. Science 1992, 255, 306. (29) Pflugrath, J. W.; Wigand, G.; Huber, R.; Vertesy, L. J. Mol. Biol. 1986, 189, 383. (30) Rardin, R. L.; Tolmann, W. B.; Lippard, S. J. New J. Chem. 1991, 15, 417. (31) Rodriguez-Cruz, S. E.; Jockusch, R. A.; Williams, E. R. J. Am. Chem. Soc. 1999, 121, 1986. (32) Bezanilla, F.; Armstrong, C. M. J. Gen. Physiol. 1972, 60, 588. (33) Katz, A. K.; Glusker, J. P.; Beebe, S. A.; Bock, C. W. J. Am. Chem. Soc. 1996, 118, 5752.

J. Phys. Chem. B, Vol. 111, No. 42, 2007 12293 (34) (a) Falke, J. J.; Drake, S. K.; Hazard, A. L.; Peersen, O. B. Q. ReV. Biophys. 1994, 27, 219. (b) Drake, S. K.; Zimmer, M. A.; Kundrot, C.; Falke, J. J. J. Gen. Physiol. 1997, 110, 173. (c) Drake, S. K.; Zimmer, M. A.; Miller, C. L.; Falke, J. J. Biochemistry 1997, 36, 9917. (35) Nonner, W.; Catacuzzeno, L.; Eisenberg, B. Biophys. J. 2000, 79, 1976. (36) Boda, D.; Busath, D.; Henderson, D.; Sokolowski, S. J. Phys. Chem. B 2000, 104, 8903. (37) (a) Corry, B.; Allen, T. W.; Kuyucak, S.; Chung, S. H. Biochim. Biophys. Acta 2000, 1509, 1. (b) Corry, B.; Allen, T. W.; Kuyucak, S.; Chung, S. H. Biophys. J. 2001, 80, 195. (c) Corry, B.; Vora, T.; Chung, S. H. Biochim. Biophys. Acta 2005, 1711, 72. (d) Corry, B.; Chung, S. H. Cell. Mol. Life Sci. 2006, 63, 301. (38) Dudev, M.; Wang, J.; Dudev, T.; Lim, C. J. Phys. Chem. B 2006, 110, 1889. (39) Dobler, M. Ionophores and Their Structures; Wiley: New York, 1981.