Paramagnetic Perturbation of the 19F NMR Chemical Shift in

Jul 16, 2009 - AlbaNova University Center. ... strategy for modeling 19F NMR chemical shifts in the presence of paramagnetic agents in an aqueous solu...
0 downloads 0 Views 244KB Size
10916

J. Phys. Chem. B 2009, 113, 10916–10922

Paramagnetic Perturbation of the 19F NMR Chemical Shift in Fluorinated Cysteine by O2: A Theoretical Study Xin Li,†,‡ Zilvinas Rinkevicius,† Yaoquan Tu,§ He Tian,‡ and Hans Ågren*,† Department of Theoretical Chemistry, Royal Institute of Technology, AlbaNoVa UniVersity Center, S-106 91 Stockholm, Sweden, Laboratory for AdVanced Materials and Institute of Fine Chemicals, East China UniVersity of Science and Technology, Shanghai 200237, People’s Republic of China, and Biophysical Chemistry, School ¨ rebro UniVersity, 701 82 O ¨ rebro, Sweden of Science and Technology, O ReceiVed: March 24, 2009; ReVised Manuscript ReceiVed: June 22, 2009

We present a combined molecular dynamics and density functional theory study of dioxygen-induced perturbation of the 19F NMR chemical shifts in an aqueous solution of fluorinated cysteine under 100 atm of O2 partial pressure. Molecular dynamics simulations are carried out to determine the dominant structures of O2 and the fluorinated cysteine complexes in water, and the collected structural information is exploited in computation of 19F chemical shifts using density functional theory. The obtained results indicate that the density redistribution of the O2 unpaired electrons between the dioxygen and fluorinated cysteine is responsible for the experimentally observed perturbation of the 19F NMR chemical shifts, where the Fermi contact interaction plays the key role. The O2-induced paramagnetic 19F chemical shift, averaged over the simulation trajectory, is comparable with the reported experimental values, proving the availability of the developed strategy for modeling 19F NMR chemical shifts in the presence of paramagnetic agents in an aqueous solution. The applicability of the combined molecular dynamics/density functional theory approach for dioxygen NMR perturbation to all resonating nuclei including 1H, 13C, 15N, and 19F is emphasized, and the ramification of this for investigations of membrane protein structures is discussed. 1. Introduction Membrane-soluble O2 can, as a paramagnetic agent, significantly influence the spin-lattice relaxation times and give rise to perturbations of the chemical shifts in nuclear magnetic resonance (NMR) spectra. As such, it has, in recent years, become a powerful agent in the study of membrane proteins1,2 since it is highly difficult to determine the structure of membrane proteins through either diffraction measurements or highresolution NMR experiments. Molecular oxygen is inhomogeneously distributed in membranes, with an increasing concentration gradient from the hydrophilic exterior to the hydrophobic interior, and the unpaired electrons on dioxygen can couple to the surrounding 1H, 13C, 15N, or 19F nuclei and result in depthdependent effects on both resonance frequencies and relaxation properties. Moreover, the O2 molecule is uncharged, apolar, less perturbing than traditional paramagnetic agents,3 and can be easily added or removed at an appropriate partial pressure. Therefore, membrane-soluble dioxygen serves as an ideal paramagnetic probe to determine the immersion depth and solvent exposure of proteins, or, conversely a targeted residue mutagenesis may be used to determine the local oxygen concentration in membranes. Molecular oxygen has also been utilized to measure the intermolecular accessibility of watersoluble proteins.4 The use of molecular oxygen can be effectively combined with 19F NMR, which offers advantages such as 100% natural abundance, high sensitivity, a wide chemical shift range (over 900 ppm), and a low background of naturally occurring resonances,5 and 19F nuclei can be introduced at a * To whom correspondence should be addressed. E-mail: agren@ TheoChem.kth.se. † AlbaNova University Center. ‡ East China University of Science and Technology. § ¨ Orebro University.

desired residue of proteins through cysteine mutagenesis and subsequent chemical modification, for example, the thioalkylation reaction with 3-bromo-1,1,1-trifluoropropanone. A series of 19F NMR studies has been reported by Prosser et al. by using dioxygen as a paramagnetic probe to determine the immersion depth of fluorine nuclei in membrane-embedded CF3(CF2)5C2H4-O-maltose,6 to reveal the secondary structure of transmembrane proteins,7 and to get insight into the solvent exposure and peptide binding behavior in proteins.8 These 19F NMR experiments have indicated a significant usefulness of paramagnetic O2 in the study of structures of membrane proteins. For instance, it has been confirmed that the mechanism of relaxation rate enhancement is dipolar rather than of Fermi contact origin, while the chemical shift perturbation mainly arises from a Fermi contact mechanism.2,9 Although the application of dioxygen as a paramagnetic agent is not limited to 19 F NMR,10 it is the sensitivity of the fluorine chemical shift that makes 19F NMR one of the most suitable techniques for the determination of membrane protein topologies. Unfortunately, tools for first-principles calculations and analysis of chemical shifts observed in paramagnetic NMR experiments are still in their infancy of development and can, at the moment, provide only limited support to the interpretations of experimental results.11-16 Unlike the theory of the NMR spin Hamiltonian parameters in closed-shell molecules,17 the firstprinciples theory and the corresponding computational approaches, in the case of open shells, still remain largely unexplored due to the difficulties arising from the necessity to treat the electronic and nuclear Zeeman effects simultaneously.16 In recent years, many efforts have been made in this field to develop a computationally tractable yet accurate methodology for predictions of nuclear shielding constants or chemical shifts in open-shell molecules. Rinkevicius et al.11 reported the first

10.1021/jp902659s CCC: $40.75  2009 American Chemical Society Published on Web 07/16/2009

Paramagnetic Perturbation of the 19F NMR Chemical Shift fully first-principles calculation of nuclear shielding tensors in paramagnetic molecules and demonstrated the accuracy and applicability of the density functional theory (DFT) methods in paramagnetic NMR calculations. This methodology was used later on to study the influence of hydrogen bonds on nuclear shielding constants of 13C and 1H nuclei in paramagnetic nitronylnitroxide radicals.12 Recently, Moon and Patchkovskii13 extended the theory and derived a general formula for doublet open-shell systems with no thermally accessible electronically excited states. Following this work, Pennanen and Vaara14 generalized the nonrelativistic formulation and incorporated the relativistic spin-orbit effect on both hyperfine and Zeeman interactions consistent up to O (R4) order, where R is the fine structure constant. Later, Hroba´rik et al.15 extended the theory to spin states higher than the doublet by taking into account the zero field splitting (ZFS) for cylindrically symmetric metallocenes. Recently, Pennanen and Vaara16 presented a general theory of nuclear shielding for molecules with a ground state of arbitrary spin multiplicity. So far, to the best of our knowledge, no first-principles investigations of paramagnetic perturbation of 19F chemical shifts induced by O2 have been reported. In this work, we take the first step toward developing a methodology for determining the perturbation of chemical shifts by paramagnetic agents by simulating the paramagnetic perturbation of 19F chemical shifts in an aqueous solution of fluorinated cysteine enriched by dioxygen and applying a 100 atm O2 partial pressure. The selection of fluorinated cysteine as a model compound is motivated by the cysteine mutagenesis and thioalkylation method used by Prosser et al.7 in the experiment of measuring the structure of a membrane protein. All simulations are performed using a combined molecular dynamics and density functional theory approach, which previously has been successfully used to investigate chemical shifts of closed-shell molecules in solution environments.18 2. Simulations Protocol and Computational Details The computational determination of the paramagnetic perturbation of the 19F chemical shift in an aqueous solution of fluorinated cysteine induced by O2 can be split into three distinct steps, (a) molecular dynamics simulations of fluorinated cysteine solution under selected conditions (partial O2 pressure, temperature, etc.) in order to determine the time scale of dioxygen residing near the -CF3 group of the fluorinated cysteine and the structure of the “O2 + fluorinated cysteine” complexes; (b) density functional theory calculations of the nuclear shielding constants of fluorine atoms in the paramagnetic O2 + fluorinated cysteine complex and stand-alone diamagnetic fluorinated cysteine using the structural information from a previous step as the input data; and (c) determination of the paramagnetic effect of O2 by computing the time-averaged change of the 19F chemical shift perturbation caused by the addition of dioxygen into an aqueous solution of fluorinated cysteine using the information gathered in the above steps. In the following subsections, we describe the computational techniques and approximations employed in each step of the outlined computational scheme. 2.1. Molecular Dynamics Simulations: Fluorinated Cysteine in Aqueous Solution. The molecular dynamics simulations have been carried out for the zwitterionic form of fluorinated cysteine (see Figure 1), which is the lowest-energy structure of fluorinated cysteine in aqueous solution. The simulation box is chosen to contain one fluorinated cysteine and two O2 molecules as the solute and 850 water molecules

J. Phys. Chem. B, Vol. 113, No. 31, 2009 10917

Figure 1. Chemical structure of the zwitterionic form of fluorinated cysteine.

as the solvent. Here, the concentration of molecular oxygen in the simulation box is given by Henry’s law at the oxygen partial pressure of 100 atm.19 An isothermal-isobaric (NPT) ensemble is used in the simulations with P ) 100 atm and T ) 298 K, maintained by using the Parrinello-Rahman barostat20 and the Nose´-Hoover thermostat,21 respectively. All of the simulations are carried out with a time step of 1 fs. The cutoff radius for the Coulomb and van der Waals short-range interactions is set to 13 Å. During the simulations, periodical boundary conditions are applied, and the long-range interactions between the atoms are calculated by the particle mesh Ewald (PME) summation.22 The system is subject first to an energy minimization and subsequently to a 1 ns simulation to reach equilibrium. After that, another 2 ns simulation is performed for data acquisition, with snapshots saved every 0.1 ps. The molecular dynamics simulations were carried out using the Gromacs program package version 3.3.323,24 with the force field tailored for systems investigated in this work. More specifically, for fluorinated cysteine, the Lennard-Jones parameters for modeling nonbonded van de Waals interactions were taken from the general Amber force field (GAFF) database,25 while the partial atomic charges used for modeling the intraand intermolecular nonbonded Coulomb interactions were derived from the restrained electrostatic potential (RESP) procedure.26 The geometry of fluorinated cysteine was optimized using the hybrid B3LYP exchange-correlation functional27 and the 6-31++G(d,p) basis set.28 The water environment was modeled by means of the polarizable continuum model (PCM).29 The partial charges for the equilibrium geometry were obtained by fitting the electrostatic potential from a single-point B3LYP/cc-pVTZ27,30 calculation, in which the solvation effects were taken into account by PCM, and thus, the RESP charges could be expected to appropriately model the interactions between the atomic partial charges in aqueous solutions. All of the above enumerated DFT calculations were carried out using the Gaussian03 program package.31 During the molecular dynamics simulations, the internal coordinates of fluorinated cysteine were fixed, except that the -CF3 group was allowed to rotate freely since we were only interested in the interactions between fluorine atoms and O2 molecules. For O2, we employed the force field model proposed by Zasetsky et al.;32 this model consists of two Lennard-Jones centers and three point charges and is designed to reproduce the experimental quadrupole moment of the oxygen molecule in the gas phase. Finally, for water molecules, the well-known TIP4P model was used.33 The selected force field parameters are expected to give reasonable structural parameters for the van der Waals bonded O2 + fluorinated cysteine complex as well as the residing time scale of O2 in the adjacent region of the fluorine atoms in cysteine. Since the perturbation of the 19F chemical shift caused by the paramagnetic agent O2 is expected to be highly dependent on the distance between O2 and the -CF3 group, the residence time of O2 in the vicinity of the fluorine atoms will directly affect the magnitude of the observed chemical shift perturbation. The cutoff radius, which determines the vicinity of O2 to the -CF3 group, is set to 4.0 Å, that is, O2 is considered to be close to the fluorinated part of cysteine if at least one of the O-F distances is shorter than 4.0 Å. The selection of this cutoff

10918

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

Figure 2. Evolution of the shortest O-F distance between O2 and the -CF3 group with respect to simulation time. The dashed line shows the cutoff radius of 4.0 Å.

radius is governed by the behavior of the Fermi contact contribution to the hyperfine coupling tensor of the fluorine atoms, which becomes negligible when the O-F distance is larger than 4.0 Å. Here, we follow the methodology by Impey et al.34 to determine the mean residence time τO2 of O2 in the adjacent region of the fluorine atoms by simply fitting to a function N(t), which decays exponentially with respect to the τO2, that is, N(t) ∼ exp(-t/τO2). The function N (t) is defined as N(t) ) 〈∑j Pj(t0,t0 + t;t*)〉t0 and can be computed directly from the trajectory obtained in molecular dynamics simulations. In the expression of N(t), the summation runs over the two O2 molecules, and the angle brackets denote averaging over all t0 times. Pj(t0,t0 + t;t*) takes the value of 1 if O2 molecule j is in the adjacent region of the -CF3 group at times t0 and t0 + t and does not leave the region for any continuous period longer than t*. In all other cases, Pj(t0,t0 + t;t*) takes the value of 0. To be accurate, the parameter t* is set to 0 ps in our calculation, and τO2 is determined to be 1.75 ps. Such a small value of τO2 indicates that the O2 molecules enter and leave the vicinity of the -CF3 group very frequently (see Figure 2), in accordance with the fact that the O2 molecules and the -CF3 group interact very weakly. To accurately determine the influence of dioxygen residence time on the resulting paramagnetic chemical shift perturbation, we checked the O-F distance between O2 and the -CF3 group in all of the 20000 snapshots saved during the 2 ns simulation and found that there are 1511 snapshots in which at least one O2 molecule is close to the -CF3 group (shortest RO-F < 4.0 Å). Neglecting the effect of distant O2 on the 19F chemical shift, it is straightforward to calculate the average chemical shift perturbation in the 1511 snapshots and then scale it by factor of 1511/20000, that is, 0.076, to obtain the final average over the whole simulation trajectory. However, DFT calculations of 1511 snapshots are too time-consuming, and in the reported combined MD/DFT calculations of NMR parameters, several tens of snapshots can give reasonable results (see Bu¨hl’s papers in ref 18). Therefore, to reduce computational costs, we selected 75 snapshots from every 20 of the 1511 ones, and calculated the average of these 75 snapshots instead, which was then scaled by 0.076 to obtain the final average value of the chemical shift perturbation. 2.2. Density Functional Theory Calculations of Nuclear Shielding Constants. The second step in our proposed scheme for determining the paramagnetic perturbation of 19F chemical shifts is the evaluation of the nuclear shielding tensors of fluorine atoms in the -CF3 group of the isolated diamagnetic fluorinated cysteine and of the paramagnetic O2 + fluorinated cysteine

Li et al. complex. In the case of diamagnetic molecules, here, the fluorinated cysteine, the calculation of the 19F nuclear shielding constant is a straightforward application of Ramsey theory, that is, only the orbital contributions to the shielding tensor must ftot be evaluated and the total shielding tensor f σ becomes equal to the sum of the so-called paramagnetic and diamagnetic forb contributions to the orbital shielding tensor f σ . Unfortunately, the Ramsey theory is not directly applicable to paramagnetic systems, like the O2 + fluorinated cysteine complex, and in addition to the above-mentioned orbital nuclear shielding fT f tensor, a new temperature-dependent contribution f σ to f σ tot must also be considered. Under these circumstances, following Pennanen and Vaara,16 the Cartesian τ components of the total shielding tensor in a paramagnetic system can be written as

tot orb T orb στ ) στ + στ ) στ -

µB γkT

∑ gaAbτ〈SaSb〉0

(1)

ab

where µB is the Bohr magneton, γ is the gyromagnetic ratio of the magnetic nucleus, k is the Boltzmann constant, T is the temperature, ga and Abτ are the Cartesian components of the electronic g-tensor and hyperfine coupling tensor, respectively, and finally, 〈SaSb〉0 is the Cartesian component of the Boltzmannaveraged 3 × 3 matrix of the products of the effective spin b S components, computed in the limit of absence of the external T can magnetic field. As shown by Pennanen and Vaara,16 the στ be further decomposed into nine distinct contributions (see Table 1 in ref 16) by taking the electronic g-tensor corrected up to O (R2) order and the hyperfine coupling tensor corrected up to O (R4) order. Therefore, the evaluation of the total nuclear ftot shielding tensor f σ in a paramagnetic system becomes computationally very complex; it involves the computation of the orbital shielding tensor, which solely defines the nuclear shielding tensor in closed-shell molecules, the electronic gtensor, the hyperfine coupling tensor, and, in the case of a system with more than one unpaired electron, the zero-field splitting tensor. However, in the case of an organic system, like the O2 + fluorinated cysteine complex, in which the spin-orbit contribution to the hyperfine coupling tensor is negligible and the electronic g-tensor only slightly deviates from the free fT electron g-factor ge, the f σ expression can be significantly simplified by taking into account only the isotropic Fermi contact contribution (see the first term in Table 1 in ref 16), that is

T στ ≈ -δτ

1 µBge S(S + 1)Acon γ 3kT

(2)

where Acon is the Fermi contact part of the isotropic hyperfine fT coupling constant. In fact, a set of sample calculations of f σ of 19F in the O2 + fluorinated cysteine complex using the complete expression given by eq 1 (see Table 1 in ref 16 for fT the list of considered contributions to f σ ) deviates from the result obtained using the approximate expression given by eq 2 only at most by 0.5%. In view of these findings, we will use exclusively eq 2 for the evaluation of the temperature-dependent part of the nuclear shielding tensor as the errors introduced by fT the simplified treatment of f σ in general are expected to be

Paramagnetic Perturbation of the 19F NMR Chemical Shift

J. Phys. Chem. B, Vol. 113, No. 31, 2009 10919

minor compared to the accuracy with which we compute the Fermi contact part of the hyperfine coupling constants. In this work, the nuclear shielding of fluorinated cysteine and its van der Waals complex with dioxygen have been determined forb using density functional theory. The orbital shielding tensor f σ has been calculated by the unrestricted DFT gauge-independent atomic orbital (GIAO) methodology,35 as implemented in the Gaussian03 program,31 and the corresponding orbital shielding orb constant σiso was obtained as the average of diagonal elements f orb of f σ . The isotropic part of the hyperfine coupling tensor of the fluorine atoms in the O2 + fluorinated cysteine complex has been determined using the density functional restricted-unrestricted approach36 based on the spin-restricted DFT formalism, which is implemented in the DALTON program.37 All DFT calculations were performed using the hybrid B3LYP exchange-correlation functional27 in combination with the ordinary Huzinaga basis set,38 which is frequently used for EPR and NMR property calculations. The B3LYP functional has proven to be able to describe the spin density properly and to give comparable results to experimental data in our previous studies of paramagnetic NMR chemical shifts.11,12 The Huz-II basis set (which is also known as “IGLOII”) has been extended to fluorine atoms by decontraction of s shells and addition of two tight s functions in order to facilitate accurate evaluation of the Fermi contact contribution to the hyperfine coupling constants of the fluorine atoms. In nuclear shielding calculations, we considered the minimal model of the solvated O2 + fluorinated cysteine complex and included only dioxygen and fluorinated cysteine molecules with geometrical arrangement obtained from the molecular dynamics snapshots. Since it is found in simulations that at most only one O2 molecule is close to the -CF3 group at a time, the input for DFT calculation contains only one O2 and one fluorinated cysteine. The choice of such a simple model leads to a neglect of the solvation effects on the hyperfine coupling constants and orbital nuclear shielding constants in each considered snapshot. However, since the hyperfine coupling constants of the fluorine atoms are mostly dependent on the dioxygen orientation and the distance with respect to the -CF3 group of fluorinated cysteine, this model is expected to describe the behavior of the nuclear shielding constant of 19F in aqueous solution at least on a semiquantitative level. Furthermore, use of this minimal model significantly reduces computational cost in the determination of nuclear shieldings of the fluorine atoms for each snapshot and consequently allows one to perform averaging of chemical shifts over extended time intervals. 2.3. Determination of Paramagnetic Perturbation of the Chemical Shift. The chemical shift of a paramagnetic molecule, like the O2 + fluorinated cysteine complex, is typically decomposed into three distinct contributions, namely, orbital, contact and pseudocontact contributions13

δtot ) δorb + δcon + δpc

(3)

where the orbital chemical shift δorb is, by definition, equal to the difference of orbital nuclear shielding constants between the reference compound and the investigated compound and the contact δcon and pseudocontact δpc shifts are equal to the corresponding shielding constants with opposite signs. The last two contributions to the chemical shift arise from the temperature-dependent part of the nuclear shielding tensor and are nonvanishing only in paramagnetic molecules. Taking this into account, the paramagnetic perturbation of the chemical shift in a diamagnetic molecule by a paramagnetic contrast agent can be expressed as

orb pc ∆δtot ) -∆σiso - σcon - σiso

(4)

where ∆σ orb iso is the change of the isotropic nuclear orbital shielding constant of the diamagnetic molecule, caused by its complexation pc are the with the paramagnetic contrast agent, and σ con and σ iso contact and pseudocontact shifts of the “diamagnetic molecule + paramagnetic contrast agent” complex, which arises due to the paramagnetic nature of the contrast agent. Therefore, determination of the chemical shift perturbation by a contrast agent requires the evaluation of the orbital nuclear shielding for the diamagnetic molecule and the full nuclear shielding tensor (see eq 1) for its complex with the contrast agent. In this work, an accurate estimation of the contact part is given by eq 2, and the pseudocontact part, which is expected to be averaged to 0 in an isotropic environment due to fast exchange and isotropic tumbling of the paramagnet,2 is disregarded. The obtained result uniquely defines the perturbation of the chemical shift of the investigated molecule by the given contrast agent and is independent of the choice of the reference compound used in the determination of the chemical shift. With the application of the outlined computational procedure, the determination of the perturbation of the 19F chemical shift of fluorinated cysteine in aqueous solution by dioxygen can be split into two separate tasks, (a) the evaluation of the averaged orbital nuclear shielding constants in fluorinated cysteine in aqueous solution under ambient pressure and (b) the evaluation of the averaged orbital and contact nuclear shielding tensors of the O2 + fluorinated cysteine complex solvated in water under the partial pressure of O2 increased to 100 atm (here, as discussed in previous paragraphs, we disregard the pseudocontact part of the nuclear shielding constant and compute only the contact part according to eq 2). In this way, one obtains the averaged value of the paramagnetic perturbation 〈∆δtot〉 characteristic of fluorinated cysteine under selected conditions, provided that the time scale of the average is sufficiently large to sample accurately the conformational space of the O2 + fluorinated cysteine and is significantly larger than the O2 residence time near the -CF3 group of fluorinated cysteine. Here, we would like to point out that the primary factor determining the quality of the obtained 〈∆δtot〉 is that of the trajectories obtained in molecular dynamics simulations, which ultimately defines the contact shift (via strong fluorine atom hyperfine coupling constants) dependence on distance between O2 and the -CF3 group of the fluorinated cysteine. 3. Results 3.1. Dynamical Motion of O2 from MD Simulations. Figure 2 shows the closest O-F approach between O2 and the -CF3 group in each snapshot, which is one of the decisive factors that affect the magnitude of the paramagnetic perturbation of the 19F chemical shift induced by dioxygen. At each time, the shortest O-F distance is compared with the cutoff radius of 4.0 Å to characterize the residence of the O2 molecule, that is, inside or outside of the adjacent region of the fluorine nuclei. It can be seen that O2 approaches and leaves the fluorine nuclei for many cycles during the simulation, and the mean residence time of O2 in the adjacent region of -CF3 is calculated to be 1.75 ps using a fitting procedure described in the previous section. Therefore, the simulation time of 2 ns is considered to be long enough to give a rational time average of the nuclear shielding constants of fluorine atoms in fluorinated cysteine. Moreover, in a combined MD and NMR investigation of the oxygen-induced paramagnetic shifts reported by Al-Abdul-

10920

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

Figure 3. Orbital shielding constant of 19F NMR (solid circles) in the absence of O2 and (squares) in the presence of O2.

Wahid et al.,39 the diffusional behavior of molecular oxygen is found to play an insignificant role in affecting the 13C NMR paramagnetic shift, and it is the local concentration of O2 that determines the depth-dependent 13C shift in membranes. Since the concentration of O2 in our simulations is given by Henry’s Law, which is applicable under a 100 atm of oxygen partial pressure, the MD simulations are expected to give reliable configurations for the fluorinated cysteine-O2 system. 3.2. Orbital and Contact 19F Shifts. The orbital nuclear shielding constants of 19F in fluorinated cysteine are calculated in the absence and presence of O2, corresponding to the experimental conditions of ambient pressure and a 100 atm of oxygen partial pressure, respectively. From these two calculations, the dioxygen-induced orbital shift can be computed as a difference between the isotropic orbital nuclear shielding orb . The time evolutions of the constants, that is, ∆δorb ) -∆σiso orbital shielding constants are shown in Figure 3. The existence of molecular oxygen in solution gives rise to a slight decrease of the shielding constant for each snapshot, and the orbital shift is a little bit enlarged accordingly. Over the 75 snapshots, the average perturbation of the orbital shift induced by O2 amounts to ∆δorb ) 0.84 ppm. Here, the long bar over ∆δorb denotes the average over the 75 snapshots, and after multiplication by a factor of 0.076, it will become the average over the whole simulation trajectory, 〈∆δorb〉, as described in section 2.1. The contact 19F shift in this work is computed as δcon ) -σcon, where σcon is determined from the isotropic hyperfine coupling constant of the fluorine atom using eq 2. Over the 75 snapshots, the average contact shift is δcon ) 43.64 ppm. Similarly, 〈δcon〉 ) 0.076 × δcon. It is undoubted that the contact shift gives the major contribution to the perturbation of the 19F chemical shift induced by the O2 molecule. The contact shift is known to be related to the unpaired spin density that arises from the spin delocalization of the unpaired electrons on dioxygen and spin polarization. Therefore, it is of importance to examine the spin density distribution in the fluorinated cysteine-O2 system. As a typical example, the spin density contour of the forth snapshot is shown in Figure 4, where positive spin density is shown in blue while the negative one is in green. Direct spin delocalization of dioxygen gives rise to positive spin density at fluorine nuclei and downfield 19F shifts, that is, positive contact shift, while spin polarization gives rise to negative spin density at fluorine nuclei and consequently upfield 19F shifts. On the basis of the shape of the spin density distribution shown in Figure 5, we conclude that the absolute value of the contact shift depends on the shortest O-F distance between the dioxygen and fluorine nuclei, while the sign is related to the relative orientation between them, that is, the F-O-O angle

Li et al. corresponding to the closet O-F approach (illustrated in the inset of Figure 5a). Although the contact shift in each snapshot mainly depends on two factors (the shortest O-F distance and the corresponding F-O-O angle), we can elucidate a single dependence in a statistical way. As shown in Figure 5a and b, the contact shifts are averaged over every 0.2 Å in the shortest O-F distance and over every 10° in the corresponding F-O-O angle, respectively. It can be seen that when the shortest O-F distance becomes close to 4.0 Å, δcon decays very fast. Therefore, it is reasonable to set 4.0 Å as the cutoff radius and neglect the paramagnetic effect of the remote O2 molecules. Moreover, when the F-O-O angle is around 90 or 180°, the fluorine nuclei are significantly affected by the negative spin densities, and the contact shifts become negative, corresponding to the upfield 19F shifts. This can be easily interpreted by the shape of the spin density distribution in Figure 4, in which the negative spin density is mainly located along the principal axis and in the plane bisecting the axis. On the basis of the shape of distributions shown in Figure 5a and b, we use a simple function f(r,φ) consisting of a quadratic polynomial of φ and an r-9 term to fit the contact shift

f(r, φ) )

c1φ2 + c2φ + c3 r9

(5)

where r denotes the shortest O-F distance in each snapshot and φ denotes the corresponding F-O-O angle, as illustrated in the inset of Figure 5a. Note that r is in angstroms and φ is in degrees. The coefficients c1, c2, and c3 are obtained as -1.8511 × 103, 4.6429 × 105, and -2.6912 × 107, respectively, through a least-squares fitting procedure. Figure 6 shows the correlation between the contact shifts obtained from DFT calculations and the fitted function f(r,φ). Although f(r,φ) is very simple, it demonstrates two facts at least qualitatively, (a) the sign of the O2-induced chemical shift perturbation is determined by the distribution of spin density around O2 molecule, and (b) the magnitude of the chemical shift perturbation decays much faster than the relaxation rate enhancement (which is of the order 1/r6) with respect to O-F distance between paramagnetic O2 and resonating fluorine atoms. 3.3. Time-Average of the O2-Induced Perturbation of the 19 F Chemical Shift. The O2-induced perturbation of the chemical shift is obtained as the sum of orbital and contact contributions according to eq 4, ∆δtot ) ∆δorb + δcon, where we disregarded pseudocontact shift. As shown in Figure 7, the magnitude of ∆δtot fluctuates with respect to the number of snapshots, while the running average, that is, the average value up to a given snapshot, converges well within (1 ppm during the last 10 snapshots, indicating that the average over the 75 snapshots is sufficient to give the time average value of the oxygen-induced paramagnetic shift. The error estimation of the average of ∆δtot over the 75 snapshots is obtained to be 7.86 ppm by using the blocking average method of Flyvbjerg and Petersen.40 Therefore ∆δtot ) 44.48 ( 7.86 ppm. As described in section 2.1, the time average of the O2-induced paramagnetic shift over the whole simulation trajectory is obtained by a further multiplication by 0.076

〈∆δtot〉 ) 0.076 × ∆δtot ) 3.38 ppm ( 0.60 ppm

(6)

Paramagnetic Perturbation of the 19F NMR Chemical Shift

J. Phys. Chem. B, Vol. 113, No. 31, 2009 10921

Figure 4. Spin density contour plotted at isovalue ) 2 × 10-5 for the forth snapshot. Blue indicates positive spin density, and green represents negative values.

Figure 6. Correlation between the contact shifts obtained from DFT calculations and the fitted function f(r,φ), where r denotes the shortest O-F distance in each snapshot and φ denotes the corresponding F-O-O angle.

Figure 5. Evolution trend of the contact shift with respect to the shortest O-F distance between O2 and the -CF3 group and the corresponding F-O-O angle; (a) average contact shift over every 0.2 Å in the O-F distance; (b) average contact shift over every 10° in the F-O-O angle. Digits on the top of the columns show the number of snapshots used in averaging.

The positive value of 〈∆δtot〉 indicates a downfield shift of the 19F NMR signal. This result is comparable to those observed for dioxygen-induced 19F shift perturbations under a 100 atm of oxygen partial pressure, as reported by Prosser et al., for the lipid-bilayer-embedded CF3(CF2)5C2H4-O-maltose (3.5-7 ppm)6 and the membrane-embedded diacylglycerol kinase (0-3 ppm).7 The agreement with experimental data proves the validity of such a combined MD/DFT approach in calculating paramagnetic NMR parameters, and the accuracy and reliability can be further improved in future studies by taking into account solvent effects and membrane/protein structures, as well as by using more snapshots in averaging. 4. Conclusion In this paper, the O2-induced perturbation of the chemical shift in 19F NMR measurements is studied for the first time

Figure 7. Time evolution of the total paramagnetic shift of the fluorinated cysteine-O2 system in the selected 75 snapshots; (squares) paramagnetic shift in each snapshot; (solid circles) average value up to this point. The dashed line shows the convergence to the final average.

through a combined molecular dynamics and first-principles approach, using an aqueous solution of fluorinated cysteine as a model system. Classical molecular dynamics simulations are carried out to produce instantaneous configurations of the fluorinated cysteine-O2 system, and the subsequent density functional calculations are performed with the hybrid B3LYP functional and the Huzinaga basis set. The calculated Fermi contact shift is in each snapshot found to depend not only on the shortest O-F distance between the O2 and fluorine atoms but also on the corresponding F-O-O angle. Such dependencies are also confirmed by the spin density distribution in the

10922

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

cysteine-O2 system, where positive spin density at fluorine nuclei gives rise to downfield 19F shifts and negative spin density gives rise to upfield 19F shifts. Moreover, it is noteworthy that the contact shift decays sharply with respect to the shortest O-F distance, making it possible to neglect the paramagnetic effect of the distant dioxygen. With the pseudocontact shift disregarded for fluorine nuclei in the cysteine-O2 system, the total perturbation of the chemical shift of 19F by dioxygen is equal to the sum of orbital and contact contributions. A downfield shift of 3.38 ( 0.60 ppm is obtained as a time average of the oxygen-induced paramagnetic shifts over the whole simulation time, which is comparable with the experimental data. It can be noted that a combined MD and DFT approach, as the one presented here, is applicable to all resonating nuclei, including 1H, 13C, 15N, and 19F. The accuracy and usefulness of further theoretical investigations can be improved by taking into account the bulk environment such as solvent effects and membrane/protein structures. The scope for investigations and use of our technology is certainly very wide, something we will capitalize on in the near future. Acknowledgment. This work was supported by a grant from the Swedish Infrastructure Committee (SNIC) for the project “Multiphysics Modeling of Molecular Materials”, SNIC 023/ 07-18. References and Notes (1) Prosser, R. S.; Evanics, F.; Kitevski, J. L.; Patel, S. Biochim. Biophys. Acta 2007, 1768, 3044. (2) Bezsonova, I.; Forman-Kay, J.; Prosser, R. S. Concepts Magn. Reson. Part A 2008, 32A, 239. (3) Altenbach, C.; Marti, T.; Khorana, H. G.; Hubbell, W. L. Science 1990, 248, 1088. (4) (a) Teng, C.-L.; Bryant, R. G. J. Am. Chem. Soc. 2000, 122, 2667. (b) Teng, C.-L.; Hinderliter, B.; Bryant, R. G. J. Phys. Chem. A 2006, 110, 580. (5) (a) Gerig, J. T. Prog. NMR Spectrosc. 1994, 26, 293. (b) Lau, E. Y.; Gerig, J. T. J. Am. Chem. Soc. 2000, 122, 4408. (6) (a) Prosser, R. S.; Luchette, P. A.; Westerman, P. W. Proc. Natl. Acad. Sci. U.S.A. 2000, 97, 9967. (b) Prosser, R. S.; Luchette, P. A.; Westerman, P. W.; Rozek, A.; Hancock, R. E. W. Biophys. J. 2001, 80, 1406. (7) Luchette, P. A.; Prosser, R. S.; Sanders, C. R. J. Am. Chem. Soc. 2002, 124, 1778. (8) Evanics, F.; Kitevski, J. L.; Bezsonova, I.; Forman-Kay, J. D.; Prosser, R. S. Biochim. Biophys. Acta 2007, 1770, 221. (9) Prosser, R. S.; Luchette, P. A. J. Magn. Reson. 2004, 171, 225. (10) (a) Evanics, F.; Hwang, P. M.; Cheng, Y.; Kay, L. E.; Prosser, R. S. J. Am. Chem. Soc. 2006, 128, 8256. (b) Bezsonova, I.; Evanics, F.; Marsh, J. A.; Forman-Kay, J. D.; Prosser, R. S. J. Am. Chem. Soc. 2007, 129, 1826. (11) Rinkevicius, Z.; Vaara, J.; Telyatnyk, L.; Vahtras, O. J. Chem. Phys. 2003, 118, 2550. (12) Telyatnyk, L.; Vaara, J.; Rinkevicius, Z.; Vahtras, O. J. Phys. Chem. B 2004, 108, 1197. (13) Moon, S.; Patchkovskii, S. In Calculation of NMR and EPR Parameters. Theory and Applications; Kaupp, M., Bu¨hl, M., Malkin, V. G., Eds.; Wiley-VCH: Weinheim, Germany, 2004; Chapter 13, p 325. (14) Pennanen, T. O.; Vaara, J. J. Chem. Phys. 2005, 123, 174102. (15) Hroba´rik, P.; Reviakine, R.; Arbuznikov, A. V.; Malkina, O. L.; Malkin, V. G.; Ko¨hler, F. H.; Kaupp, M. J. Chem. Phys. 2007, 126, 024107. (16) Pennanen, T. O.; Vaara, J. Phys. ReV. Lett. 2008, 100, 133002. (17) (a) Ramsey, N. F. Phys. ReV. 1950, 78, 699. (b) Helgaker, T.; Jaszun´ski, M.; Ruud, K. Chem. ReV. 1999, 99, 293. (18) For example, see: (a) Bu¨hl, M.; Parrinello, M. Chem.sEur. J. 2001, 7, 4487. (b) Bu¨hl, M.; Mauschick, F. T. Phys. Chem. Chem. Phys. 2002, 4, 5508. (c) Bu¨hl, M. J. Phys. Chem. A 2002, 106, 10505. (d) Pennanen, T. S.; Vaara, J.; Lantto, P.; Sillanpa¨a¨, A. J.; Laasonen, K.; Jokisaari, J. J. Am. Chem. Soc. 2004, 126, 11093. (e) Yazyev, O. V.; Helm, L.; Malkin, V. G.;

Li et al. Malkina, O. L. J. Phys. Chem. A 2005, 109, 10997. (f) Pennanen, T. S.; Lantto, P.; Sillanpa¨a¨, A. J.; Vaara, J. J. Phys. Chem. A 2007, 111, 182. (g) Li, X.; Rinkevicius, Z.; Tu, Y.; Tian, H.; Ågren, H. J. Phys. Chem. B 2008, 112, 11347. (19) The concentration of O2 in aqueous solution is given by CO2 ) kHPO2, where kH is Henry’s Law constant (kH ) 1.3 × 10-3 mol · L-1 · atm-1 at 298 K). (20) (a) Parrinello, M.; Rahman, A. J. Appl. Phys. 1981, 52, 7182. (b) Nose´, S.; Klein, M. L. Mol. Phys. 1983, 50, 1055. (21) (a) Nose´, S. Mol. Phys. 1984, 52, 255. (b) Hoover, W. G. Phys. ReV. A 1985, 31, 1695. (22) (a) Darden, T.; York, D.; Pedersen, L. J. Chem. Phys. 1993, 98, 10089. (b) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. J. Chem. Phys. 1995, 103, 8577. (23) (a) van der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. J. Comput. Chem. 2005, 26, 1701. (b) Lindahl, E.; Hess, B.; van der Spoel, D. J. Mol. Model. 2001, 7, 306. (c) Berendsen, H. J. C.; van der Spoel, D.; van Drunen, R. Comput. Phys. Commun. 1995, 91, 43. (24) van der Spoel, D.; Lindahl, E.; Hess, B.; van Buuren, A. R.; Apol, E.; Meulenhoff, P. J.; Tieleman, D. P.; Sijbers, A. L. T. M.; Feenstra, K. A.; van Drunen, R.; Berendsen, H. J. C. Gromacs User Manual, version 3.3; www.gromacs.org (2005). (25) (a) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. J. Comput. Chem. 2004, 25, 1157. (b) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. J. Comput. Chem. 2005, 26, 114. (26) Bayly, C. I.; Cieplak, P.; Cornell, W. D.; Kollman, P. A. J. Phys. Chem. 1993, 97, 10269. (27) (a) Becke, A. D. J. Chem. Phys. 1993, 98, 5648. (b) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. J. Phys. Chem. 1994, 98, 11623. (28) (a) Hehre, W. J.; Ditchfield, R.; Pople, J. A. J. Chem. Phys. 1972, 56, 2257. (b) Hariharan, P. C.; Pople, J. A. Theor. Chim. Acta 1973, 28, 213. (c) Francl, M. M.; Petro, W. J.; Hehre, W. J.; Binkley, J. S.; Gordon, M. S.; DeFrees, D. J.; Pople, J. A. J. Chem. Phys. 1982, 77, 3654. (29) (a) Cance`s, M. T.; Mennucci, B.; Tomasi, J. J. Chem. Phys. 1997, 107, 3032. (b) Cossi, M.; Scalmani, G.; Rega, N.; Barone, V. J. Chem. Phys. 2002, 117, 43. (30) (a) Dunning, T. H., Jr. J. Chem. Phys. 1989, 90, 1007. (b) Woon, D. E.; Dunning, T. H., Jr. J. Chem. Phys. 1993, 98, 1358. (31) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A., Jr.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian 03, revision D.01; Gaussian, Inc.: Wallingford, CT, 2004. (32) Zasetsky, A. Y.; Svishchev, I. M. Chem. Phys. Lett. 2001, 334, 107. (33) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (34) Impey, R. W.; Madden, P. A.; Mcdonald, I. R. J. Phys. Chem. 1983, 87, 5071. (35) Wolinski, K.; Hilton, J. F.; Pulay, P. J. Am. Chem. Soc. 1990, 112, 8251. (36) Rinkevicius, Z.; Telyatnyk, L.; Vahtras, O.; Ågren, H. J. Chem. Phys. 2004, 121, 7614. (37) DALTON, Release 2.0; http://www.kjemi.uio.no/software/dalton/ dalton.html (2005). (38) (a) Huzinaga, S. J. Chem. Phys. 1965, 42, 1293. (b) Helgaker, T.; Jaszun´ski, M.; Ruud, K.; Goˇrska, A. Theor. Chem. Acc. 1998, 99, 175. (39) Al-Abdul-Wahid, M. S.; Yu, C.-H.; Batruch, I.; Evanics, F.; Poms, R.; Prosser, R. S. Biochemistry 2006, 45, 10719. (40) Flyvbjerg, H.; Petersen, H. G. J. Chem. Phys. 1989, 91, 461.

JP902659S