Analysis of local molecular motions of aromatic sidechains in proteins by 2D and 3D fast MAS NMR spectroscopy and quantum mechanical calculations

We report a new multidimensional magic angle spinning NMR methodology, which provides an accurate and detailed probe of molecular motions occurring on timescales of nanoto microseconds, in sidechains of proteins. The approach is based on a 3D CPVC-RFDR correlation experiment recorded under fast MAS conditions (nR = 62 kHz), where C–H CPVC dipolar lineshapes are recorded in a chemical shift resolved manner. The power of the technique is demonstrated in model tripeptide Tyr–(D)Ala–Phe and two nanocrystalline proteins, GB1 and LC8. We demonstrate that, through numerical simulations of dipolar lineshapes of aromatic sidechains, their detailed dynamic profile, i.e., the motional modes, is obtained. In GB1 and LC8 the results unequivocally indicate that a number of aromatic residues are dynamic, and using quantum mechanical calculations, we correlate the molecular motions of aromatic groups to their local environment in the crystal lattice. The approach presented here is general and can be readily extended to other biological systems.


Introduction
In recent years solid-state Nuclear Magnetic Resonance (NMR) spectroscopy has become one of the most powerful techniques in structural studies of condensed matter, in particular for synthetic polymers and biological systems. 1 A number of highly exciting applications have been reported for the latter group of samples including characterization of proteins and other biopolymers that are inaccessible to X-ray crystallography and solution-state NMR spectroscopy. 2 Today the palette of available methodological NMR approaches revealing general and specific information about local geometry, intra-and intermolecular interactions, and higher order structure is very rich. 3 The full assignment of resonances and determination of molecular structure can be done by employing different multi-dimensional solid-state NMR protocols. This issue was exhaustively discussed in recent review articles by us and by others. [4][5][6][7] To understand in detail the biological functions of proteins, knowledge of the molecular dynamics is crucial. It is well known that protein motions occurring on the pico-to millisecond time scales play important functional roles, and that dynamics modulates interactions between physiological partners. The unique feature of NMR spectroscopy is the possibility of the analysis of the dynamics of biomolecules at atomic resolution. 8 Once signals are assigned, measurements can be made of relaxation rates, residual anisotropic interactions (dipolar, quadrupolar and chemical shift), relaxation dispersion profiles, magnetizationand hydrogen exchange rates; all yielding site-specific information about dynamics. However, one of the biggest achievements in the field of solid-state NMR spectroscopy is the introduction of fast magic angle spinning (MAS) probes presently delivering rotation frequencies higher than n R = 110 kHz. 9 These very fast rates have opened new possibilities in studies of proteins. 10 On the other hand, anisotropic tensorial interactions, e.g. dipolar and chemical shift, which are the source of crucial information about the structure and dynamics of natural products are eliminated or averaged to their isotropic values due to very fast sample spinning. Hence, a number of techniques have been developed during recent years to reintroduce those interactions and a large variety of pulse schemes have been published. [11][12][13][14][15] Recoupling of anisotropic spin interactions is achieved by application of pulse sequences that are often based on symmetry rules introduced by Levitt. 16 These sequences are called CN n n and RN n n , where the symbols N, n, and n are small symmetry integer numbers that depend on the rotation properties of the spin angular moment during the rotor-synchronized train of radio frequency (RF) pulses. In this strategy the relationship between spins, space, rotation and radiofrequency pulses is used to amplify the desired and suppress the unwanted interactions. The recently published RN-CSA symmetry schemes are suitable for CSA recoupling under fast spinning speed (n R Z 40 kHz). 17 Symmetry-based approaches were also employed by several research groups to study hetero-nuclear dipolar couplings in proteins. 18 These sequences are highly effective in suppressing 1 H-1 H homonuclear dipolar couplings, which make them suitable in fully protonated systems. In these symmetry-based experiments the RF pulses have to be synchronized with the sample rotation, which is not always easy to achieve in case of very fast sample spinning (e.g. n R Z 100 kHz) due to the proportionality between the required RF fields and n R , and also to the transients which are not negligible with such small rotor periods.
In a recent article employing small molecules as model samples, we have shown that a simple Cross Polarization with Variable Contact (CPVC) experiment performed under fast MAS conditions is free from the discussed supra-limitations and provides very accurate values of hetero-nuclear ( 13 C-1 H, 15 N-1 H) dipolar couplings, D, which reflect the inter-atomic distances and/or molecular motions. 19 We have also shown that CPVC methods are very robust with respect to RF-inhomogeneity, CSAs and offsets. Recent papers, from us and others, have fully confirmed the applicability of this approach. [20][21][22] In the current work we show that this methodology is a starting point for more advanced multi-dimensional techniques, which can be applied to study large biological systems. The power of this methodology is illustrated on 13 C enriched models: tripeptide Tyr-(D)Ala-Phe and two proteins, GB1 and LC8. We demonstrate that using our approach, we can site-specifically assign the side-chain resonances of the aromatic residues, which is typically a challenge in proteins. With these assignments in hand, we infer the molecular dynamics of the aromatic residues. To do so, we employ an average Hamiltonian approach to simulate the dipolar lineshapes. This approach permits unequivocal identification of the motional modes of the aromatic sidechains undergoing nano-to microsecond time scale motions. Furthermore, by combining our experimental NMR results with theoretical calculations, we show that the motional modes correlate with the location of a particular amino acid in the context of the protein structure. The methodology reported here is general and readily extendable to other biological systems of interest, including non-crystalline proteins, protein assemblies, and other bio-macromolecules. At this stage of hardware development, the 3D CPVC approach (as most of multi-dimensional solid-state NMR techniques) can only be applied to isotopically enriched samples. However, introducing new NMR concepts and techniques (e.g. very fast DNP CP/MAS) may shortly extend the applicability of the method to unlabeled samples.

Materials
U-13 C Tyr-(D)Ala-Phe (YAF) tripeptide was synthesized using labeled amino acids purchased from Sigma-Aldrich (D-alanine-13 C 3 -99%) and Cambridge Isotope Laboratories, Inc. (L-tyrosine-13 C 9 -99%, 15 N-99% and L-phenylalanine-13 C 9 -99%, 15 N-99%). Its synthesis was carried out in solution according to procedure reported in our previous paper. 23 U-13 C, 15 N-GB1 was bought from Giotto Biototech Company as ''ready-to-use'' sample packed in 1.3 mm rotor. Dynein light chain 8 (LC8) protein was prepared according to the protocol reported by us previously. 24 3 mg of U-13 C, 15 N-LC8 prepared by controlled precipitation from polyethylene glycol were packed into a 1.3 mm rotor and used for subsequent MAS NMR experiments.

NMR spectroscopy
All experiments were measured on a Bruker Avance III spectrometer equipped with a 1.3 mm 1 H/ 31 P-15 N probehead, spinning at n R = 62 kHz, and operated at frequencies of 600.1 and 150.9 MHz for 1 H and 13 C, respectively. In all cases the p/2 pulse durations were of 2.5 ms for 1 H and 1.5 ms for 13 C. We used a low power (RF = 10 kHz) swftppm decoupling 25

Simulations of CPVC dipolar lineshapes
The simulations of CPVC dipolar lineshapes were performed in the time-domain mode. Based on the following considerations we could do it by application of the simple algebraic expressions quoted below instead of having to go through timeconsuming full quantum mechanical calculations.
Under MAS, the CP magnetization transfer from the 1 H spin polarization to the to-be-detected 13 C spin polarization can be performed either in zero-quantum (ZQ) or double-quantum (DQ) mode. ZQ-CP is optimized under a so-called Hartmann-Hahn-sideband condition requiring that the difference of the RF amplitudes applied to the two nuclei is equal to once or twice the spinning speed. 28 DQ-CP requires that the RF-amplitude sum equals once or twice the spinning speed. 29 The CPVC method uses a sum or difference of RF-fields equal to the spinning rate as it reintroduces the maximum of dipolar interactions, hence leading to the maximum splitting of D = D/O2 in the case of static moieties. Moreover, it has been shown that ZQ-CP transfers perform better than DQ-CP ones, due to their larger possible RF fields.
In our recent publication 20 it was reported that spectra obtained under standard experimental conditions closely follow the results predicted under the assumption of ideal Hartmann-Hahn matching. It was also shown that non-ideal conditions, such as RF inhomogeneity, CSAs, resonance offsets, and homonuclear proton interactions, have very limited effects on the CPVC lineshape. Globally, the distortions are mainly observed in the region near zero frequency, leading us to exclude the central portion of the spectrum when performing the lineshape fittings. Another effect is a small increase of the apparent D value (less than 1 kHz), which is inconsequential and is neglected in our analysis.
These observations assured us that the current CPVC spectra can be analyzed in terms of the spin-evolution formalism of the idealized first-sideband Hartmann-Hahn situation, according to the intended RF-amplitude configuration. For a single 13 C-1 H spin pair, the 13 C signal after a CP contact time, t, is given by In the case of an immobile spin pair, the o CP frequency is given by o D being the dipolar coupling constant (DCC) in frequency unit and b the angle between the rotor axis and the CH vector. 30,31 In the case of a powder sample, two horns are observed in the dipolar spectrum, and their splitting is equal to D = o D /O2. As a result, the inter-nuclear distance is then directly determined from the splitting value without any lineshape fitting. In the absence of molecular motions, the dipolar-interaction tensor is axially symmetric, having components Ào D /2, Ào D /2, o D . When the CH bond undergoes very fast orientational motions, the size of the effective dipolar tensor is reduced to one having principal components D xx , D yy , D zz , thereby possibly losing its axial symmetry. 32 In such a case, we identify the largest of the absolute values of the three components as 'the' DCC, represented by o D and the difference of the two other components divided by o D , as the asymmetry parameter Z D .
The o CP frequency is then given by where a and b are the azimuthal and polar angles of the rotor-axis direction with respect to the dipolar-tensor principal axes.
We simulated the CPVC signal in a home-written Fortran-95 program by summation of S(t) over 2000 powder-averaged tensor orientations. This signal was then shifted by the experimental dead time, followed by Lorentzian apodization, and FT. Since the term 1 in the expression for S(t) creates the large time-independent signal component that in the experimental spectra is eliminated by appropriate baseline correction prior to the FT, we simulated S(t) simply as cos(o CP t). The experimental CPVC spectra were obtained by FT of the time-domain signals extracted from the original 3D datasets for each well-resolved aromatic peak. The best-fit parameters for the dipolar coupling constant (o D ), dipolar asymmetry parameter (Z D ), and line broadening parameter (LB) were obtained by simplex leastsquares minimization. Since experimental artifacts such as RF inhomogeneity and chemical-shift offsets often cause lineshape distortions in the central region, we disregard those regions in the fitting procedure. The width of this region is determined by the particular lineshape features in each individual case. We note that dipolar lineshapes corresponding to the desired carbon atom can be obtained using any cross peak with the frequency corresponding to this atom in the indirect or direct 13 C dimension, and, as discussed in the Results section, the dipolar parameters derived from the different cross peaks corresponding to the same carbon atom are consistent to within the experimental error.

Quantum mechanical calculations
The energy of the infinite crystal of YAF tripeptide was calculated by the CASTEP program, 33 which is a DFT-based code. Electron-correlation effects were modeled using the generalized gradient approximation of Perdew, Burke, and Ernzerhof. 34 The calculations were performed with 'on-the-fly' pseudopotentials, and a planewave cutoff energy of 600 eV with integrals taken over the Brillouin zone using a Monkhorst-Pack 35 grid with a minimum k-point sampling of 0.05 Å À1 . The atomic coordinates were derived from the Cambridge Crystallographic Database 36 (refcode IWANID 23 ). The positions of all atoms were optimized at the same computational level prior to the aromatic ring rotations. The energy dependence on the PHE/TYR rotation was calculated on a 2 Â 2 Â 1 supercell so that the rotation of one phenyl ring did not result in the same rotation of the neighboring rings under the periodic conditions. The energy calculations were performed with and without the dispersion correction. 37 For the cluster calculations, an aromatic ring was selected, then the selection was expanded to all atoms closer than 5 Å (including the atoms of the periodic images in the neighboring unit cells), and then the selection was expanded to include complete amino acid residues. Free C-termini amino acids were ended with amide (CONH 2 ) groups. The energy dependence on the phenyl ring rotation around the C b -C g bond was calculated at Hartree-Fock level using a standard 6-31G* basis set within the Gaussian-09 program package. 38 The orientation of the OH group of TYR residues was kept fixed so that the hydrogen bond network of the crystal was not disturbed.
The same procedure for cluster constructions was used for the GB1 and LC8 proteins. Crystal structures 2QMT and 3DVT from Protein Data Bank 39 served as the starting structures. 40 Hydrogen atoms were added manually and their positions were optimized at the same computational level prior to the ring rotations and energy calculations.

Results and discussion
NMR spectroscopy 2D CPVC spectra. In the previous reports the 2D CPVC NMR approach was applied to study structural constraints in relatively simple model samples. 19,20 In these cases, very well resolved 13 C isotropic signals have allowed unambiguous assignment of specific sites and analysis of C-H dipolar splittings. In the case of proteins, recording dipolar interactions site specifically requires at least one additional chemical shift dimension, to alleviate the signal overlap. This problem is illustrated here using two proteins under investigation, GB1 and LC8. GB1 is a 56-residue B1 binding domain of protein G, an immunoglobulinbinding protein expressed in group C and G streptococcal bacteria. LC8 is an 89-residue microtubule-associated protein, dynein light chain 8 (LC8), which is an integral subunit of cytoplasmic dynein, a multi-subunit microtubule-based retrograde motor.
The 2D CPVC spectra of GB1 recorded with the spinning frequency of o R = 62 kHz are shown in Fig. 1a (full spectrum) and Fig. 1b and c (expansions around aromatic and aliphatic regions). The appropriate spectra for LC8 protein are attached in the ESI † (Fig. S2).
On the first glance, the spectral resolution appears to be sufficiently high as to permit assignment of almost all 13 C resonances. However, upon further inspection of the spectrum it becomes clear that we can only determine coarse values of dipolar splittings in the aromatic and aliphatic regions. Inspection of F1 slices indicates that peaks corresponding to the individual amino acid residues overlap, precluding extraction of site-specific information concerning the different motional modes. This overlap is further exacerbated in LC8, due to the larger number of residues. As one can see in Fig. 1b and c, both aromatic and aliphatic residues are overcrowded, and so the dipolar lineshapes for the individual carbon sites overlap.
3D CPVC-based experiments for analysis of aromatic sidechain dynamics: validation in YAF tripeptide. To overcome the spectral overlap problem, we introduced two new 3D experiments. As shown in Fig. 2, we have combined CPVC with homonuclear 13 C-13 C sequences, RFDR (Fig. 2b) and SHANGHAI (Fig. 2c), which allow for the determination of short and long 13 C-13 C distances, respectively. In the current paper we focus on the 3D CPVC-RFDR approach.
The applicability and accuracy of the CPVC-RFDR method were first validated using a simple model compound, the U-13 C labeled YAF tripeptide. This sample is an excellent test case due to the well-defined and distinct dynamics of the two aromatic residues. 41 Indeed, the phenyl ring of Tyr is rigid, while that of Phe is mobile and undergoes p jumps at the picosecond rate at room temperature. 41 These fast (on the NMR time scale) motions are expected to give rise to partial averaging of the Phe C-H dipolar tensors associated with the Cd1, Cd2 and Ce1, Ce2 carbon atoms (the two pairs of resonances overlap), resulting in a reduced dipolar frequency of o D = 11.25 kHz, which is ca. half of the rigid-limit value of 22.5 kHz, and a non-zero dipolar asymmetry parameter, Z D . Fig. 3 shows the 2D RFDR plane from the 3D CPVC-RFDR spectrum of the YAF peptide and five selected CPVC lineshapes extracted from the 3D spectrum. The dipolar interaction parameters,  Table S2 of the ESI. † The dipolar parameters for the Tyr and Phe residues are very different; the reduced dipolar frequency and non-zero dipolar asymmetry parameter for Phe indicate the presence of dynamics. The reduction of the dipolar splitting due to the molecular motions of phenyl rings in organic polymers and biological systems is well known and was discussed in prior reports including our previous paper. [41][42][43][44] The Cd1, Cd2 and Ce1, Ce2 carbons of the Phe residue exhibit motionally reduced dipolar interaction parameters: o D = 12.3/12.2 kHz and Z D = 0.42/0.45. These dipolar parameters suggest composite motions consisting of the p-jumps and additional motional modes, vide infra. In contrast, the Cz carbon, whose dipolar tensor does not experience such motional averaging due to symmetry (Fig. 3A), exhibits dipolar interaction parameters very close to the rigid limit values, o D = 22.9 kHz and Z D = 0.0. As discussed later in the text, this finding provides an important clue in defining the motional models for this residue. For the tyrosine residue, the dipolar interaction parameters for Cd1, Cd2 and Ce1, Ce2 carbons indicate an almost rigid phenyl ring, with o D = 23.2/22.8 kHz and Z D = 0.20/0.16, which is close (but not identical) to the rigid-limit values of 23 kHz and 0.0, respectively. That tyrosine's Cd1, Cd2 and Ce1, Ce2 as well as the phenylalanine's Cz show slightly depressed dipolar couplings and non-zero dipolar asymmetry parameters (the deviations are well outside of the experimental errors) suggest the presence of some mobility in the ns-ms regime, likely small-amplitude out-of-plane wobbling motions.
Analysis of aromatic sidechain dynamics in proteins: GB1 and LC8. With the 3D CPVC-RFDR methodology validated on the model YAF tripeptide, we have next turned our attention to the analysis of aromatic sidechain dynamics in nanocrystalline GB1. The starting point of the analysis is the assignment of the 13 C signals. [45][46][47] We note that GB1 has been extensively studied by MAS NMR, and resonance assignments have been established by the Rienstra's group previously. 48 In this work, we focused attention on the dynamics of aromatic amino acids side groups. The molecular motion of other residues, e.g. Ca and Cb carbons, can be probed using our approach (vide supra) as well as other 3D MAS NMR methodologies discussed elsewhere. 49,50 To assign the aromatic sidechain carbons, we employed the 13 C-13 C 2D RFDR correlation, which yields one-bond C-C connectivities for all amino acid residues. Fig. 4 shows the expanded aromatic region with labels corresponding to selected quaternary atoms. GB1 contains six aromatic residues: three tyrosines (Y3 and Y45 in b-sheets and Y33 in an a-helix), two phenylalanines (F30 in an a-helix and F52 in a b-sheet), and one tryptophan (W43, in a b-sheet). The quaternary carbons, which are separated from other signals, were assigned first, and the chemical shifts of Cz carbons are d = 158.5 (Y33Cz), 158.0 (Y3Cz), and 157.5 ppm (Y45Cz). Using these Cz resonances, the assignments of the aromatic sidechain spin systems for the three Tyr residues were straightforward.
The assignment of Phe-52 is also straightforward. The chemical shift of F52Cg is 139.9 ppm, and it correlates with its Ce/Cz whose shift is 131.0 ppm. Other correlations for Phe appear as intense overlapping signals found in the region of 130-132 ppm, and we could not unequivocally assign Phe-30 resonances. The assignment of resonances for Trp-43 is unambiguous. The correlation peaks between quaternary C43e carbon (d = 138.2 ppm) and other carbons are clearly present. Fig. 5 shows the CPVC dipolar lineshapes for the aromatic carbons in GB1 extracted from the 3D CPVC-RFDR spectrum. Indeed, the dipolar parameters were obtained as best fits by the AHT-based simulations of the lineshapes, and are summarized in Table S2 of the ESI. † The results indicate that W43 does not undergo any significant dynamics. The dipolar interaction parameters extracted for all carbon atoms that have attached protons are very close to the rigid limit values, o D = 23 kHz and Z D = 0. The observed deviations from these parameters for GB1 are attributed to experimental errors, which are larger than in the YAF tripeptide because of lower signal-to-noise ratio in the spectra. The fact that the CPVC lineshapes for W43 and the corresponding dipolar interaction parameters extracted from 14 different diagonal and cross peaks are consistent is an important validation of the robustness of the CPVC-RFDR experiment. We next discuss the analysis of the molecular dynamics of the aromatic residues in LC8. Fig. 6 shows the expansion around the aromatic region of the 13 C-13 C 2D RFDR correlation spectrum. The full spectrum is displayed in Fig. S1 of the ESI. † LC8 contains six tyrosines (Y32, Y41, Y50, Y65, Y75, and Y77), five phenylalanines (F46, F62, F73, F76, and F86), three histidines (H55, H68, and H72), and one tryptophan (W54). We have reported the residue-specific assignments of the aliphatic and carbonyl carbons of LC8 previously. 51,52 However, in those earlier data sets aromatic sidechain assignments could not be obtained. As shown in Fig. 6, the quality of the 2D RFDR spectrum is such that the assignments of the majority of the aromatic sidechains were readily accomplished. The chemical shifts are summarized in Table S1 of the ESI. † With these assignments in hand, the analysis of the dipolar lineshapes from the 3D CPVC-RFDR spectra was straightforward for many of the aromatic residues in LC8. Fig. 7 shows the experimental and simulated dipolar lineshapes for representative residues that yield resolved cross peaks. The corresponding dipolar parameters are summarized in Table S2 of the ESI. † The only tryptophan residue in LC8, W54, exhibits o D = 21.8-22.0 kHz and Z D = 0.06-0.20. These values indicate that W54 is nearly rigid in LC8, which is consistent with the fact that it is located in a b-sheet and packed against several other aromatic residues: F46, Y50, and H70. Similarly, in H55 the field values of o D = 21.9 AE 0.2 kHz and Z D = 0.11 AE 0.05 indicate the lack of large-amplitude motions of the aromatic ring. This is not very surprising given the fact that H55 (i) belongs to a b-sheet, (ii) comprises the interface of the LC8 homodimer, and (iii) is packed against three aromatic residues: F86 located in the same LC8 monomer as well as F86 and H55 of the neighboring molecule.
F76, the only phenylalanine with non-overlapping signals whose dipolar lineshapes we have been able to analyze, is dynamic. The dipolar parameters recorded for Cd1/Cd2, o D = 14.6 kHz and Z D = 0.05, indicate complex motions similar to those related to our observations for Y33 of GB1, and the model that will be discussed in the next section apply.
For tyrosines, spectral overlap precludes analysis of dipolar lineshapes for Y32 and Y50. For the other four tyrosines, Y41, Y65, Y75, and Y77, the results unequivocally indicate the presence of dynamics. The dipolar parameters for all tyrosines suggest complex motional modes, the models for which are discussed in the next section.
Motional models for aromatic groups. As discussed above, for several residues in GB1 and LC8 and for F3 in YAF we observed motionally reduced dipolar interaction parameters, namely o D smaller than the rigid limit value of 22.5 kHz and a non-zero Z D value. Non-zero dipolar asymmetry arising in the presence of dynamic averaging has a profound effect on the lineshapes as it introduces an apparent broadening and it shifts the positions of the peak splittings, as illustrated in Fig. 8A.
Under such circumstances, it is not possible to determine the dipolar interaction parameters from the peak splittings alone, as has been suggested in prior studies, 19 and lineshape simulations are required. As discussed earlier in the text, accurate dipolar interaction parameters are derived by simulations of CPVC spectra using the average Hamiltonian approach.
Interestingly, the dipolar interaction parameters for many carbon sites in our experiments are indicative of motions more complex than p-jumps around the Cg-Cz axis of the aromatic ring, where o D = 14.0 kHz and Z D = 0.6. To describe the possible motional models giving rise to such dipolar interaction parameters, we define an xyz axis system in the aromatic ring illustrated in Fig. 8B: z -along the Cg-Cz axis; x -perpendicular to z in the aromatic plane; and y -perpendicular to the aromatic plane.
The models are kept as simple as possible by restricting them to rapid jumps of the dipolar tensor between two equally occupied orientations that are rotated away over angles AEf i from an unoccupied central orientation, whereby the rotation axis is parallel to one of the principal tensor axes (i = x, y or z).
If the rotation is about the z axis, then the conversion of the dipolar tensor from the static central situation to the average of the rotated situations is given by A similar averaging is obtained for rotation about x or y axes. The main motional process of the tyrosine ring is 1801 flips about the z-axis. For the four carbons that are not along the z-axis, this can be viewed as rotational jumps about a central position consisting of an axially symmetric dipolar tensor with its main axis along z, making f y = AE601 jumps about y. These parameters deviate from the experimentally observed values by more than the experimental error.
Moreover, we cannot ascribe the deviations to p-ring flips that occur at correlation times between the extreme-narrowing limit of 10 À7 s and the rigid limit of 10 À4 s, thereby creating interesting transition lineshapes. Indeed, we have performed such simulations (not shown), and not only do the lineshapes poorly fit the experiment, but these intermediate correlation times also give rise to very short T 1r relaxation times that produce huge negative center peaks in the CPVC spectra. Those are not observed in our experimental data.
In the case of Y33 (GB1) and F76 (LC8), o D value found by the fit equals that of the two-site flip model, but Z D is nearly zero instead of the 0.6 predicated by the model. This suggests that there is an additional rotation about the x-axis (the orientation where the dipolar component is not affected). Such upand-down flaps of varying angles give the following average dipolar tensors: Keeping the line broadening equal to the fitted value of LB = 1.7 kHz, the resulting simulated spectra assuming the above flip angles, plotted on the experimental spectrum for Y33 of GB1, are illustrated in Fig. 9.
A particular observation is that the peak splitting increases from 8.0 to 8.9 kHz when the rotation angle f x increases from 0 to 401. This contradicts the conventional assumption that a smaller splitting implies a smaller ''order parameter'', i.e., a larger amplitude of molecular motion. This result, counterintuitive on the first glance, underscores the importance of fitting the lineshapes for deriving accurate dipolar tensor parameters.
The above result suggests that, in addition to the 1801 flips, the phenyl rings in Y33 (GB1) and F76 (LC8) also flap up and down over an angle of about 2 Â 401. This could be related to trans-gauche-type transformations at C b and be accommodated by bond-angle changes of the OH group and its hydrogenbonding. We note that the relatively large line broadening parameter required to fit the spectra may either be a relaxation effect caused by these large-amplitude motions, or a distribution of the f x angles.
For F3 residue in YAF and for Y41, Y65, Y75, and Y77 in LC8, the CPVC lineshapes deviate in a different way from a simple p-jump model. In this case Z D is close to the predicted 0.6 but o D is reduced to 10.4-12.4 kHz. This can be reproduced by applying one additional small-angle two-site jump. In the case of F3 in YAF it is clear that the second rotation has to be about the z-axis, because the Cz carbon exhibits rigid-limit dipolar interaction values, which would be inconsistent with rotations about either the x-or the y-axis. Incorporating the second rotation about the z-axis into our motional model for the YAF's F3 residue, we obtain the best fit to the experimental results when the ring undergoes additional rotations by AE151 about the z-axis. So, we may conclude that when the phenyl ring makes p flips about the z-axis it has a wide range of possible orientations around its energy minimum, ca. 151 to each side, in agreement with the very shallow minimum around 0 and 1801 found in the quantum mechanical calculations, vide infra.  For tyrosine residues in LC8, we cannot constrain the motional model in a similar way because Cz carbons do not carry a proton. In the absence of such constraints, more than one motional profile is consistent with the experimental data.
Theoretical calculations of the energy barriers for mobile aromatic groups. To corroborate the motional models observed experimentally (vide supra) we have turned our attention to the quantum chemical calculations of rotational energy barriers. We have conducted these calculations on GB1, as the 3D structure of the polymorph used in our NMR experiments is known, in contrast to that of the LC8 sample used in the current work.
An alternative computational approach for the assessment of the side-chain motions could be molecular dynamics simulations with their tremendous advances in performance and applicability in recent years. [53][54][55] However, the size of the observed dipolar couplings suggests that micro-to millisecond simulation times would be required for a proper capturing of the molecular motions, which is still beyond the capability of conventional simulation techniques. The below described calculations of rotational barriers thus represent a computationally cheap and widely applicable approach for a qualitative link between observed NMR data and side-chain mobility.
In explaining the differences in the C-H dipolar splittings of aromatic groups in the crystal lattice we assumed that even subtle distinctions in the local environment will influence the energy barrier of free rotation around the 1-4 axis of phenyl rings. The aim of the study reported in this section is an analysis of energetic parameters employing the theoretical approach and their correlation with experimental values of the DCC. For testing the methodology we have used X-ray data for YAF peptide (refcode IWANID) constructing the appropriate cluster surrounding each phenyl residue. Fig. 10 shows the crystal lattice motif with clearly visible phenyl groups of Phe inside a cage. The aromatic residue has enough space for free rotation. The molecular motion of phenyl ring of Tyr is more restricted due to presence of the neighboring methyl groups and hydrogen bonds.
Two different computational approaches were used to estimate the mobility of the aromatic rings in YAF peptide. First, crystal calculations with fully periodic boundary conditions of the optimized YAF structure and of that with changed orientation of Phe or Tyr aromatic rings were used to estimate the barrier of aromatic ring rotation. Significant differences between the barriers of rotation of aromatic rings of Tyr and Phe are clearly visible in Fig. 11.
Indeed, the rotation barrier of Phe is much lower than that of Tyr, which is in agreement with previous observations of significant differences in the mobility of these phenyl rings in  the YAF peptide. 23,41 The rotation barriers in a real system will be probably slightly lower than the calculated ones for both rings because the surrounding molecules can adapt to the transition geometry.
The second approach of the rotation barrier estimation consists in energy calculations for a cluster of molecules cut out of the infinite crystal structure with the studied aromatic ring in the center of the cluster. This approach is commonly used for prediction of NMR parameters of crystalline samples. 56 The Hartree-Fock energies calculated for the cluster models of YAF peptide overestimate the rotation barriers slightly (Fig. 12), but for a qualitative insight into the mobility of aromatic residues within a polypeptide structure these calculations are sufficient. The advantage of the cluster approach is its applicability to protein structures where the calculations with periodic boundary conditions would be prohibitively computationally expensive. The comparison with the energies calculated with the periodic approach validates the cluster model, which is then applied for the GB1 protein.
Analyzing the surrounding of aromatic residues we have found that the ring of Tyr33 faces the area between a helix and a sheet (Fig. 13). In this position the phenyl ring has enough space for molecular motion, which is confirmed by the cluster calculations. The energy barrier for Tyr33 rotation was found to be less than half of that for Tyr45 and Phe52 (Fig. 14). Tyr45 is located in the neighborhood of Phe52. The distance between the aromatic groups is below 4 Å and favors the p-p contacts. The rings under discussion are not aligned in one plane but form a small angle; and it is known from literature that a small offset strengthens such interactions. 57     Thus we can assume that molecular motions of Tyr45 and Phe52 are limited to small-amplitude wobbling around the 1-4 axis. Phe52 has an additional T-shape short-range contact (3.3 Å) with the adjacent tryptophan residue (Fig. 13). The order of the calculated barrier heights agrees well with the order of experimental dipolar couplings (Fig. 14) for these three residues. Tryptophan contains a bulky aromatic system, which cannot rotate freely in the GB1 protein.
We have tried to apply the same procedure for the analysis of LC8 protein trying to find correlations between NMR measurements and quantum mechanical calculations of rotational barriers. Unfortunately, in the case of LC8 protein, the calculations are complicated by the fact that several polymorphic structures of LC8 are known. We have chosen the 3DVT structure for the calculations, but the calculated rotational barriers were not fully consistent with the experiment, which indicates that the 3DVT structure was not a good representative of the experimental structure.
A knowledge of the 3D structure of the studied molecule is crucial for the understanding of its mobility. If, however, the structural data are missing, the calculations of rotational barriers might be possibly used to provide structural constraints or even to differentiate between several structural models.

Conclusions
In conclusion, under fast MAS conditions, the measurements of heteronuclear ( 15 N-1 H and 13 C-1 H) dipolar couplings in proteins can be obtained in a site-resolved way in a 3D experiment containing the CPVC dipolar dimension and two isotropic chemical shift dimensions. We have demonstrated that, for aromatic sidechains, 13 C-1 H dipolar couplings can be accurately recorded and that they provide information on dynamic processes occurring on the timescales of dipolar interaction, as we have illustrated for two proteins, GB1 and LC8. Remarkably, by the analysis of dipolar lineshapes, specific motional modes can be inferred. Furthermore, in YAF, GB1, and LC8, the observed dynamic behavior could be readily accounted for by the analysis of tertiary contacts formed between the aromatic sidechains and other amino acids, and these observations are generally consistent with quantum mechanical calculations of the rotational energy barriers. Broadly, we expect that the methods introduced in this work will be useful for probing nano-to microsecond time scale motions in a variety of biological systems.