Quantum Chemical Studies on Stability and Chemical Activities in

Oct 30, 2015 - Our analysis indicates a general way of interpreting the electronic energy spectra to understand stability and functions of large biomo...
0 downloads 0 Views 2MB Size
Article pubs.acs.org/JPCB

Quantum Chemical Studies on Stability and Chemical Activities in Calcium Ion Bound Calmodulin Loops Samapan Sikdar,† Mahua Ghosh,*,† Molly De Raychaudhury,‡ and J. Chakrabarti*,†,§ †

Department of Chemical, Biological and Macromolecular Sciences and §Unit of Nanoscience and Technology-II and The Thematic Unit of Excellence on Computational Materials Science, S. N. Bose National Centre for Basic Sciences, Sector III, Block JD, Salt Lake, Kolkata 700098, India ‡ Department of Physics, West Bengal State University, Barasat, Kolkata 700126, India S Supporting Information *

ABSTRACT: Quantum chemical (QC) calculations for macromolecules require truncation of the molecule, highlighting the portion of interest due to heavy computation cost. As a result, an estimation of the effects of truncation is important to interpret the energy spectrum of such calculations. We perform density functional theory based QC calculations on calcium ion bound EF-hand loops of Calmodulin isolated from the crystal structure in an implicit solvent. We find that the terminal contributions of neutral capping are negligible across the entire ground-state energy spectrum. The coordination energy range and the nature of hybridization of the coordination state molecular orbitals remain qualitatively similar across these loops. While the HOMO and LUMO of loops in the N-terminal domain are dominated by the acidic aspartates, and the polar/hydrophobic residues, respectively, these levels of the C-terminal domain loops show strong localized electron density on the phenyl rings of the tyrosines. The Fukui index calculation identifies the hydroxyl oxygen in the phenyl ring of Y99 as a potent nucleophile. Our analysis indicates a general way of interpreting the electronic energy spectra to understand stability and functions of large biomolecules where the truncation of the molecule and, hence, the terminal capping effects are inevitable.



INTRODUCTION

While we have considered only the CSMOs of loop1 in the earlier study,5 we consider the entire energy spectrum of all the CaM loops in the present work. The loops 1 and 2 form the Nterminal domain (NTD) and the loops 3 and 4 form the Cterminal domain (CTD), the domains being connected by a linker region. Table 1 shows the amino acid sequences of 12residue-long EF hand loops of CaM, where the Ist, IIIrd, Vth, VIIth, and XIIth residues of the EF hand loops are involved in Ca2+ coordination along with a water molecule, leading to pentagonal bipyramidal geometry. Mutational studies11−15 highlight the roles of both ligand and nonligand residues in maintaining the coordination geometry and structural stability of the Ca2+ bound loops of CaM. Although the Ca2+ ion binding affinities of isolated CaM loops have been probed by the grafting approach,16 the loops show cooperativity in Ca2+ ion binding.17−21 However, the free energy contribution due to cooperativity is small compared to the overall binding energy.17 This leads us to consider the Ca2+ binding loops in isolation. In general, the polar amino acids, like tyrosine, serine, and threonine, are functionally significant through participation in phosphorylation, a post-translational chemical modification (PTM), at their phenyl ring hydroxyl Oh atom.22 The reversible phosphorylation of proteins plays essential role in cell signaling and various cellular processes.23−26 Phosphorylation involves

The energy spectrum of an interacting many-electron system determines a host of properties including the stability and chemical reactivity of the system.1 Quantum chemical (QC) calculations are routinely performed to derive the electronic energy spectrum. However, the QC calculations on large molecular systems of biological interest pose a formidable task due to the involvement of a large number of electrons in the constituent atoms. In these cases, the calculations are often executed on truncated systems2−5 considering only the portion of direct interest. The effects of truncated systems on the outcome of the QC calculations have been shown to be nontrivial in our earlier work5 on density functional theory (DFT)-based vacuum QC calculations of calcium (Ca2+) ion coordination to an isolated loop (loop 1) of Calmodulin (CaM), an EF-hand protein with four Ca2+ binding loops involved in various processes in eukaryotic cells.6 These calculations show that the coordination state molecular orbitals (CSMOs), derived from contributions of the valence orbitals of the Ca2+ ion and the loop atoms, yield a robust microscopic picture of coordination stability, largely independent of the conditions imposed on the terminal residues. However, the energy levels near the HOMO and the LUMO are affected by such details of the system.5,7−9 Furthermore, earlier studies8,10 show that the treatment of solvent, at least at an implicit level, is essential to calculate the HOMO and the LUMO levels with reasonable accuracy. © XXXX American Chemical Society

Received: October 5, 2015 Revised: October 28, 2015

A

DOI: 10.1021/acs.jpcb.5b09713 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Table 1. Amino Acid Sequence of the 12-Residue-Long EF Hand Loops of CaM, Where the Coordinating Residues Are Marked in Bold loop1 loop2 loop3 loop4

I

II

III

IV

V

VI

VII

VIII

IX

X

XI

XII

D20 D56 D93 D129

K21 A57 K94 I130

D22 D58 D95 D131

G23 G59 G96 G132

D24 N60 N97 D133

G25 G61 G98 G134

T26 T62 Y99 Q135

I27 I63 I100 V136

T28 D64 S101 N137

T29 F65 A102 Y138

K30 P66 A103 E139

E31 E67 E104 E140

and the C-terminus by N-methylamide (CH3NH) group. This resembles neutral terminal capping in order to satisfy the required valency of the terminal atoms and has been used in earlier works as well.4,31−33 We further add the missing hydrogen atoms in the X-ray crystal structure of each of these loops. Both are performed using VMD based on CHARMM27 topopologies.40 The resulting structures are then solvated using the TIP3P water model and neutralized with counterions. The four solvated Ca2+ binding loops are then subjected to minimization involving only the hydrogen and the terminal capping atoms keeping the other heavy atoms of the protein fixed in NAMD41 using the CHARMM27 force-field parametes.40 These optimized geometries of these loops are then considered for QC calculation. The QC-DFT calculations are done by variationally minimizing the many-body electronic Hamiltonian for fixed nuclear coordinates where the electronic correlations are approximately taken into account via general functionals of the electron density. The surrounding solvent effect has been implicitly considered through the polarizable continuum model (PCM).34 The PCM considers the molecule as embedded in a cavity surrounded by infinite dielectric medium of water (with dielectric constant of 78.39) and has been reported from earlier studies to successfully reproduce electronic properties of oligopeptides in solution.42 We perform the QC calculation using the Gaussian03 package43 with the B3LYP functional for the electron exchange-correlation in the 6-31G(d,p) basis set.44 We calculate the total ground-state energy by the Self Consistent Field method and compute the MOs. The MOs are linear combinations of atomic orbitals, the coefficients of which are referred to as the MO coefficients. In our work we perform natural population analysis (NPA) to define the partial atomic charges. We also perform all atom QC geometry optimization at the same level of theory and using the same basis set in polarizable solvent for the Ca2+ bound loop 1 from NTD and loop 3 from CTD. Fukui Index. The atom condensed Fukui index have the N − qNk and f −k = qNk − qN−1 following form: f +k = qN+1 k k , where qk represents the charge on atom k for a system of N electrons. denotes the same for a system of N − 1 electrons upon qN−1 k for a removal of an electron from HOMO and similarly qN+1 k system of N + 1 electrons upon addition of an electron to the LUMO. Here, f −k measures the capacity of atom k to provide a nucleophilic attack. Similarly, f +k denotes the tendency to act as an electrophile.45 Multiple Sequence Alignment (MSA). A majority of the Ca2+ binding proteins possess EF-hand motif,46,47 characterized by the helix−loop−helix structure, where the Ca2+ ion binds to the 12-residue-long loop region. We have considered 100 EFhand loop sequences from 36 distinct subfamilies of EF-hand Ca2+ binding proteins from “The EF-hand Calcium Binding Proteins Data Library” (http://structbio.vanderbilt.edu/cabp_ database) and MSA is performed using Clustal Omega.48 The degree of conservation of a residue ‘x’ at loop position ‘i’ is

reversible covalent linkage of phosphate group from adenosine triphosphate (ATP) to specific protein sites, catalyzed by specific enzymes like kinases. The side chain Oh atom from serine, threonine, or tyrosine acts as a nucleophile attacking the phosphorus (Pγ) of the ATP in the presence of the kinase22 to get reduced to ADP. A recent DFT based Gaussian 03 study shows that the Oh atom can act as a nucleophile for transfer of phosphoryl group from ATP to tyrosine in insulin receptor tyrosine kinase via a transition state with barrier height of ∼14 kcal/mol.27 Phosphorylation of Ca2+ bound CaM leads to differential regulation of CaM-dependent enzymes.28 The residues Y99 and Y138 are potential in vivo phosphorylation sites28 in CaM involving their phenyl ring Oh atom as found in the mass spectrometric analysis. Although serine, threonine, and tyrosine are all present in CaM loops, neither serine nor threonine undergoes phosphorylation. In this study, we carry out implicit solvent QC calculations on all the four Ca2+ bound EF hand loops of CaM29,30 isolated from the crystal structure (PDB ID: 1CLL)30 to extract the ground-state energy levels including the HOMO and the LUMO levels. The isolated Ca2+ binding EF-hand loops from the crystal structure are truncated with neutral capping by acetyl (COCH3) group at the N-terminus and N-methylamide (CH3NH) group at the C-terminus as used in earlier studies.4,31−33 The QC-DFT calculations are performed in polarizable continuum solvent34 using the B3LYP functional and the 6-31G(d,p) basis set, extensively used in various QC studies on biologically relevant metal binding sites including Ca2+ ions.35−38 Our objectives in the present study are manyfold: (1) We can precisely estimate the terminal capping contribution in different parts of the spectra. (2) We aim to extract the robust features of the stability of the coordination geometry across the loops, emphasizing the roles of ligand and nonligand residues. This is in sharp contrast to earlier studies which focus only on ligand atoms.4 (3) We interpret the HOMO and the LUMO levels in terms of chemical reactivity in the system.1,3 We elucidate the functional implications of different loop residues from the HOMO. In particular, we aim to provide the rationale behind different behaviors of the polar loop residues in phosphorylation. We substantiate the functional implications through the calculation of the Fukui index. It is a measure of the chemical reactivity of a system and is defined as the change in electronic density upon removal of an electron from the HOMO or upon addition of an electron to the LUMO, with the molecular geometry being held fixed.39



METHODS System Details and Computation. The 1.7 Å X-ray crystal structure of Ca2+ loaded (holo) CaM (PDB ID: 1CLL)30 is used in our study. We isolate the Ca2+ binding EF-hand loops: loop 1 (D20 to E31), loop 2 (D56 to E67), loop 3 (D93 to E104), and loop 4 (D129 to E140) from the crystal structure retaining the bound water molecules and the Ca2+ ions. We cap the N-terminus by acetyl (COCH3) group B

DOI: 10.1021/acs.jpcb.5b09713 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B calculated as the ratio of number of sequences having ‘x’ at position ‘i’ to the total number of EF-hand loop sequences and expressed in percentage.

The HOMO−LUMO gaps in the N-terminal domain (NTD) loops are 5.01 eV (Ca2+-loop 1) and 5.37 eV (Ca2+loop 2). Similarly, HOMO−LUMO gaps of 5.71 eV (Ca2+-loop 3) and 5.19 eV (Ca2+-loop 4) are observed for the C-terminal domain (CTD) loops. The coordination levels of the Ca2+-loop 1 complex occupy an energy range of −22.82 eV to −19.52 eV. The coordination energy range in implicit solvent is similar to that calculated in vacuum.5 The other three loops also occupy similar coordination energy range: Ca2+-loop 2 (−23.13 eV to −19.34 eV), Ca2+-loop 3 (−22.84 eV to −19.42 eV) and Ca2+loop 4 (−23.26 eV to −19.70 eV). The CSMOs are interspersed by a few noncoordinating energy levels characterized by negligible contribution from the Ca2+ ion. CSMO. The valence orbitals of the Ca2+ ion contributing to CSMOs are spd hybridized through mixing of 3s, 3p, 3d, and 4s levels of Ca2+ ion, while the coordinating oxygen atoms contribute to CSMOs through their valence 2s-2p hybrid orbitals, similar to the earlier vacuum calculation.5 The contribution of the ligand and nonligand atoms from the four different loops across the coordination energy range is illustrated in Figure 2. The carboxylate oxygen, O1 of the coordinating aspartates at the Ist position of all the loops shows maximum contribution (C(α) I,O1) between E ∼ −20.5 eV and −19.5 eV along with secondary peaks in a lower energy range −23 eV to −21 eV (Figure 2a). The contributions (C(α) III,O1, (α) (α) (α) , C ,C ,C ) of ligand oxygen from the other C(α) V,O1 VII,O XII,O1 XII,O2 residues show similar feature in comparable energy ranges



RESULTS First, we illustrate the effects of neutral terminal capping on the ground-state energy spectrum. Then, we provide detailed insight into the CSMO, the HOMO, and the LUMO levels. All the energy levels are reported with respect to the HOMO energy level. The absolute contribution of N-terminal capping is taken as the sum of squares of MO coefficients, C(α) N−cap of the constituent atoms of the COCH3 group in the αth loop and similarly C(α) C−cap for the C-terminal CH3NH group. We also calculate the contribution of the loop, C(α) pro , as the sum of squares of MO (α) coefficients of its constituent atoms. We show C(α) N−cap/Cpro and 2+ (α) (α) CC−cap/Cpro for the different Ca bound loops of CaM in (1) (1) Figure 1. The relative terminal contributions CN−cap /Cpro

Figure 1. Relative contribution of terminal neutral capping: (a) (1) (1) (1) (2) (2) 2+ C(1) N−cap/Cpro and (b) CC−cap/Cpro for Ca -loop 1, (c) CN−cap/Cpro and (2) (2) (3) (3) (3) 2+ (d) CC−cap/Cpro for Ca -loop 2, (e) CN−cap/Cpro and (f) CC−cap/C(3) pro for (4) (4) (4) 2+ Ca2+-loop 3 and (g) C(4) N−cap/Cpro and (h) CC−cap/Cpro for Ca -loop 4.

(1) 2+ (Figure 1a) and C(1) C−cap/Cpro (Figure 1b) in Ca -loop 1 are negligible at the HOMO and LUMO levels. The terminal capping contributions exist at discrete energy levels below the HOMO, albeit negligible compared to contribution from protein atoms. The terminal capping contributions of Ca2+(2) (2) (2) (2) loop 2 complex, CN−cap /Cpro (Figure 1c) and CC−cap /Cpro (Figure 1d) exhibit similar behavior across the energy range. (3) (3) (3) The C(3) N−cap/Cpro (Figure 1e) and CC−cap/Cpro (Figure 1f) of 2+ (4) (4) (4) Ca bound loop 3 and CN−cap/Cpro (Figure 1g) and C(4) C−cap/Cpro 2+ (Figure 1h) of Ca -loop 4 complex indicate negligible capping effect on the ground-state energy spectra.

Figure 2. CSMO energy range for Ca2+-loop 1 (circle), Ca2+-loop 2 (square), Ca2+-loop 3 (diamond), and Ca2+-loop 4 (triangle) corresponding to the QC calculation on X-ray crystal structure geometry. C

DOI: 10.1021/acs.jpcb.5b09713 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

significant contribution from the polar coordinating O1 atom of N60 (loop position V). The side chain coordinating oxygen atoms of E67 (loop position XII), N60 (loop position V), and the backbone ligand O of T62 (loop position VII) contribute to the immediately higher CSMO (Figure 3d) within ∼0.1 eV. Similar to loop 2, the LCMO of Ca2+-loop 3 (Figure 3e) complex at E = −22.84 eV is strongly characterized by the ligand carbonyl O atom of Y99 (loop position VII) and the nonligand hydrophobic I100 (loop position VIII). The ligand O1 atom of N97 (loop position V) also participates in the LCMO. The other energetically low CSMO (Figure 3f) of Ca2+ bound loop 3 at E = −22.75 eV is dominated by the side chain coordinating atoms of E104 (loop position XII) and N97 (loop position V) along with the backbone carbonyl O of Y99 (loop position VII). The energetically lowest CSMOs of Ca2+-loop 4 are similar to that of Ca2+ bound loop 1. The coordinating carbonyl O atom of polar Q135 (loop position VII) and the adjacent noncoordinating hydrophobic V136 (loop position VIII) dominate the LCMO at E = −22.75 eV (Figure 3g), while the coordinating O1 and O2 atoms of E140 (loop position XII) and D133 (loop position V) contribute predominantly to the immediately higher CSMO (Figure 3h). HOMO and LUMO. We show the HOMO and LUMO levels corresponding to NTD loops (Ca2+-loop 1 and Ca2+-loop 2) in Figure 4 and for the CTD loops (Ca2+-loop 3 and Ca2+loop 4) in Figure 5. NTD Loops. The HOMO of Ca2+-loop 1 (Figure 4a) complex is predominantly derived from the side chain

(Figure 2b−f). The coordinating water molecule (Figure 2g) shows diffused contribution (C(α) W,O) at E ∼ −23 eV to −21 eV. The nonligand contribution C(α) VIII,N (Figure 2h) predominates at E ∼ −20 eV and −19.5 eV along with small contribution at energetically lower CSMOs with E ∼ −23 eV. Figure 3 illustrates the two energetically lowest CSMOs of different Ca2+ bound loops of CaM. The lowest coordination

Figure 3. (a) LCMO and (b) the immediately higher CSMO of Ca2+loop 1. (c) LCMO and (d) the immediately higher CSMO of Ca2+loop 2. (e) LCMO and (f) the immediate higher CSMO of Ca2+-loop 3. (g) LCMO and (h) the immediate higher CSMO of Ca2+-loop 4. The acidic (red), polar (magenta), and hydrophobic (green) residues are color-coded.

molecular orbital (LCMO) of Ca2+-loop 1 (Figure 3a) complex at E = −22.82 eV finds predominant contribution from the ligand O atom of polar coordinating T26 (loop position VII) and the adjacent noncoordinating hydrophobic I27 (loop position VIII). The coordinating carboxylate O1 and O2 atoms of acidic E31 also contributes to LCMO, albeit weak compared to the formers. The immediate higher CSMO (Figure 3b) at E = −22.77 eV is characterized by strong participation of the coordinating oxygen atoms of the acidic residues, E31 (loop position XII) and D24 (loop position V). The LCMO of the Ca2+-loop 2 (Figure 3c) complex at E = −23.13 eV is predominantly derived from the coordinating carbonyl O of T62 (loop position VII) and the noncoordinating hydrophobic residue, I63 (loop position VIII) along with

Figure 4. (a) HOMO and (b) LUMO of the Ca2+-loop 1 complex. (c) HOMO and (d) LUMO of the Ca2+-loop 2 complex. (e) Hybridization of O1 atom of D22 at CSMO and HOMO. (f) Hybridization of O1 atom of D24 at CSMO and HOMO. D

DOI: 10.1021/acs.jpcb.5b09713 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

CTD Loops. The HOMO of the Ca2+ bound loop 3 (Figure 5a) shows strong localized electron density on the side chain phenyl ring of Y99. The LUMO (Figure 5b) at E = 5.71 eV is also characterized through strong contribution from the Y99 side chain. The HOMO of Ca2+-loop 4 complex (Figure 5c) shows localized charge density on the side chain carboxylate group of E139. On the other hand, the electronic density at the corresponding LUMO (Figure 5d) is essentially localized on the phenyl ring of Y138 and the side chain carboxy-amide group of N137. The orbital specific contributions of the Y99 ring atoms reflect strong hybridization through the 2p electrons of the ring atoms at both these levels in the Ca2+-loop 3 complex. Similarly, the carboxylate oxygen atoms of E139 and the phenyl ring atoms of Y138 show hybridization with dominating 2p character at the HOMO and the LUMO levels, respectively, of Ca2+ bound loop 4. Fukui Index. The HOMO of Ca2+ bound loop 3 has strong contribution from Y99, where the side chain phenyl ring has strong localized electron density. Such sites being rich in electrons can act as nucleophiles. Since earlier QC calculations indicate phosphorylation of tyrosine through a transition state involving Oh of phenyl ring as a nucleophile, we ascertain the nucleophilic character for all non-hydrogen atoms of Y99 through the computation of atom-condensed Fukui index.49 Figure 5e illustrates the nucleophilicity index, f −k , for the atoms of Y99. The Oh attached to the phenyl ring shows maximum nucleophilicity followed by the Cγ atom of the ring itself. We illustrate the electrophilic index, f +k , of the Y99 atoms in Figure 5f. The phenyl ring carbon atoms Cβ, Cγ, and Cδ2 show significantly large electrophilicity compared to other atoms. That the phenyl ring Oh does not bear electrophilic character, but possess the capacity to act as a nucleophile, is evident from Figure 5e and f. The nucleophilicity and electrophilicity indices are compared for different polar residues like serine, threonine, and tyrosine in SI, Table S1. The table clearly shows nucleophilicity character of Oh in Y99 orders of magnitude larger than corresponding indices in any other polar residues.

Figure 5. (a) HOMO and (b) LUMO of the Ca2+-loop 3 complex. (c) HOMO and (d) LUMO of the Ca 2+-loop 4 complex. (e) Nucleophilicity index of Y99 atoms in the Ca2+-loop 3 complex belonging to unoptimized system. (f) Electrophilicity index of Y99 atoms in the Ca2+-loop 3 complex belonging to unoptimized system.

carboxylate oxygen, O1 and O2 atoms of D22, D24, and the main chain carbonyl oxygen O atom of G23. The LUMO of Ca2+ bound loop 1 (Figure 4b) at E = 5.01 eV finds major contribution through both backbone and side chain atoms of T28 and T29. The HOMO in the Ca2+-loop 2 complex (Figure 4c) is characterized by strong contributions of the side chain carboxylate oxygen atoms of D64. The contribution to the corresponding LUMO (Figure 4d) originates from the strong localized electron density on the side chain phenyl ring of F65. Both D22 and D24 contribute to the low-lying CSMOs as well the HOMO. However, the nature of hybridization is different in these two levels. We illustrate the orbital (2s, 2px, 2py, and 2pz) specific contributions in the global coordinate (1) system of D22-O1 atom (CIII,D22,O1 ) at the CSMO (E = −21.65 eV) and the HOMO (Figure 4e). The CSMO finds predominant contribution through the 2s orbital, while the 2p orbitals show minimal participation in hybridization. At the HOMO the 2p orbitals contribute significantly to hybridization. Similarly, the orbital specific contributions (C(1) V,D24,O1) of the D24-O1 atom at the CSMO (E = −21.5 eV) and the HOMO (Figure 4f) are shown. The CSMO is sp hybrid with predominant s-character, while the HOMO is characterized by a huge contribution of the 2py orbital. The hybridization at LUMO of Ca2+-loop 1 involving the backbone and the side chain oxygen atoms of T28 and T29 also exhibit a predominant 2p character. A similar hybridization pattern is observed for both HOMO and LUMO of the solvated Ca2+-loop 2 complex involving the 2p electrons of the carboxylate oxygen atoms of D64 and the phenyl ring atoms of F65.



DISCUSSION The hierarchy and nature of hybridization of the CSMOs remain largely similar across the CaM loops. However, subtle differences arise owing to the nature of coordinating residue at position V. The Ca2+ bound loops 1 and 4 are similar as they harbor an acidic residue (D24, D133) at position V, while Ca2+ bound loops 2 and 3 share similar features as they harbor a polar residue (N60, N97) at the equivalent position. The details of some low-lying CSMOs highlighting both ligand and nonligand contributions are given in Table S2 (SI). Our electronic structure calculation indicates that the phenyl ring Oh atom of Y99 can indeed perform a nucleophilic attack which is a prerequisite for phosphorylation, consistent with experimental evidence.28 Since the Ca2+ bound loops are truncated from the full CaM, there may be structural perturbations. This leads us to compare the SPE results obtained from the X-ray crystal structure unoptimized geometry with those of QC optimized geometry for Ca2+ bound loop 1 and loop 3 complexes. The root-meansquared displacement (RMSD) of the representative loops between the crystal structure and the QC optimized structure is ≤0.2 Å, so that the overall geometries of the loops remain similar. The coordination distances in Ca2+-loop 1 and Ca2+loop 3 complexes for both optimized and unoptimized E

DOI: 10.1021/acs.jpcb.5b09713 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

few EF-hand loop sequences from different proteins including CaM loop 3 which have (1) coordinating tyrosine at position VII; (2) flanking residues glycine at position VI and isoleucine at position VIII; and (3) other coordinating residues remaining the same. Such proteins are Ca2+/CaM dependent protein kinase (CDPK), Calbindin D28k, Calretinin, Calcineurin, and Hippocalcin (see Figure 6). Our QC results indicate that the EF-hand loop tyrosines from these proteins may pose as a putative phosphorylation site owing to sequence conservation.

calculations exhibit no appreciable differences (see SI, Table S3). The relative contributions of N-terminal acetyl and Cterminal N-methylamide capping for both these loops are as negligible as the unoptimized system (see SI, Figure S1). The HOMO−LUMO gaps for the optimized Ca2+-loop 1 and Ca2+loop 3 complexes are 5.09 and 5.76 eV, respectively, similar to the unoptimized results. The CSMOs of Ca2+-loop 1 and Ca2+loop 3 complexes after optimization occupy an energy range ∼ −23.0 eV to −19.0 eV similar to that of the unoptimized system. The energetically lowest CSMOs (Figure S2), HOMO, and LUMO levels (Figure S3) of these loops from the optimized QC calculation are illustrated in SI. The LCMO of optimized Ca2+-loop 1 complex at E = −22.52 eV finds a major contribution from the coordinating oxygen atoms of E31, D24, and T26. The immediate higher CSMO within 0.04 eV is predominantly derived from the ligand oxygen atom of T26 and the adjacent noncoordinating I27, along with weak contributions from the coordinating carboxylate oxygens of E31. The contributions to these closely spaced CSMOs get flipped between the unoptimized X-ray and the optimized geometries, albeit their qualitative similarity. The LCMO of optimized Ca2+-loop 3 complex at E = −22.79 eV shows strong hybridization involving the coordinating carbonyl O of Y99 and the adjacent noncoordinating I100, while the other energetically low CSMO at E = −22.72 eV is predominantly derived from the coordinating oxygen atoms of E104 and N97 similar to the X-ray structure calculation (see Figure 3). The HOMO of Ca2+ bound loop 1 upon optimization shows strong hybridization involving the side chain carboxylate groups of D20, D22, and D24. Both G23 and G25 also contribute significantly through their backbone atoms. On the other hand, the corresponding LUMO at E = 5.09 eV is derived from T28 and T29. Thus, the HOMO and the LUMO remain similar before and after optimization of the Ca2+-loop 1 complex. The HOMO and the LUMO (E = 5.76 eV) of optimized Ca2+ bound loop 3 show strong localized electron density over the side chain phenyl ring of Y99 similar to the unoptimized calculation. We also calculate the nucleophilic index of Y99 atoms corresponding to the optimized geometry. A similar picture is obtained after optimization (see SI, Figure S4). The Oh atom shows maximum nucleophilic reactivity compared to the other ring atoms as obtained from the previous calculation. Thus, neither the structure nor the functional implications change under fully optimized calculations. Next we discuss the general implications of our results in the context of EF-hand proteins. The degree of conservation (see Methods for calculation) of the residues contributing to the CSMOs in the EF-hand50 family of proteins is estimated from the Multiple Sequence Alignment (MSA) of EF-hand loop sequences as shown in SI Figure S5. Among the coordinating residues, the axial aspartate (position I, 83%) and the bidentate glutamate (position XII, 82%) are strongly conserved across the family. However, the backbone coordinating polar residue (position VII) is moderately conserved (45%). The nonligand residue types are found to be conserved as well. For instance, the VIIIth position of the loop is a hydrophobic residue (isoleucine, valine, leucine, and methionine) in 95% of the cases. Owing to such strong conservation of different loop residues, the robustness in the hierarchy of contributions from different loop positions toward the CSMOs can be generalized to other EF-hand Ca2+ binding proteins, elucidating a conserved mechanism of Ca2+ ion coordination. There are a

Figure 6. Multiple sequence alignment of EF-hand loop sequences with respect to CaM loop 3 having tyrosine at position VII (in bold face), which we predict to serve as putative phosphorylation sites.



CONCLUSION In conclusion, we have computed the electronic energy spectra of the X-ray structure geometries of isolated Ca2+ ion bound loops of CaM through implicit solvent QC calculation based on DFT approach. We show that the contributions of the neutral capping at the terminal residues of the loops are negligible in the ground state. Our calculation indicates general features of coordination geometry of the CaM loops: (i) The coordinating acidic and polar residues play crucial role in stabilizing the Ca2+ ion within the EF-hand loops. (ii) The noncoordinating residues, especially the hydrophobic residue at position VIII show significant contribution to stabilization of coordination geometry. The HOMO of the NTD loops 1 and 2 are derived from the acidic aspartates, whereas the LUMO are characterized by polar residues in Ca2+ ion bound loop 1 and hydrophobic residue in loop 2. The HOMO and the LUMO of the CTD Ca2+-loop 3 complex exhibits strong localized charge density on the phenyl ring of polar Y99. Similar localized charge density is observed on the acidic E139 at the HOMO and on the phenyl ring of Y138 at the LUMO of Ca2+ ion bound loop 4. We further show through the calculation of Fukui index that the HOMO density on the side chain phenyl ring Oh of Y99 can provide a site for nucleophilic attack, which is a prerequisite to its phosphorylation. We show a systematic way to interpret the ground-state energy levels to understand structural and functional aspects of EF-hand50 proteins. This can be extended to large biomolecules in general where truncation of the molecule is required to carry out QC calculations.



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.5b09713. Tables showing nucleophilicity and electrophilicity indices of different polar residues, contribution of ligand and nonligand residues to CSMOs of all the loops in the unoptimized calculation, and the coordination distances of Ca 2+ -loop 1 and Ca 2+ -loop 3 complexes in unoptimized and optimized systems. Figures illustrating the effects of terminal capping in the ground state after QC optimization, representative CSMOs, HOMO and F

DOI: 10.1021/acs.jpcb.5b09713 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B



LUMO levels of optimized geometries of Ca2+-loop 1 and Ca2+-loop 3, along with the MSA of EF-hand loop sequences from different EF-hand subfamilies. (PDF)

(16) Ye, Y.; Lee, H. W.; Yang, W.; Shealy, S.; Yang, J. J. Probing sitespecific calmodulin calcium and lanthanide affinity by grafting. J. Am. Chem. Soc. 2005, 127, 3743−3750. (17) Linse, S.; Helmersson, A.; Forsen, S. Calcium binding to calmodulin and its globular domains. J. Biol. Chem. 1991, 266, 8050− 8054. (18) Martin, S. R.; Andersson Teleman, A.; Bayley, P. M.; Drakenberg, T.; Forsen, S. Kinetics of calcium dissociation from calmodulin and its tryptic fragments. A stopped-flow fluorescence study using Quin 2 reveals a two-domain structure. Eur. J. Biochem. 1985, 151, 543−550. (19) Pedigo, S.; Shea, M. A. Discontinuous equilibrium titrations of cooperative calcium binding to calmodulin monitored by 1-D 1Hnuclear magnetic resonance spectroscopy. Biochemistry 1995, 34, 10676−10689. (20) Pedigo, S.; Shea, M. A. Quantitative endoproteinase GluC footprinting of cooperative Ca2+ binding to calmodulin: proteolytic susceptibility of E31 and E87 indicates interdomain interactions. Biochemistry 1995, 34, 1179−1196. (21) Waltersson, Y.; Linse, S.; Brodin, P.; Grundstroem, T. Mutational effects on the cooperativity of Ca2+ binding in calmodulin. Biochemistry 1993, 32, 7866−7871. (22) Lassila, J. K.; Zalatan, J. G.; Herschlag, D. Biological phosphoryltransfer reactions: understanding mechanism and catalysis. Annu. Rev. Biochem. 2011, 80, 669−702. (23) Canagarajah, B. J.; Khokhlatchev, A.; Cobb, M. H.; Goldsmith, E. J. Activation mechanism of the MAP kinase ERK2 by dual phosphorylation. Cell 1997, 90, 859−869. (24) Hurley, J. H.; Dean, A. M.; Thorsness, P. E.; Koshland, D. E., Jr.; Stroud, R. M. Regulation of isocitrate dehydrogenase by phosphorylation involves no long-range conformational change in the free enzyme. J. Biol. Chem. 1990, 265, 3599−3602. (25) Johnson, L. N. The regulation of protein phosphorylation. Biochem. Soc. Trans. 2009, 37, 627−641. (26) Johnson, L. N.; Barford, D. The effects of phosphorylation on the structure and function of proteins. Annu. Rev. Biophys. Biomol. Struct. 1993, 22, 199−232. (27) Zhou, B.; Wong, C. F. A computational study of the phosphorylation mechanism of the insulin receptor tyrosine kinase. J. Phys. Chem. A 2009, 113, 5144−5150. (28) Corti, C.; Leclerc L’Hostis, E.; Quadroni, M.; Schmid, H.; Durussel, I.; Cox, J.; Dainese Hatt, P.; James, P.; Carafoli, E. Tyrosine phosphorylation modulates the interaction of calmodulin with its target proteins. Eur. J. Biochem. 1999, 262, 790−802. (29) Babu, Y. S.; Bugg, C. E.; Cook, W. J. Structure of calmodulin refined at 2.2 A resolution. J. Mol. Biol. 1988, 204, 191−204. (30) Chattopadhyaya, R.; Meador, W. E.; et al. Calmodulin structure refined at 1.7 A resolution. J. Mol. Biol. 1992, 228, 1177−1192. (31) Wu, E. L.; Mei, Y.; Han, K.; Zhang, J. Z. Quantum and molecular dynamics study for binding of macrocyclic inhibitors to human alpha-thrombin. Biophys. J. 2007, 92, 4244−4253. (32) He, X.; Zhang, J. Z. The generalized molecular fractionation with conjugate caps/molecular mechanics method for direct calculation of protein energy. J. Chem. Phys. 2006, 124, 184703(1−7). (33) Zhang, D. W.; Chen, X. H.; Zhang, J. Z. Molecular caps for full quantum mechanical computation of peptide-water interaction energy. J. Comput. Chem. 2003, 24, 1846−1852. (34) Tomasi, J.; Mennucci, B.; Cammi, R. Quantum mechanical continuum solvation models. Chem. Rev. 2005, 105, 2999−3093. (35) Siegbahn, P. E.; Borowski, T. Modeling enzymatic reactions involving transition metals. Acc. Chem. Res. 2006, 39, 729−738. (36) Pavlov, M.; Siegbahn, P. E. M.; Sandström, M. Hydration of beryllium, magnesium, calcium, and zinc ions using density functional theory. J. Phys. Chem. A 1998, 102, 219−228. (37) Mehandzhiyski, A. Y.; Riccardi, E.; van Erp, T. S.; Koch, H.; Åstrand, P.-O.; Trinh, T. T.; Grimes, B. A. Density functional theory study on the Interactions of metal ions with long chain deprotonated carboxylic acids. J. Phys. Chem. A 2015, 119, 10195−10203.

AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected]. *E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS S.S. thanks UGC, INDIA for fellowship; M.G. and J.C. thank DST for funding, and M.D.R. thanks the Associateship program of the SNBNCBS.



REFERENCES

(1) Fukui, K. Role of frontier orbitals in chemical reactions. Science 1982, 218, 747−754. (2) Blomberg, M. R. A.; Siegbahn, P. E. M. A quantum chemical approach to the study of reaction mechanisms of redox-active metalloenzymes. J. Phys. Chem. B 2001, 105, 9375−9386. (3) Kagawa, H.; Mori, K. Molecular orbital study of the interaction between MgATP and the myosin motor domain: The highest occupied molecular orbitals indicate the reaction site of ATP hydrolysis. J. Phys. Chem. B 1999, 103, 7346−7352. (4) Lepsik, M.; Field, M. J. Binding of calcium and other metal ions to the EF-hand loops of calmodulin studied by quantum chemical calculations and molecular dynamics simulations. J. Phys. Chem. B 2007, 111, 10012−10022. (5) Sikdar, S.; Ghosh, M.; De Raychaudhury, M.; Chakrabarti, J. Quantum chemical studies on the role of residues in calcium ion binding to Calmodulin. Chem. Phys. Lett. 2014, 605, 103−107. (6) Kretsinger, R. H.; Nockolds, C. E. Carp muscle calcium-binding protein. II. Structure determination and general description. J. Biol. Chem. 1973, 248, 3313−3326. (7) Cohen, A. J.; Mori-Sanchez, P.; Yang, W. Insights into current limitations of density functional theory. Science 2008, 321 (5890), 792−794. (8) Lever, G.; Cole, D. J.; et al. Electrostatic considerations affecting the calculated HOMO-LUMO gap in protein molecules. J. Phys.: Condens. Matter 2013, 25, 152101−152107. (9) Mori-Sanchez, P.; Cohen, A. J.; Yang, W. Localization and delocalization errors in density functional theory and implications for band-gap prediction. Phys. Rev. Lett. 2008, 100, 146401−146404. (10) Rudberg, E. Difficulties in applying pure Kohn-Sham density functional theory electronic structure methods to protein molecules. J. Phys.: Condens. Matter 2012, 24, 72202−72208. (11) Beckingham, K. Use of site-directed mutations in the individual Ca2(+)-binding sites of calmodulin to examine Ca2(+)-induced conformational changes. J. Biol. Chem. 1991, 266, 6027−6030. (12) Browne, J. P.; Strom, M.; Martin, S. R.; Bayley, P. M. The role of beta-sheet interactions in domain stability, folding, and target recognition reactions of calmodulin. Biochemistry 1997, 36, 9550− 9561. (13) Drake, S. K.; Lee, K. L.; Falke, J. J. Tuning the equilibrium ion affinity and selectivity of the EF-hand calcium binding motif: substitutions at the gateway position. Biochemistry 1996, 35, 6697− 6705. (14) Evenas, J.; Thulin, E.; Malmendal, A.; Forsen, S.; Carlstrom, G. NMR studies of the E140Q mutant of the carboxy-terminal domain of calmodulin reveal global conformational exchange in the Ca2+saturated state. Biochemistry 1997, 36, 3448−3457. (15) Yang, J. J.; Gawthrop, A.; Ye, Y. Obtaining site-specific calciumbinding affinities of calmodulin. Protein Pept. Lett. 2003, 10, 331−345. G

DOI: 10.1021/acs.jpcb.5b09713 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B (38) Blomberg, M. R.; Siegbahn, P. E. Quantum chemistry applied to the mechanisms of transition metal containing enzymes – cytochrome c oxidase, a particularly challenging case. J. Comput. Chem. 2006, 27, 1373−1384. (39) North, M. A.; Bhattacharyya, S.; Truhlar, D. G. Improved density functional description of the electrochemistry and structureproperty descriptors of substituted flavins. J. Phys. Chem. B 2010, 114, 14907−14915. (40) Brooks, B. R.; Brooks, C. L., 3rd; Mackerell, A. D., Jr.; Nilsson, L.; Petrella, R. J.; Roux, B.; Won, Y.; Archontis, G.; Bartels, C.; Boresch, S.; et al. CHARMM: the biomolecular simulation program. J. Comput. Chem. 2009, 30, 1545−1614. (41) Phillips, J. C.; Braun, R.; et al. Scalable molecular dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781−1802. (42) Improta, R.; Vitagliano, L.; Esposito, L. Peptide bond distortions from planarity: new insights from quantum mechanical calculations and peptide/protein crystal structures. PLoS One 2011, 6, e24533. (43) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A.; Vreven, T.; Kudin, K. N.; Burant, J. C., et al. Gaussian 03, Revision C.02; Gaussian, Inc., Wallingford, CT, 2004. (44) Dudev, T.; Lim, C. Metal binding affinity and selectivity in metalloproteins: insights from computational studies. Annu. Rev. Biophys. 2008, 37, 97−116. (45) Chattaraj, P. K.; Sarkar, U.; Roy, D. R. Electrophilicity index. Chem. Rev. 2006, 106, 2065−2091. (46) Gifford, J. L.; Walsh, M. P.; Vogel, H. J. Structures and metalion-binding properties of the Ca2+-binding helix-loop-helix EF-hand motifs. Biochem. J. 2007, 405, 199−221. (47) Lewit-Bentley, A.; Rety, S. EF-hand calcium-binding proteins. Curr. Opin. Struct. Biol. 2000, 10, 637−643. (48) Sievers, F.; Wilm, A.; Dineen, D.; Gibson, T. J.; Karplus, K.; Li, W.; Lopez, R.; McWilliam, H.; Remmert, M.; Soding, J.; et al. Fast, scalable generation of high-quality protein multiple sequence alignments using Clustal Omega. Mol. Syst. Biol. 2011, 7, 539−544. (49) Kohn, W.; Becke, A. D.; Parr, R. G. Density functional theory of electronic structure. J. Phys. Chem. 1996, 100, 12974−12980. (50) Grabarek, Z. Structural basis for diversity of the EF-hand calcium-binding proteins. J. Mol. Biol. 2006, 359, 509−525.

H

DOI: 10.1021/acs.jpcb.5b09713 J. Phys. Chem. B XXXX, XXX, XXX−XXX