|Home | About | Journals | Submit | Contact Us | Français|
We performed density functional calculations to examine the effects of solvation, hydrogen bonding, backbone conformation, and the side chain on 15N chemical shielding in proteins. We used N-methylacetamide (NMA) and N-formyl-alanyl-X (with X being one of the 19 naturally occurring amino acids excluding proline) as model systems. In addition, calculations were performed for selected fragments from protein GB3. The conducting polarizable continuum model was employed to include the effect of solvent in the density functional calculations. Our calculations for NMA show that the augmentation of the polarizable continuum model with the explicit water molecules in the first solvation shell has a significant influence on isotropic 15N chemical shift but not as much on the chemical shift anisotropy. The difference in the isotropic chemical shift between the standard β-sheet and α-helical conformations ranges from 0.8 ppm to 6.2 ppm depending on the residue type, with the mean of 2.7 ppm. This is in good agreement with the experimental chemical shifts averaged over a database of 36 proteins containing >6100 amino acid residues. The orientation of the 15N chemical shielding tensor as well as its anisotropy and asymmetry are also in the range of values experimentally observed for peptides and proteins.
Chemical shielding reflects electronic environment of nuclei under observation and therefore contains important information about molecular structure and conformational dynamics. Detailed understanding of the sources of various contributions to chemical shielding in proteins is not only critical for our ability to predict chemical shifts and thus facilitate NMR signal assignment, but also potentially important for improvement in structure characterization of proteins (Cornilescu et al. 1999; Lipsitz and Tjandra 2003) and analysis of protein dynamics (Hall and Fushman 2006).
The past several decades of research have led to the development of accurate theoretical methods and computational schemes for chemical shift calculations in peptides and proteins (reviewed in (Shen and Bax 2007)). The computational methods range from those that base on empirical shielding surface (Wishart and Nip 1998) or sequence homology (Wishart et al. 1997; Shen and Bax 2007) to those that use ab initio quantum mechanical (QM) calculations (de Dios et al. 1993; Oldfield 1995; Xu and Case 2001; Xu and Case 2002). The first-principle QM calculations are very accurate, but their use requires significant computational efforts and, therefore, is limited to small peptides. Nevertheless, accurate QM calculations for small peptides are very important since they elucidate the role of various chemical factors which determine and control chemical shielding in proteins. QM chemical shift predictions for 13Cα, 13Cβ, and 13C′ in proteins are typically more accurate than for 15N (Xu and Case 2002). The main reason is that 15N chemical shift is influenced by numerous factors such as backbone torsion angles and ψ, side-chain torsion angle χ1, hydrogen bonds, neighboring residue types, and electrostatic interactions. This complexity indeed makes accurate prediction of 15N shifts very challenging (Le and Oldfield 1996). An earlier quantum chemical calculation (de Dios et al. 1993) of the 15N shielding has been performed using protein fragments and including electrostatic effects from the rest of the protein in addition to intramolecular hydrogen bonding. However, to our knowledge, bulk solvent effects have not been considered in chemical shielding calculations for proteins. It is acknowledged that since solution NMR measurements in proteins are typically conducted in the presence of water, a highly polar solvent, the effects of bulk solvent molecules can be crucial. To study solvent effects on chemical shifts, one can utilize continuum models. A recent review (Tomasi et al. 2005) points out two aspects important for the applications of continuum models to QM studies of chemical shifts. The first aspect is the perturbation effect of the solvent on the electronic wave function of the solute and the geometric distortion of the solute molecule. The second aspect concerns the importance of both short-range and long-range solute-solvent interactions in determining the solvent effect on the nuclear shieldings. The current concept is that short-range interaction can be effectively handled by a number of explicitly treated solvent molecules from the first solvation shell, while the long-range effects can be described effectively by continuum methods. Tomasi et al. also noted that because the characteristic time scale in NMR spectroscopy is milliseconds and longer, whenever explicit solvent molecules are used, it is necessary to correctly account for the statistical picture inherent in the dynamic nature of the solvation shell (picosecond time scale). Recent years have seen continuum methods being applied to study the solvent effect on chemical shifts of small molecules (Mennucci et al. 2001; Mennucci and Martinez 2005; Aidas et al. 2007). For example, Mennucci and Martinez (Mennucci and Martinez 2005) compared continuum-only description, discrete description in terms of solute-solvent clusters, and mixed discrete/continuum description in order to identify and characterize different aspects of solvation.
In this paper we report density-functional theory calculations for N-methylacetamide (NMA) to analyze in detail the specific and bulk effects of the solvent water on 15N chemical shielding. Figure 1a shows the chemical structure of NMA. We then apply the continuum model in both structure optimization and chemical shielding calculations of N-formyl-alanyl-X amides, where X is one of the 19 naturally occurring amino acids excluding proline. For each compound, calculations were carried out for two backbone conformations, corresponding to a standard α-helix (= −58°, ψ= −47°) and a standard β-sheet ( = −139°, ψ = 135°) respectively, with the latter illustrated for N-formyl-alanyl-Ala in Fig. 1b. The results are compared with the available experimental data as well as with the previous approach that does not account for solvent effects. We also compare the effect of polarizable continuum and structure-related hydrogen bonding on 15N chemical shielding for several residues from protein GB3.
All calculations reported in this paper were performed using the GAUSSIAN03 suite of programs (Frisch et al. 2004). We use density-functional theory (DFT) with three-parameter Becke–Lee–Yang–Parr (B3LYP) exchange-correlation functional (Becke 1988; Lee et al. 1988; Becke 1993). The solvent is taken into account approximately, by employing the conducting polarized continuum model (CPCM) (Barone and Cossi 1998; Cossi et al. 2003) where the solute molecule is placed into a cavity surrounded by the solvent considered as a continuum medium with certain dielectric constant. The charge distribution of the solute polarizes the dielectric medium, which generates surface charges around the cavity and hence in turn polarizes the solute. In our calculations, the dielectric constant of water (78.39) was used and the cavity was chosen to be built up by the simple united atom topological model (UA0), in which the van der Waals surface was built by placing a sphere around each solute heavy atom while hydrogen atoms were enclosed in the sphere of the atom to which they are bonded. The number of surface elements for each sphere was 60, and an area of 0.2 Å2 was set for each surface element. Gauge-invariant atomic orbitals (GIAO) were employed to compute NMR properties (Ditchfield 1974; Wolinski et al. 1990) as implemented in Gaussian03. We used the 6-311+G(2d,p) basis set for all our calculations except for the geometry optimization of the dipeptides which was performed with the 6-31+G(d) basis set.
NMA serves as a simple model representing the amide linkage in proteins. It enables us to eliminate the conformational complexity of peptides and to concentrate only on the solvent’s influence on 15N chemical shielding. To this end we performed a series of calculations for gas phase, continuum-only model, and cluster-continuum model. Figure 1c shows NMA with the first water coordination shell represented by three water molecules that make direct hydrogen bonds with NMA. The parameters used for these models are shown in Table 1, and the results are presented in Table 2.
To elucidate the role of various solvent effects we define the following differences of the various chemical shieldings:
Here the subscripts in chemical shieldings consist of two parts separated by a slash. The first part refers to the geometry of the model and the second part indicates if the chemical shift calculation for the given geometry is performed in gas phase (vac) or in polarizable continuum (cont). Three geometries have been employed in this analysis: NMA_vac refers to the NMA molecule structure optimized in vacuo; NMA+3w refers to the optimized structure of NMA and three water molecules from the first coordination shell within the polarized continuum model; NMA* refers to the molecular geometry obtained by extracting NMA atoms from the NMA+3w structure. These Δσ terms partition the solvent effects into several important contributions. The first difference, , shows how the shielding constant calculated in vacuum is changed due to the distortions in the geometry of NMA caused by its aggregation with three water molecules in continuum. represents the effect of solvation of NMA by the hydrogen-bonded water molecules and by more distant water surrounding modeled by a continuum reaction field. accounts for the total shift due to both the geometrical distortions and the solvation. Our results (Fig. 2) show that , being −0.65 ppm as a deshielding effect, is small compared to the deshielding effect of solvation, , which came out to be −19.38 ppm. Thus we conclude that the solvation is the dominating solvent effect in 15N chemical shielding in NMA. To further analyze the roles played by the three hydrogen-bonded water molecules and the more distant water reaction field, we tried to partition the solvation shift, , as , with the superscripts “1” and “2” representing two artificial paths to account for the two sources of contributions. The two paths differ in the order in which the contributions from explicit and continuum waters are taken into account, see Eqs. 3–7. For either path, the two contributions are additive. As calculated, it is not possible to quantitatively separate the contributions this way as and are not equal, and by definition neither are and . Although takes on different values depending on which water molecules (hydrogen-bonded or distant) are considered first, it is qualitatively clear that the effect of the bound water molecules on 15N chemical shielding in NMA is at least as important as, if not more important than the effect of the more distant water molecules. This is understandable since hydrogen bonding is expected to influence the electronic environment of 15N greatly in the case of NMA with saturated hydrogen bonds.
We also applied a similar analysis to 15N chemical shielding anisotropy (CSA) values obtained from NMA calculations, namely:
with being −23.24 ppm, −3.29 ppm, −15.85 ppm, and −10.69 ppm. This clearly shows that the continuum-only model calculation for NMA* yields 15N CSA value that is already very close to that for NMA clustered with three explicit waters in continuum, since a further correction of −3.29 ppm appears minor compared to −23.24 ppm. This suggests that the continuum-only model will produce a smaller relative error in the CSA values than in the isotropic chemical shift.
We notice that in solvent, the NH bond (in the peptide plane) and CO bond tend to stretch while the CN bond becomes shorter (Table 1). This agrees with recent ab initio and DFT calculations (Selvarengan and Kolandaivel 2004). When the NH bond length was varied between its value in vacuum and in NMA+3w structure, with the other structural parameters being optimized in the continuum-only model, the 15N chemical shielding decreased as the bond stretches (Fig. 3a), indicating a deshielding effect of bulk water which tends to correlate well with the NH bond length (Pearson’s |r|≈1.00). A similar linear correlation exists between the 15N chemical shielding and CO bond length in the continuum-only model (Fig. 3b). This is similar to previously found correlations between 13C and 17O chemical shielding and bond length (Oldfield 2002; Aidas et al. 2007). Our continuum-only calculations for NMA showed that elongation of the NH bond by 0.02 Å results in reduction of 15N chemical shielding by 2.33 ppm (~1.8%), and the elongation of the CO bond by 0.03 Å reduces this shielding by 3.92 ppm (~3.1%). The variation of 15N chemical shielding with the NH and CO bond lengths (accompanied by corresponding changes in the optimized geometry) is small compared to other nuclei (Aidas et al. 2007) or 15N in a different chemical environment (Manalo and de Dios 2002).
All dipeptide structures, where appropriate, adopted χ1 angles close to 180° after the geometry optimization, provided that the optimization started with such conformation. This might not represent the global energy minimum of the dipeptide though, since the energy barriers between the rotameric conformations could hinder the side chain’s rotation to its energy minimum during optimization. The calculated 15N shielding tensors are shown in Tables 3 and and44.
We assume that the differences between the DFT-calculated 15N chemical shift and its true value are systematic and depend only on the local electron density around the nitrogen. Therefore these deviations should be the same for all nitrogens in the similar chemical environment. We chose methylamine (CH3NH2) as a reference compound, because it is the simplest molecule that has similar chemical bonding structure for the amide nitrogen, and reliable experimental 15N chemical shielding in methylamine is available. Then to cancel out possible systematic errors, the chemical shift δsample is computed as (Benzi et al. 2004)
where δsample, refer to the N-formyl-alanyl-X amides under study; σ0= 244.6 ppm is the absolute 15N chemical shielding of liquid ammonia at 25°C (Jameson et al. 1981); is the experimental 15N chemical shielding for methylamine, reported to be 249.5 ppm (Cramer 2004) and is the corresponding theoretical chemical shielding computed at B3LYP/6-311+G(2d,p) level of theory.
Figure 4a–b shows the values of 15N chemical shift as a function of the residue type for the α-helical and β-sheet conformations of N-formyl-alanyl-X. The profile of 15N chemical shift as a function of residue type (Fig. 4a,b) agrees well with the statistically averaged experimental 15N chemical shifts in the α-helices and β-sheets in proteins (Wang and Jardetzky 2002). This indicates that a geometry-optimized structure in a water continuum provides a good model for chemical shifts in solution and the dynamics of local waters are effectively averaged out. If we assume that the available protein chemical shift database samples all possible configurations, then a particular residue type in a particular secondary structure, averaged over the database, ought to produce a reasonable mean isotropic chemical shift for that residue in that secondary structure. Since our calculation used standard α-helical and β-sheet conformations, the isotropic chemical shifts we obtained may represent a mean as well. In fact, by comparison, our calculations overestimate the isotropic chemical shift by 2–3 ppm (on average) for the α-helical conformation for most N-formyl-alanyl-X amides (Fig. 4a), but no over- or underestimation is obvious for the β-sheet conformation (Fig. 4b). Considering the standard error associated with the statistical average, this overestimation may be even less significant.
Now we turn our attention to the overall difference between 15N chemical shifts for α-helical and β-sheet conformations. We define the difference in 15N chemical shifts between the standard β-sheet and α-helical conformations as
where X refers to one of the 19 amino acid residues. This difference in our calculation ranges from 0.76 ppm to 6.19 ppm with an average of 2.74 ppm. Now we would like to compare our results with the gas phase calculations of Poon et al. (Poon et al. 2004) who used the same dipeptide model but without structure optimization and solvent. Their gas phase calculations indicated that the difference in 15N chemical shifts between the standard β-sheet and α-helix varies between 13.2 ppm and 24 ppm, with the average value of 15.8 ppm. Our calculation resulted in a greatly reduced chemical shielding difference between the β-sheet and α-helix, (Fig. 4c), which is now in agreement with the results from the statistical analysis of a database containing more than 6100 amino acid residues in proteins, where this difference ranged from 2.66 ppm to 5.80 ppm with an average of 3.74 ppm (Wang and Jardetzky 2002). A close examination shows that the biggest discrepancy between our values and the experimental data is for threonine (about 3.4 ppm, see Fig. 4c). However, it is well within the range of statistical uncertainty in the protein database, and we note that a variation in chemical shift due to side-chain conformation can be up to ~5 ppm according to statistical analysis for amino acids Val, Ile, Thr, Phe, His, Tyr, and Trp (Wang and Jardetzky 2004), see also our calculations below. This discrepancy may be well due to the fact that side-chain configuration averaging is not considered in our calculation. Figure 4c demonstrates that by applying the continuum-only model the magnitude of is reduced by ~13 ppm from that of Poon et al. (Poon et al. 2004). This reduction in the difference in shielding between the two secondary structures is due to the different deshielding effects the bulk water has for the two backbone conformations: the bulk water deshields 15N by ~18 ppm in the α-helical conformation but only by ~5 ppm in the extended conformation of the β-sheet. This can be qualitatively understood based upon the following consideration. We considered the solvated dipeptide as a molecule embedded in a cavity in bulk dielectric with the dielectric constant of liquid water. The polarization of the surrounding dielectric continuum by the electrostatic potential of the dipeptide induces electric charges, which are distributed on the surface of the cavity. The dipeptide in the α-helical conformation is more “globular” and compact than the dipeptide in the extended β-sheet conformation. Therefore, there are more cavity surface charges in proximity to the amide nitrogen in the helical conformation of the dipeptide than in the β-sheet conformation.
As noted above in continuum-only calculation of NMA, neglecting close-contact solvent can result in underestimation by about 9 ppm of the deshielding effect of solvent in the NMA model with saturated hydrogen bonds. Since secondary structure elements in proteins are often hydrogen-bonded, it is important to know the magnitude of possible underestimation for the 15N chemical shielding in dipeptide caused by ignoring these (specific) interactions. In order to assess this effect, we performed 15N chemical shielding calculations for selected residues from protein GB3. In these calculations, a fragment, C-Y-X-N, containing the residue of interest (X) and its preceding residue (Y) was taken directly from the crystal structure of GB3 (PDB code 1IGD)(Derrick and Wigley 1994) and modified by replacing the end atoms to become N-formyl-Y-X-NH2 (see Fig. 1d). In addition, in the case of α-helix, the side chain of the residue preceding X was replaced with CH3. The hydrogen-bonded residues (through the NH group of X and, where applicable, CO group of the preceding residue) were also taken from the GB3 structure (Fig. 1d) and modified to become either CH3-CO-NH-CH(CH3)-CO-NH-CH3 or CH3-CO-NH-CH(CH3)-CO-NH-CH(CH3)-CO-NH-CH3 (as detailed in Table 5). We then performed chemical shift calculations in vacuum and in continuum for the dipeptide alone and for this hydrogen-bonded cluster. The results (Table 5) show that the dipeptide in continuum model yields deshielding up to 15.3 ppm for α-helix and 9.5 ppm for β-sheet. Depending on the hydrogen bonding geometry, the cluster in continuum calculation can further deshield 15N by about 1 to 4 ppm. This suggests that the polarizable continuum model can account for hydrogen bonding in a realistic protein secondary structure and, therefore, might be a reasonable first approximation for computing 15N chemical shielding.
The chemical shielding tensor contains a wealth of potentially useful structural information, which could be lost when the tensor is reduced to isotropic shielding. Knowledge of individual components and orientation of 15N shielding tensor could be important for many NMR applications, including accurate analysis of protein dynamics from 15N relaxation data (e.g. (Hall and Fushman 2006; Fushman and Cowburn 2001)), TROSY-based experiments and cross-correlation effects involving 15N CSA (Fushman and Cowburn 1999), and the use of residual 15N chemical shift anisotropy upon molecular alignment (e.g. (Lipsitz and Tjandra 2003)) as restraints for structure refinement. It is therefore important to understand the dependence of particular components of the 15N chemical shielding tensor on solvent, side chain, and conformation. In our calculations, the 15N CSA values (Tables 3, ,4)4) fall in the range of experimental values reported for proteins (Tjandra et al. 1996; Fushman et al. 1998; Fushman et al. 1999; Kroenke et al. 1999; Cornilescu and Bax 2000; Kurita et al. 2003; Loth et al. 2005; Wylie et al. 2006; Wylie and Rienstra 2008). The difference in CSA values between the α-helix and β-sheet (mean CSAs of −164.3 ppm and −160.0 ppm, respectively) is consistent with the observations for ubiquitin (Cornilescu and Bax 2000) and GB1 (Wylie and Rienstra 2008). It is worth pointing out that in our data this difference arises primarily from σ22, which is systematically higher in α-helix (by 7.8 ppm on average), while the other two components of the 15N shielding tensor (particularly σ11) show a considerably smaller and less systematic difference between the β-sheet and α-helix conformations (see Tables 3, ,4).4). The calculated 15N CSA values also agree with the solid state NMR measurements in short peptides (Hartzell et al. 1987; Oas et al. 1987; Hiyama et al. 1988; Shoji et al. 1989; Mai et al. 1993; Wu et al. 1995). A good agreement with the experimental data is also found for the angle β between the least shielded component (σ11) of the 15N shielding tensor and the NH bond. The values of β obtained here, from 13.5° to 19.8°, are well in the range of the experimental values (12°–24°) obtained by different NMR techniques, both solution and solid-state (Hartzell et al. 1987; Oas et al. 1987; Hiyama et al. 1988; Shoji et al. 1989; Mai et al. 1993; Fushman et al. 1998; Cornilescu and Bax 2000; Kurita et al. 2003; Loth et al. 2005; Hall and Fushman 2006; Vasos et al. 2006). There seems to be a weak correlation between the β angle and secondary structure, with slightly smaller angles for the β-sheet than for the α-helix (mean β angles are 14.8° and 16.5°, respectively). This also agrees with the experimental findings in GB3 (Hall and Fushman 2006) and ubiquitin (Fushman et al. 1998).
Our calculations show a considerable spread in 15N CSA values, from −154.7 ppm to −168.4 ppm, depending on the residue type and the backbone conformation. This range, however, is smaller than the 15N CSA dispersion observed by solution NMR in ubiquitin and GB3 (Fushman et al. 1998; Fushman et al. 1999; Kover and Batta 2001; Hall and Fushman 2006) and by solid-state NMR in GB1(Wylie et al. 2006; Wylie and Rienstra 2008). This likely reflects the fact that these calculations do not take into account the complexity of local electronic environment in proteins, including interactions with neighboring atoms (e.g. hydrogen bonding, charge and ring-current effects), deviations of the backbone and side chain conformations from those considered here, averaging by anisotropic dynamics etc. Note, for example, that while the type of amino acid residue X varied in our calculations, the torsion angle χ1 was close to 180° (where applicable) for both backbone conformations. In order to explore the effect of side chain’s rotameric state on the 15N chemical shielding tensor, we performed a set of calculations for glutamate (X=Glu) in N-formyl-alanyl-X in the β-sheet conformation, in which the angle χ1 was fixed at −180°, −150°, 60°, and −60° in the geometry optimization. The results showed a significant variation in the anisotropy of the shielding tensor (CSA = −161.38, −165.29, −164.94, and −153.98 ppm, respectively), which is bigger than for the isotropic shielding (107.96, 105.91, 114.05, and 111.18 ppm, respectively). These results emphasize the importance of side-chain conformation for the calculation of the 15N shielding tensor, which could be one of the reasons for the wider spread in the experimental 15N CSA values than that calculated here.
Also we notice that in the α-helical conformation, valine, isoleucine, and threonine have particularly low absolute 15N CSA values. This can be due to the presence of branched side chains in these residues that disturbed the structure by tilting the NH bond out of the peptide plane by about 13°. The effect could result from a particular combination of the torsion angles (, ψ, χ1) as it is not observed in the β-sheet conformation.
It is instructive to discuss the orientation of the individual components of the 15N chemical shielding tensor. The expectation from solid state NMR measurements and planar symmetry of the peptide bond (see e.g. (Oas et al. 1987)) is that the least shielded component (σ11) is lying in the peptide plane and tilted by a small angle (see Fig. 5) from the NH-bond, while the intermediate component, σ22, is orthogonal to the peptide plane. Interestingly, while the least shielded component (σ11) of our calculated 15N chemical shielding tensor lies almost in the peptide plane for all residues, independent of the backbone conformation, the orientations of the other two components differ between the two conformations (Tables 3,,4).4). Our calculations show that, in the β-sheet, σ22 is almost orthogonal to the peptide plane (the biggest tilt is ~9° for valine), which automatically places the most shielded component, σ33, close to the peptide plane. The deviations from “ideal” picture are more dramatic for the α-helix. Here the σ22 component is tilted by as much as 32° for glycine and 24° for valine and isoleucine, and the σ33 component is also significantly tilted away from the peptide plane. These results demonstrate the effect of the backbone conformation on the orientation of the 15N shielding tensor.
Another important aspect is the asymmetry of the shielding tensor, as it is often assumed (e.g. in 15N relaxation analysis) that the tensor is axially symmetric, although solid state NMR data on short peptides indicated that deviations from axial symmetry could be substantial (Oas et al. 1987; Hiyama et al. 1988). For the N-formyl-alanyl-X examples considered here the asymmetry of the 15N shielding tensor, defined as
ranges from 0.43 to 0.56 for the β-sheet and from 0.22 to 0.43 for the α-helix. The higher 15N shielding asymmetry in the β-sheet conformation is consistent with the data reported for ubiquitin (Cornilescu and Bax 2000). These differences in the asymmetry of 15N shielding between the two backbone conformations could be related to the differences in the orientation of the tensor. The absolute values of the asymmetry are somewhat higher than the experimentally observed in solution (Cornilescu and Bax 2000; Loth et al. 2005) but comparable to solid-state NMR data (Wylie et al. 2006), which likely reflects motional averaging expected to be more pronounced in proteins in solution.
To examine the effects of solvation, backbone conformation, and the side chain on 15N chemical shielding in proteins, we performed the density-functional theory calculations with the polarizable continuum solvent model for NMA and N-formyl-alanyl-X amides, where X is one of the 19 naturally occurring amino acids excluding proline. The main results of our calculations can be summarized as follows
Supported by NIH grant GM065334 to DF and by American Chemical Society Petroleum Research Fund (44481-G6) to DSK.