Metal-Binding Affinity and Selectivity of Nonstandard Natural Amino

11754

J. Phys. Chem. B 2009, 113, 11754–11764

Metal-Binding Affinity and Selectivity of Nonstandard Natural Amino Acid Residues from DFT/CDM Calculations Todor Dudev† and Carmay Lim*,†,‡ Institute of Biomedical Sciences, Academia Sinica, Taipei 115, Taiwan, and the Department of Chemistry, National Tsing Hua UniVersity, Hsinchu 300, Taiwan ReceiVed: May 7, 2009; ReVised Manuscript ReceiVed: July 4, 2009

Unnatural amino acid residues are increasingly being used in metalloprotein design and engineering to expand the repertoire of protein structures/folds and functions. However, natural but nonstandard amino acid residues (not in the basic set of 20) possessing metal-ligating groups such as selenocysteine (Sec), pyrrolysine (Pyl), and γ-carboxyglutamic acid (Gla) have attracted little attention, and their potential as metal-binding entities in metalloprotein engineering has not been assessed. In particular, the metal-binding affinity/selectivity of these three rare residues remains unclear. Herein, the metal-binding affinity/selectivity of the Gla, Pyl, and Sec side chains have been systematically studied using a combined density functional theory and continuum dielectric method. The calculations reveal an advantage of using these noncanonical protein building blocks instead of the standard 20 amino acid residues. Gla2-, Pyl0, and Sec- have greater potential in trapping the metal cation than their standard amino acid counterparts. They prefer binding to Zn2+ rather than to Mg2+ or Ca2+ in a protein cavity due to the better electron-accepting ability and lower coordination number preference of Zn2+, as compared to Mg2+ and Ca2+. Between Ca2+ and Mg2+, Gla2- prefers Ca2+, whereas Pyl0 and Sec- poorly discriminate between the two metal cations. The results herein suggest that Gla2-, Pyl0, and Sec- could be employed as very efficient metal-binding entities in engineering metalloproteins with preprogrammed properties. Introduction Metal cofactors, which are found in ∼40% of all known proteins,1 are engaged in performing a plethora of tasks ranging from protein structure stabilization to enzyme catalysis, signal transduction, nitrogen fixation, photosynthesis, and respiration.2-5 Backbone peptide groups and the side chains of several standard amino acid residues such as Asp, Glu, Cys, His, Asn, Gln, Ser, and Thr, are most frequently found to coordinate the metal cation(s) in metalloproteins.6 Amino acid residues with negatively charged deprotonated side chains (Asp, Glu, and Cys) are often used to sequester the metal ion from the cellular fluids, stabilizing the protein structure and/or attaining proper protein fold.6-11 Moreover, as a carboxylate can bind to a metal ion monodentately (via one of its oxygen) or bidentately (via both oxygen atoms), Asp/Glu can be used as a molecular “switch” to regulate protein function.12 Amino acid residues with neutral coordinating groups (His, Asn, Gln, Ser, Thr), characterized by weaker electrostatic interactions with the metal cation than the charged amino acid residues,10 enhance the binding affinity and/ or selectivity of the active site and help fine-tune the geometrical and electrostatic parameters of the binding pocket for proper catalytic action.4,11,13 Because of the specific ligating properties of the metalcoordinating amino acid residues (Asp, Glu, Cys, His, Asn, Gln, Ser, and Thr), they have been extensively used in the design and engineering of new metal-binding sites in proteins with preprogrammed properties. As summarized in recent reviews,14-17 these endeavors pursue the following various goals: (i) increas* To whom correspondence should be addressed. E-mail: carmay@ gate.sinica.edu.tw. † Academia Sinica. ‡ National Tsing Hua University.

ing the protein affinity for a given metal cation, (ii) altering the selectivity of the metal-binding site so that it can bind another (very often heavy and toxic) metal type, (iii) enhancing the mechanical and thermal stability of the protein, (iv) modulating the catalytic activity of the active site, (v) altering the structure/ fold/function of the host protein, and (vi) influencing the protein-protein or protein-DNA interactions. To expand the spectrum of properties of the engineered or de novo designed metalloproteins, unnatural amino acid residues that can bind metal cations have been studied experimentally.18,19 These include bipyridylalanine,19 4-β-(pyridyl)-alanine, 1-methylhistidine, and p-aminophenylalanine.18 Interestingly, naturally occurring but nonstandard amino acid residues (i.e., not in the basic set of 20) possessing metal-ligating groups such as selenocysteine (Sec), pyrrolysine (Pyl), and γ-carboxyglutamic acid (Gla) have attracted little attention,20-23 and their potential as metal-binding entities in metalloprotein engineering has not been fully assessed. Selenocysteine is recognized as the 21st amino acid residue that is coded indirectly in the genetic code.24-26 It is an analogue of cysteine (Figure 1a) where the sulfur atom is replaced by the next element down group VIA, selenium (Figure 1b). Several non-metalloenzymes exploit the specific properties of selenium in performing their catalytic tasks. These include glutathione peroxidases scavenging for destructive peroxides,27,28 thioredoxin reductases responsible for cell growth and survival in higher eukaryotes,29 and tetraiodothyronine deiodinases involved in the synthesis and regulation of essential thyroid hormones.26 As an example, using selenium instead of sulfur improves the catalytic rate of formate dehydrogenases by 300 times.26 A density functional theory (DFT) study on the Co2+ binding to L-cysteine and L-selenocysteine in the gas phase30 found that replacing sulfur by selenium results in a very small difference in the Co2+

10.1021/jp904249s CCC: $40.75  2009 American Chemical Society Published on Web 07/30/2009

Metal-Binding Affinity of Rare Natural Amino Acids

J. Phys. Chem. B, Vol. 113, No. 34, 2009 11755 by a second carboxylate group (Figure 1f). Proteins rich in Gla are coagulation factors II (prothrombin), VII, IX, and X, protein S, osteocalcin (noncollagenous protein found in bone and dentin), and conotoxins (peptide toxins found in Conus snail venom; e.g., conantokin-G and conantokin-T).35-40 The presence of a second carboxylate in Gla (Figure 1f) increases the affinity of the side chain for metal cations so, not surprisingly, Gla is utilized by Nature mostly for binding divalent (rather than monovalent) cations, Ca2+ and Mg2+, during the blood coagulation process. In vitro experiments have shown that Gla-rich proteins/peptides can bind other divalent or trivalent metal cations as well, such as Ba2+, Mn2+, Zn2+, Cd2+, La3+, Eu3+, Gd3+ and Tb3+.41-43 Some early studies44,45 have studied malonate (modeling Gla) interacting with Ca2+ and Mg2+, using the Hartree-Fock method with the 3-21G* basis. However, the major determinants of the metal cation selectivity of Gla remain elusive. In this work, we have systematically studied the metal-binding properties of the Gla, Pyl, and Sec side chains using DFT to treat the metal cations and the first-shell ligands, in combination with a continuum dielectric representation of the rest of the protein (see Methods section). Gla, Pyl, and Sec with side chain pKa values of 3.2/4.8,46 ∼7,47 and 5.4,48 respectively, are expected to have a net charge of -2, 0, and -1, respectively, in the vicinity of a metal cation around physiological pH. The free energies of Gla2-, Pyl0, and Sec- side chains binding to Mg2+, Ca2+, and Zn2+ were evaluated in the gas phase and in a protein environment. The effects on the metal-binding affinity/ selectivity of (1) the coordinating group’s denticity, protonation state, and charge-donating properties, (2) the metal’s electronaccepting ability, and (3) the binding cavity’s dielectric constant were assessed. The metal-binding affinities of Gla2-, Pyl0, and Sec- were compared with those of the respective standard amino acid side chains. The results reveal the key factors determining the cation binding affinity/selectivity in these systems. Methods

Figure 1. Standard amino acid residues (left) and their nonstandard counterparts (right).

binding affinity. Apart from this work, no other computational studies on the metal-binding affinity/selectivity of selenocysteine have been reported (to the best of our knowledge). In 2002, the 22nd amino acid residue, pyrrolysine, was discovered.31-33 This nonstandard naturally occurring amino acid residue, like Sec, is coded indirectly in the genetic code.32,33 Pyrrolysine is present in the active site of enzymes employed by some methanogenic archaea in their methane-generating metabolism. It is a derivative of lysine (Figure 1c) where one of the side chain amino hydrogen atoms is replaced by 3-methyl1-pyrrolinylcarbonyl (Figure 1d). The side chain of Pyl has two potential metal-binding centers: the carbonyl oxygen and the pyrrolinyl ring nitrogen. However, their metal-coordinating properties have not been explored experimentally or theoretically. Unlike selenocysteine and pyrrolysine, the γ-carboxyglutamic acid, another rare amino acid residue, is not genetically encoded. Instead, it is produced by post-translational carboxylation of specific Glu residues (Figure 1e) by a vitamin K-stimulated reaction,34 in which a side chain methylene hydrogen is replaced

The metal cation and ligands forming the first coordination shell were treated quantum mechanically to account for electronic effects such as polarization of the participating entities and charge transfer from the ligands to the metal cation. The rest of the protein was considered as a dielectric continuum characterized by a dielectric constant ranging from 4 for totally buried metal-binding pockets to 80 for fully solvent exposed sites. Models Used. In modeling the amino acid residues, the side chains of Cys, Sec, His, Asp, Glu, and Gla were truncated at the CR position and the H2N-HCR-COOH moiety was converted into a CH3 group, whereas the side chain of Pyl was truncated at the Cδ position and capped with a CH3 group (Figure 2). The backbone peptide group was modeled by CH3CONHCH3. Asp, Glu, Gla and Sec with side chain pKa values of 3.9,49 4.3,49 3.2/4.8,46 and 5.4,48 respectively, are expected to be deprotonated at physiological pH, so the corresponding models were assumed to be monoanionic except Gla, which is dianionic. At neutral pH, the metal-free Cys side chain with a typical pKa value between 8 and 950-52 would be protonated, so it was modeled as CH3CH2SH. Binding to a metal cation (acting as a Lewis acid) causes the cysteine pKa to drop,53 thus facilitating sulfhydryl group deprotonation under physiological conditions. Hence, Cys was ionized in the respective metal complexes. His and Pyl side chains, whose imidazole and 1-pyrrolinyl heterocyclic rings, respectively, are characterized by a pKa of ∼7,47,54 were assumed to be neutral at ambient pH in the protein interior.

11756

J. Phys. Chem. B, Vol. 113, No. 34, 2009

Dudev and Lim TABLE 1: Comparison between Computed and Experimental Average Metal-Ligand (M-L) Bond Distances (in Å) molecule +c

[Li(H2O)4] [Mg(H2O)6]2+ c [Ca(H2O)7]2+ d [Zn(H2O)6]2+ c [Zn(NH3)4]2+ d bis(ethane-1,2dithiolato-S-S′)-zinc(II)e bis(tetraselenidoSe1,Se4)-zinc(II) d [Cd(H2O)6]2+ d [Cd(NH3)6]2+ d [Hg(SCH3)2]0 d [La(H2O)9]3+ c [ZrCl6]2- c

M-L

RM-L (calc)a

RM-L (expt)b

Li-O Mg-O Ca-O Zn-O Zn-N Zn-S

1.93 2.05 2.42 2.04 2.01 2.34

1.94 ( 0.05 2.07 ( 0.03 2.40 ( 0.04 2.08 ( 0.03 2.02 ( 0.02 2.34 ( 0.02

Zn-Se

2.46

2.47 ( 0.02

Cd-O Cd-N Hg-S La-O Zr-Cl

2.25 2.38 2.32 2.55 2.48

2.27 ( 0.04 2.37 ( 0.03 2.36 ( 0.05 2.55 ( 0.04 2.47 ( 0.02

a

Calculated at the S-VWN/6-31+G* level with the SDD basis set for Se, Cd, Hg, La, and Zr. b From Cambridge Structural Database (CSD) analysis. c From Dudev and Lim.58 d This work. e From Dudev and Lim.59

Figure 2. Molecules modeling side chains of (a) Cys0, (b) Sec-, (d) His0, (e) Pyl0, (f) Asp-, (g) Glu- and (h) Gla-2, and (c) the backbone group.

In proteins, Mg2+ is predominantly hexacoordinated, as in aqueous solution, but Ca2+ is usually heptacoordinated.55,56 Thus complexes containing Mg2+ and Ca2+ were modeled as MgL6 and CaL7 (L ) H2O or an amino acid residue), respectively, in aqueous solution and in proteins. On the other hand, Zn2+ is more variable in its coordination preference and is usually hexacoordinated in solution but tetracoordinated in proteins.55,56 Accordingly, Zn2+ complexes were modeled as [Zn (H2O)6]2+ in aqueous solution and ZnL4 in proteins. DFT Calculations. (a) Geometries. Full geometry optimization for each ligand/metal complex studied was carried out using the Gaussian 03 program57 employing the S-VWN functional. The 6-31+G* basis set was used for all the atoms except Se where the SDD basis set was employed. This functional/basis set combination was chosen as it reproduces the experimentally observed metal-ligand bond distances in a number of metal-ligand complexes within experimental error (see Table 1). For each fully optimized structure, S-VWN/(6-31+G*,SDD) vibrational frequencies were computed to verify that the molecule was at the minimum of its potential energy surface. No imaginary frequency was found in any of the metal complexes. Among several possible conformations of a given metal complex, the one with the lowest energy was selected for further evaluation (see Figures S1 and S2 in the Supporting Information and the Results section). (b) Charge Transfer. The difference in the amount of charge transferred from a given ligand to two different ions or from two different ligands to the same ion was assessed by evaluating S-VWN/(6-31+G*,SDD) atomic charges obtained from MerzKollman (MK) population analysis.60,61 To ensure that the trend in the charge transfer difference is independent of the charge scheme employed, the difference in the amount of charge transferred from F to hydrogen in HF and to a metal ion in LiF or NaF (denoted by ∆CTmetal) was computed using various charge analyses; viz., MK, CHelpG, NBO, Hirshfeld, and Voronoi deformation density (VDD).62 All the charge transfer schemes show the same trend in the ∆CT, which is positive,

TABLE 2: Formation Energies, ∆Eform, for Mg2+ + 6H2O f [Mg (H2O)6]2+ Evaluated at Different Levels of Theory (in kcal/mol) method

∆Eform (in kcal/mol)

B3-LYP/6-311++G(2df,2p) MP2/6-311++G(2df,2p) MP4SDTQ/6-311++G(2df,2p) CCSD(T)/6-311++G(2df,2p) QCISD(T)/6-311++G(2df,2p)

-317.6 -318.8 -319.1 -319.8 -319.7

indicating that the charge transfer from F to H is greater than that to Li or Na: The ∆CTLi is 0.38e using MK or Hirshfeld, 0.37e using CHelpG, and 0.34e using NBO or VDD, whereas the ∆CTNa is 0.42e using MK or Hirshfeld or VDD, 0.39e using CHelpG, and 0.37e using NBO (see Table S1 in the Supporting Information). Furthermore, all the charge transfer schemes show that the charge transfer from F to Li is greater than that to Na (by 0.02 to 0.08e). (c) Gas-Phase Free Energies. Based on the fully optimized S-VWN/(6-31+G*,SDD) geometries, the electronic energies, Eelec, were evaluated using the B3-LYP functional with a large 6-311++G(2df,2p) basis set. The latter was chosen from among several other basis sets as the gas-phase formation energy of [Mg(H2O)6]2+, ∆Eform, has been found to converge at this level of theory.8 Moreover, the B3-LYP/6-311++G(2df,2p) method produces ∆Eform that is close (within 0.7%) to ∆Eform evaluated by more accurate but computationally expensive post-HF methods such as second-order Moller-Plesset MP2, fourth-order Moller-Plesset with single, double, triple and quadruple substitutions MP4SDTQ, coupled cluster with single, double and triple substitutions CCSD(T), and quadratic configuration interaction with single, double and triple substitutions QCISD(T) (see Table 2). The basis set superposition errors were shown not to be significant for this large basis set63 and were thus omitted in this work. The thermal energy, including zero point energy (ET), and entropy (S) corrections, were evaluated according to standard statistical mechanical formulas64 using S-VWN/(6-31+G*,SDD) vibrational frequencies scaled by an empirical factor of 0.9833.65 The differences ∆Eelec, ∆ET, ∆PV (work term), and ∆S between the products and reactants were used to compute the reaction free energy in the gas phase at room temperature, T ) 298.15 K, according to the following expression:

Metal-Binding Affinity of Rare Natural Amino Acids

∆G1 ) ∆Eelec + ∆ET + ∆PV - T∆S

J. Phys. Chem. B, Vol. 113, No. 34, 2009 11757

(1)

Continuum Dielectric Calculations. The reaction free energy in a given environment characterized by a dielectric constant ε ) x can be calculated according to the following thermodynamic cycle:

∆G1, the gas-phase free energy, was computed using eq 1, as described above. ∆Gsolvx, the free energy for transferring a molecule in the gas phase to a continuous solvent medium characterized by a dielectric constant, x, varying from 4 to 80, was estimated by solving Poisson’s equation using finite difference methods.66,67 Thus, the reaction free energy in an environment modeled by dielectric constant x, ∆Gx, can be computed from

∆Gx ) ∆G1 + ∆Gsolvx(products) - ∆Gsolvx(reactants) (2) The solvation free energies were evaluated using the MEAD (Macroscopic Electrostatics with Atomic Detail) program,68 as described in previous works.58 The effective solute radii, which were obtained by adjusting the CHARMM (version 22)69 van der Waals radii to reproduce the experimental hydration free energies of Mg2+, Ca2+, Zn2+, and model ligand molecules, are as follows (in Å): RMg ) 1.50, RCa ) 1.75 RZn ) 1.40, RC ) 1.88, RS(SH) ) 2.5, RS(S-) ) 2.19, RSe(Se-) ) 2.55, RN(CONH) ) 1.75, RN(His/Pyl) ) 1.80, RO(CONH) ) 1.78, RO(COO-) ) 1.575, RO(H2O) ) 1.78, RO(H2O-Mg/Ca/Zn) ) 1.70, RH(C/S/ N) ) 1.468, RH(H2O-Mg) ) 1.16, RH(H2O-Ca) ) 0.99, RH(H2O-Zn) ) 1.09. These effective solute radii reproduce (within 1 kcal/mol) the experimental hydration free energies of the cations and model protein ligands.9,58 Results γ-Carboxyglutamic Acid (Gla). (a) Preferred Binding Mode. Each of the two carboxylate groups from the Gla2- side chain can bind the metal cation either monodentately or bidentately. Combining these binding modes for the two carboxylate groups yields several possibilities for coordinating the metal cation (see Figure S1 in the Supporting Information): (a) monodentate binding of one of the carboxylates, (b) bidentate binding of one of the carboxylates, (c) monodentate binding of both carboxylates (so-called chelation bidentate mode), (d) monodentate binding of one carboxylate and bidentate binding of the other, and (e) bidentate binding of both carboxylates. To determine which of these binding mode possibilities is favored for a given metal ion, geometry optimization of the metal cation bound to Gla2- in all 5 possible modes was performed, as shown in Figure S1 in the Supporting Information. For the Mg2+ and Ca2+ complexes, fully optimized geometries could be obtained starting from the single-carboxylate monodentate (Figure 3a) and chelation bidentate (Figure 3b) binding modes, but not for the other three binding modes, which converged to the singlecarboxylate monodentate or chelation bidentate structure during geometry optimization (see Figure S1b,d,e in the Supporting

Figure 3. Fully optimized S-VWN/6-31+G* geometries of Mg2+ (yellow) and Ca2+ (blue) bound to Gla2- (a) via a single carboxylate group (monodentate) or (b) via both carboxylates (chelation bidentate) and the formation free energies (in kcal/mol) of the complexes relative to those of the respective monodentate complexes in the gas phase and in a protein environment characterized by a dielectric constant ε ranging from 4 to 20.

Information). For the Zn2+ complexes, a stationary point could be found for the chelation bidentate structure only. To determine the preferred binding mode of the Gla2- side chain for a given metal cation, the formation free energies of hexahydrated Mg2+ or heptahydrated Ca2+ binding to Gla2- in the chelation bidentate mode were compared to those of Mg2+ or Ca2+ binding to Gla2- monodentately via one of its carboxylate groups. The relative formation free energies, ∆∆Gx (x ) 1, 4, 10, and 20) in Figure 3 show that the chelation bidentate structures (with negative ∆∆Gx) are more stable than their reference monodentate counterparts in the gas phase and in a protein cavity. This finding is in agreement with earlier gasphase theoretical results,44 which show that the chelation

11758

J. Phys. Chem. B, Vol. 113, No. 34, 2009

Dudev and Lim

Figure 4. Free energies (in kcal/mol) of Gla2- binding in the chelation bidentate mode to (a) Mg2+ (yellow), (b) Ca2+ (blue) and (c) Zn2+ (black) in the gas phase and in a protein cavity.

bidentate binding in Mg2+/Ca2+-malonate complexes (with malonate modeling the Gla side chain) is preferred over monodentate coordination. It is also in accord with surveys of crystal structures of metal-malonate44 and metal-Gla complexes (this work; PDB entries 1NL1, 1NL2, 1J34, 1J35, and 1VZM) showing that the favored mode of Gla binding to a single Mg2+ or Ca2+ ion is the chelation bidentate mode. Three factors mainly contribute to the higher stability of the chelation bidentate complexes relative to the respective monodentate structures: (1) Replacing a water molecule with a second carboxylate ligand, which can donate more charge than H2O to the metal cation,59 stabilizes the chelation bidentate complexes. (2) Releasing an extra water molecule upon chelation bidentate binding yields a more favorable gas-phase entropy and solvation free energy, as compared to single-carboxylate monodentate binding. (3) As the solvation free energy depends inversely on the size of the molecule, the more compact chelation bidentate structure is better solvated than the single-carboxylate monodentate counterpart. (b) Metal SelectiWity. To assess the metal selectivity of Gla2-, its binding free energies to Mg2+, Ca2+, and Zn2+ in the stable chelation bidentate mode were computed in the gas phase and in a condensed medium (Figure 4). Due to the attractive electrostatic interactions between the dication and dianionic Gla2-, the gas-phase free energies for all three dications are very negative but attenuate quickly upon solvation such that binding of Gla2- to Mg2+ or Ca2+ in a partially solvent-exposed cavity (ε g 20) becomes unfavorable (positive ∆G20 in reactions a and b in Figure 4). Among the three dications, Zn2+, which is known to be the better Lewis acid,59,70 seems to be the preferred

cation for Gla2-. The ∆Gx (x ) 1-20) values for the Zn2+ complex (Figure 4c) are more favorable (more negative) than the corresponding values for the Mg2+ (Figure 4a) and Ca2+ (Figure 4b) complexes. That Zn2+ is a better electron-acceptor than Mg2+ and Ca2+ is evidenced by its lower charge (1.09e) compared to the charge on Mg2+ (1.27e) or Ca2+ (1.55e) in the respective complexes. In addition, the release of two extra water molecules upon Zn2+ binding compared to Mg2+ or Ca2+ binding also contributes favorably to the Zn2+-Gla interaction. Furthermore, the tetracoordinated Zn2+ complex, being more compact, is better solvated than the corresponding hexacoordinated Mg2+ or heptacoordinated Ca2+ complex. Between Mg2+ and Ca2+, two of the physiologically most important cations, which ion does Gla2- prefer to bind to? Mg2+, possessing higher charge density and better charge-accepting capacity than Ca2+,70 generally binds with higher affinity to a given ligand than Ca2+.10,76 In line with this expectation, the charge on Mg2+ is less than that on Ca2+ (see above), indicating more charge transfer from Gla2- to Mg2+ than to Ca2+ (by 0.28e). Yet, the γ-carboxyglutamic acid side chain interacts more favorably with Ca2+ than with Mg2+ in the gas phase and in a buried protein cavity. The ∆Gx (x ) 1-10) values for Gla2binding to Ca2+ (Figure 4b) are more negative than those for Gla2- binding to Mg2+ (Figure 4a) by ∼3 kcal/mol in the gas phase and by ∼9 kcal/mol in a buried protein cavity. This result is in line with experimental observations showing that Glacontaining proteins are tuned to bind preferentially Ca2+ in vivo.36,41,77 To unravel the major determinants of the Ca2+/Mg2+ selectivity of Gla2-, the structural parameters of the free and metal-

Metal-Binding Affinity of Rare Natural Amino Acids

J. Phys. Chem. B, Vol. 113, No. 34, 2009 11759

Figure 5. Free energies (in kcal/mol) of heptahydrated Ca2+ (blue) binding to (a) two Glu- side chains and (b) a single Gla2- side chain in the chelation bidentate mode in the gas phase and in a protein cavity.

bound Gla ligand as well as the solvation free energies of the Ca2+/Mg2+ complexes were compared. The structural comparison shows that the angle R between the two COO- groups in the [Mg (H2O)4 Gla]0 complex (100.7°) or in the [Ca (H2O)5 Gla]0 complex (101.8°) is significantly smaller than that in the free Gla2- (115.0°), implying some degree of strain in the chelation bidentate complexes. Along the same vein, the distance d between the two metal-bound oxygen atoms in the [Mg(H2O)4Gla]0 complex (2.82 Å) or in the [Ca(H2O)5Gla]0 complex (2.96 Å) is significantly shorter than the corresponding distance (3.67 Å) in the free Gla2- ligand. Notably, relative to the [Mg(H2O)4Gla]0 complex, the corresponding Ca2+ complex has a larger R and a longer d, and would thus be less strained. Furthermore, the Ca2+-Gla complex is better solvated than its Mg2+ counterpart. Hence, the calculations imply that the different degree of strain in the 6-membered rings formed upon chelation bidentate Gla2- binding to Mg2+ and Ca2+ and the different degree of solvation of the Ca2+ and Mg2+ complexes determines the competition between these two cations for the ligand. (c) Comparison with Glu- Binding. How does the metalbinding affinity of Gla2- compare with that of its common amino acid counterpart Glu-? To answer this question, the free energies of Ca2+ binding to a single Gla2- in the gas phase and in a condensed medium (Figure 5b) were compared with the respective values for Ca2+ binding simultaneously to two Gluligands with the same net charge of -2 (Figure 5a). In buried protein cavities (characterized by ε e 10), Ca2+ binding to Gla2is more favorable (more negative ∆Gx) than that of the two Glu- residues combined (see Figure 5). Hence, converting Gluto Gla2- in deeply buried sites of certain proteins enhances the affinity of these sites for metal cations, in particular, Ca2+. However, the higher metal affinity for Gla2- instead of two Gluresidues attenuates with increasing solvent exposure of the metal site owing to the greater desolvation penalty of Gla2- compared to that of two Glu-. In solvent-exposed protein cavities (with ε g 20), Ca2+ binding to Gla2- or Glu- becomes thermodynamically unfavorable (positive ∆G20). Thus, the greater stability of [Ca(H2O)5Gla]0 with respect to [Ca(H2O)5(Glu)2]0 is due, not to solvent effects, but to the following two ligand features: (1) The dianionic Gla2- moiety is a better charge donor than two monoanionic Glu- residues, as evidenced by Merz-Kollman population analysis revealing larger charge transfer from the ligands to the metal in [Ca(H2O)5Gla]0 compared to [Ca(H2O)5(Glu)2]0 (by 0.25e). (2) Coordination of one Gla2- ligand to the metal cation instead of

Figure 6. Fully optimized S-VWN/6-31+G* geometries of Mg2+ (yellow) bound to Pyl0 via (a) the side chain carbonyl oxygen (monodentate O), (b) the ring nitrogen (monodentate N), and (c) both the side chain carbonyl oxygen and ring nitrogen (chelation bidentate) and their relative formation free energies (in kcal/mol).

two Glu- ligands reduces the gas-phase entropy penalty, thus further lowering ∆G1 for reaction b compared to reaction a in Figure 5. Pyrrolysine (Pyl). (a) Preferred Binding Mode. Three metalbinding modes are possible for the Pyl0 side chain: (1) monodentate via the carbonyl oxygen (Figure 6a); (2) monodentate via the 1-pyrrolinyl ring nitrogen (Figure 6b); and (3) chelation bidentate via both oxygen and nitrogen atoms (Figure 6c). To determine the preferred binding mode of the Pyl0 side

11760

J. Phys. Chem. B, Vol. 113, No. 34, 2009

Dudev and Lim

Figure 7. Free energies (in kcal/mol) of Pyl0 binding in the chelation bidentate mode to (a) Mg2+ (yellow), (b) Ca2+ (blue) and (c) Zn2+ (black) in the gas phase and in a protein cavity.

chain for a given metal cation, the relative formation free energies of Mg2+ binding to Pyl0 in the three different modes were compared (Figure 6). Such a comparison shows that the chelation bidentate binding mode is preferred to the monodentate modes in a fully or partially buried protein cavity: the ∆Gx (x ) 1, 4, 10, and 20) values for 6c are more negative than those for 6a or 6b. In contrast, the two monodentate complexes are nearly equally stable in a protein cavity, differing by e1.6 kcal/ mol, although coordination to the 1-pyrrolinyl ring nitrogen appears more favorable than coordination to the carbonyl oxygen in the gas phase (∆∆G1 ) -3.2 kcal/mol; Figure 6b). Unfortunately, no appropriate PDB or CSD structures are available (to the best of our knowledge) to verify these findings. As for the Mg/Ca-Gla complexes, the greater stability of the chelation bidentate complex, 6c, relative to the monodentate complexes, 6a and 6b, is largely due to the following two factors: (1) The carbonyl oxygen or ring nitrogen in the bidentate complex donates more charge to the metal cation than the respective substituted water molecule in the monodentate complexes.10 (2) Chelation bidentate binding results in a more compact structure and an extra free water molecule, which in turn yields more favorable gas-phase entropy and solvation free energy contributions than the corresponding monodentate carbonyl oxygen or 1-pyrrolinyl ring nitrogen binding. (b) Metal SelectiWity. To assess the metal selectivity of Pyl0, its binding free energies to Mg2+, Ca2+, and Zn2+ in the stable chelation bidentate mode were computed in the gas phase and in a condensed medium (Figure 7). In analogy to Gla2-, pyrrolysine prefers binding to Zn2+ rather than to Mg2+ or Ca2+ in a protein cavity due to (i) the better electron-accepting ability of Zn2+, (ii) the release of two extra water molecules upon Zn2+ binding, and (iii) the better solvation of the smaller tetracoordinated Zn2+ complex as compared to the respective hexaco-

ordinated Mg2+ and heptacoordinated Ca2+ complexes. The ∆Gx (x ) 4-20) values for the Zn2+ complex (Figure 7c) are more favorable (more negative) than the corresponding values for the Mg2+ (Figure 7a) and Ca2+ (Figure 7b) complexes. In contrast to a fully or partially buried protein cavity, Pyl0 prefers binding to Mg2+ rather than Zn2+ in the gas phase even though Zn2+ is known to be a better Lewis acid than Mg2+. To elucidate why the gas-phase binding free energy of Pyl0 to Zn2+ (∆G1 ) -32.3 kcal/mol) is less favorable than that to Mg2+ (∆G1 ) -34.9 kcal/mol), the formation free energies of Mg2+ (Figure 8a) and Zn2+ (Figure 8c) bidentately bound to Pyl0 were compared with those of Mg2+ and Zn2+ monodentately bound to two separate ligands, each of which contains one of the metalligating entities comprising the Pyl0 side chain (see Figures 8b and 8d); viz., the carbonyl oxygen in CH3CONHCH3 and the ring nitrogen in 3-methyl-1-pyrroline. Since the metal-ligating groups belong to separate molecules, the resulting structures are strain-free and may be used as a reference for the respective strained bidentate Pyl0 complexes. Such a comparison shows that, in the gas phase, binding of Pyl0 to Zn2+ is less favorable than that to Mg2+ because of the greater ring strain in the Zn2+ complex than in the corresponding Mg2+ complex. The loss of free energy upon Pyl0 bidentate coordination due to ring strain, ∆∆G1(M))∆G1[M(bidentate-Pyl)2+]-∆G1[M(CH3CONHCH3)(3methyl-1-pyrroline)2+] (M ) Zn or Mg), is 39.1 kcal/mol for Zn2+ binding but 29.8 kcal/mol for Mg2+ binding. This difference outweighs the greater charge transfer from Pyl0 to Zn2+ than to Mg2+ (by 0.27e), thus the gas-phase free energy for Pyl0 binding to Zn2+ is less negative than that for Pyl0 binding to Mg2+. Whereas, in proteins, Gla2- prefers binding to Ca2+ than to Mg2+ (see above), Pyl0 is not effective in discriminating between Mg2+ and Ca2+, as the respective ∆Gx (x ) 1-20) values differ

Metal-Binding Affinity of Rare Natural Amino Acids

Figure 8. Fully optimized SVWN/6-31+G* geometries of (a) Mg2+ (yellow) bound to Pyl0 in the chelation bidentate mode, (b) Mg2+ bound to CH3CONHCH3 and methyl-1-pyrroline, (c) Zn2+ (black) bound to Pyl0 in the chelation bidentate mode, and (d) Zn2+ bound to CH3CONHCH3 and methyl-1-pyrroline.

by e1 kcal/mol (Figures 7a and 7b). The larger charge transfer from Pyl0 to Mg2+ than to Ca2+ (by 0.34e) is offset by the greater ring strain in the Mg2+-Pyl complex than in the corresponding Ca2+ complex: The loss of free energy due to ring strain in the Ca2+-Pyl complex, ∆∆G1(Ca), is smaller (by 14.8 kcal/mol) than ∆∆G1(Mg) in the corresponding Mg2+-Pyl complex (see above and Figure 8). Owing to the similar degree of solvation of the two metal complexes, the trends observed in the gas phase remain unaltered, thus Pyl0 binding to Mg2+ and Ca2+ in a protein environment appears almost equally favorable, as in the gas phase. (c) Comparison with Backbone and His Binding. To compare the metal-binding affinity of Pyl0 with that of its common amino acid counterparts, the free energies of Zn2+ binding to Pyl0 in a chelation bidentate mode (Figure 9b) were compared with those of Zn2+ coordinating simultaneously to CH3CONHCH3, corresponding to the Pyl’s CdO group, and the His side chain, the closest relative to the Pyl’s 1-pyrrolinyl heterocycle (Figure 9a). In a protein cavity, Zn2+ prefers to bind to Pyl0 than to the CH3CONHCH3 backbone and the His side chain due to solvation effects: The more compact [Zn(H2O)2Pyl]2+ complex is better solvated than its [Zn(H2O)2(BKB)(His)]2+ counterpart. Furthermore, the desolvation penalty for the single Pyl0 ligand (reaction b in Figure 9) is less than that for the BKB and that for His combined (reaction a in Figure 9). Thus, in a protein cavity (with ε ) 4-20), Zn2+ coordinating to Pyl0 appears more favorable than Zn2+ binding to His and backbone ligands. However, in the absence of the protein and water, Zn2+ binding to Pyl0 is less favorable than that to the CH3CONHCH3 backbone and the His side chain: The ∆G1 for reaction b in Figure 9 is less negative than that for reaction a in Figure 9 (by 8 kcal/mol). This is mainly due to two factors: (i) larger charge transfer from the ligands to Zn2+ in [Zn(H2O)2(BKB)(His)]2+ complex (by 0.30e) than in its [Zn(H2O)2Pyl]2+ counterpart, and (ii) lack of ring strain in the former structure. Selenocysteine. In contrast to Gla2- and Pyl0, which can bind to a given metal ion in various modes, Sec- has only one donor atom, the selenium atom, and can therefore bind to the metal ion only monodentately. (a) Metal SelectiWity. To assess the metal selectivity of Sec-, its binding free energies to Mg2+, Ca2+, and Zn2+ were computed in the gas phase and in a condensed medium (Figure 10). As for Gla2- and Pyl0, the ionized Sec- prefers binding to Zn2+

J. Phys. Chem. B, Vol. 113, No. 34, 2009 11761 rather than to Mg2+ or Ca2+ in the gas phase and in a protein cavity because (i) Sec- transfers more charge to Zn2+ (qZn ) 0.51e) than to Mg2+ (qMg ) 1.20e) or Ca2+ (qCa ) 1.39e), (ii) two more water molecules are released upon Zn2+ binding, as compared to Mg2+ or Ca2+ binding, and (iii) the Zn2+ complex is better solvated than the respective Mg2+ and Ca2+ complexes. The ∆Gx (x ) 1-20) values for Sec- binding to Zn2+ (Figure 10c) are much more negative than those for Sec- binding to Mg2+ (Figure 10a) or Ca2+ (Figure 10b). Moreover, Secbinding to Zn2+ can occur even in a partially solvent-exposed binding site (negative ∆G20 for reaction c in Figure 10), whereas Sec- binding to Mg2+ or Ca2+ under the same condition seems unlikely (positive ∆G20 for reaction a or b in Figure 10). These findings are in line with general expectations that the “soft” selenolate anion will have greater affinity for “borderline” Zn2+ than for “hard” Mg2+ and Ca2+. The low affinity of the “soft” selenolate anion for “hard” Mg2+ and Ca2+ also accounts for the fact that Sec- is even more selective for Zn2+, as compared to Mg2+ and Ca2+, than Gla2and Pyl0 in both fully and partially buried sites. The absolute free energy difference between the Zn-Sec and Mg/Ca-Sec complexes (Figure 10) is greater than that between the respective Zn-Gla and Mg/Ca-Gla complexes (Figure 4) or between the respective Zn-Pyl and Mg/Ca-Pyl complexes (Figure 7). For example, the ∆G10 difference between the Zn-Sec and Ca-Sec complexes (-41.5 kcal/mol, Figure 10) is greater than that between the respective Zn-Gla and Ca-Gla complexes (-20.3 kcal/mol, Figure 4) or between the respective Zn-Pyl and Ca-Pyl complexes (-11.4 kcal/mol, Figure 7). Like Pyl0, Secis not effective in discriminating between Mg2+ and Ca2+, as the respective ∆Gx (x < 10) values differ by