Model-Free or Not?

Introduction

Study of biomolecular function requires understanding the dynamics of the biological system. Nuclear magnetic resonance (NMR), despite many recent technological advances in other techniques, remains a premier method for detailed dynamics characterization. In NMR, one may measure a variety of site-specific relaxation experiments, which provide timescale sensitive information about the motion. By varying the type of experiment (T₁, T_1ρ, NOE, etc.) or experimental conditions (external magnetic field, applied field strength, magic-angle spinning (MAS) frequency, etc.), the timescale sensitivity of the measurement is modified. Then, one may resolve the dynamics both in space, via site resolution, and in timescale, via multiple experiments (Palmer, 2004; Schanda and Ernst, 2016).

However, is it possible to fully characterize the motions leading to the observed relaxation behavior? Many relaxation experiments in NMR are sensitive to the reorientational motion of anisotropic NMR interaction tensors (NMR relaxation can also be sensitive to change in scalar terms, e.g., isotropic chemical shift). For a given spin, relaxation is usually dominated by only one to two interactions. For example, relaxation of ¹⁵N in a protein backbone is determined almost entirely by the reorientation of the one-bond ¹H–¹⁵N dipole coupling and the ¹⁵N chemical shift anisotropy (CSA). But, multiple sources of motion lead to reorientation of the bond. For example, if we suppose the H–N bond to be in a protein, within a helix, then we would have local distortion of the peptide plane (one-bond libration), motion of the peptide plane within the helix, motion of the helix within the protein, and motion of the protein either in solution, in a crystal, a fibril, a membrane, etc.

This degree of complexity is illustrated in Figure 1. For a given bond in a molecule, and a given motion acting on that bond, a distribution of orientations is sampled as illustrated in Figure 1A. The orientational distribution determines the contribution of that motion to the total order parameter, $S^2$. However, not only are there many orientations sampled by a bond due to a motion, but those orientations are sampled at some rate, such that the motion has an associated correlation time or distribution of correlation times (denoted $\left(1-{\mathit{S}}^2\right)\theta \left(\mathit{z}\right)$). We illustrate this in Figure 1B; note that not only is the width of such a distribution variable, but also the functional form of the distribution itself. This results in a correlation function that decays from 1 to $S^2$, where integrating over the distribution of correlation times yields the total amplitude of the decay. Already, a single bond with just one motion acting on it yields potentially a high degree of complexity; however, we must still consider that multiple motions act on each bond, where the total correlation function is the product of the correlation functions of each individual motion (if those motions are independent from one another, Figure 1C). Finally, motion varies throughout a molecule, as a function of position, resulting in a complex, multi-dimensional description as illustrated in Figure 1D.

FIGURE 1

FIGURE 1. Complexity of reorientational dynamics. For each bond in a molecule, multiple types of motion result in orientational sampling, where the distribution of angles for each motion result in a generalized order parameter, S². Therefore, in (A) we plot a possible distribution of Euler angles for a single type of motion (population is plotted as a function of angles β and γ, where α is not required for a symmetric interaction tensor). A single motion is furthermore described by a correlation time, and may be distributed over a range of correlation times. In (B) we plot a possible distribution of correlation times $\left(1-{\mathit{S}}^2\right)\theta \left(\mathit{z}\right)$, that is, amplitude of motion as a function of the log-correlation time, $\mathit{z}=\mathrm{lo}{\mathrm{g}}_{10}\left({\tau }_c/\mathrm{s}\right)$. Each distribution is characterized by an amplitude, center, and width. Note that the integral of the distribution is $\left(1-{\mathit{S}}^2\right)$, S² being determined by the distribution of angles in (A). While (A,B) illustrate aspects of a single motion, multiple motions influence a given bond, where the total correlation function is the product of individual correlation functions. In (C), we plot four distributions of motion (color). Above each motion, we plot the distribution resulting from the product of that motion and all motions below it (black), eventually resulting in the total distribution seen at the top. Finally, we note that the total distribution varies as a function of position in the molecule, resulting in the 3D plot of the distribution as a function of correlation time and position in the molecule observed in (D). While this is just an illustration, one could imagine that motion in (D) results from three α-helices in a protein, each having a slightly different behavior, and varying dynamics as one approaches the end of each helix.

While NMR is powerful, obtaining a complete description of the complex dynamics stretches beyond the limit of what is possible based on experimental data alone, especially for large molecules such as proteins. This problem is a familiar one, addressed almost 40 years ago by Lipari and Szabo (Lipari and Szabo, 1982a), who developed a method known as the model-free approach. While we will discuss the details of this approach below, the name tells us a critical advantage of such an approach: model-free analysis allows the extraction of dynamics information from NMR relaxation data without having knowledge of the specific model of motion. Furthermore, the resulting parameters have a well-defined relationship to the distribution of orientations sampled and the distribution of correlation times.

Lipari and Szabo described the internal motion of a molecule with just two parameters: a generalized order parameter related to the amplitude of motion, $S^2$, and a mean effective correlation time, $\langle {\tau }_e\rangle$ (a third parameter, ${\tau }_{\text{M}}$, gives the correlation time of the molecule tumbling in solution). While only two parameters suggests a simple analysis, it is important to note that Lipari and Szabo did not intend to only describe simple motions having just a single correlation time and amplitude: theoretical tests of their model were performed on a wobbling-on-a-cone model (Kinosita et al., 1977) that results in a weighted sum of correlation times, and experimental work was performed on methyl groups in a protein, for which the total motion is determined by the product of methyl rotation and by reorientation of the methyl group’s C–C bond. Rather, the two parameters contain the aggregated information describing all motions that is available from the set of relaxation experiments alone.

The advantage of model-free analysis is that it does not require knowing the model of motion. For example, for relatively low fields (∼90 MHz, as used by Lipari and Szabo), all distributions of orientations and correlation times shown in Figure 2 should yield identical relaxation rate constants for the set of experiments. If we do not know which model is the correct model, the best we can do is to parameterize the results in a way that does not depend on the model of motion, as can be done with the model-free parameters $S^2$ and $\langle {\tau }_e\rangle$.

FIGURE 2

FIGURE 2. Five distributions of orientations and correlation times that yield the same model-free parameters ($\left(1-{\mathit{S}}^2\right)$ = 0.3, $\langle {\tau }_e\rangle$ = 0.1 ns). In (A–E), we plot a distribution of orientations (sphere, right); on the axes, we plot the distribution of correlation times resulting from exchange among that set of orientations. Models of motion are wobbling-on-a-cone (${\theta }_{\text{cone}}={19}^{°}$), wobbling-in-a-cone (${\theta }_{\mathrm{cone}}={28}^{°}$), symmetric two-site hop (${\theta }_{\mathrm{hop}}={39}^{°}$), asymmetric two-site hop (${\theta }_{\mathrm{hop}}={70}^{°}$), and 6-site asymmetric exchange. Insets in (A,B) show correlation times with small amplitudes. $S_{\mathrm{resid}\mathrm{.}}^2$ refers to the order parameter from residual couplings (see Supplementary Section S3), which deviates from the generalized order parameter for asymmetric motion.

When analyzing data, a model provides a framework for understanding the data, and by using a model we are always adding some information to the experimental data. In some cases, we add further information depending on how we interpret a model. A model is advantageous if that information is correct, and disadvantageous if that information is wrong. Suppose, for example, we know that the correct model in Figure 2 is a symmetric two-site hop, shown in Figure 2C; then we may extract the hop angle and exchange rates from $S^2$ and $\langle {\tau }_e\rangle$, resulting in ${\theta }_{\text{hop}}$ = 39° and $k_{\text{ex}}^{1\to 2}$ = $k_{\text{ex}}^{2\to 1}$ = 5 × 10⁹/s. However, if the true model is an asymmetric two-site hop, shown in Figure 2D, the true angle and exchange rates may be significantly different (for Figure 2D, these are ${\theta }_{\text{hop}}$ = 70° with $k_{\text{ex}}^{1\to 2}$ = 1.3 × 10⁹/s and $k_{\text{ex}}^{2\to 1}$ = 8.7 × 10⁹/s). Then, a model-free approach is the more reliable method when the correct model cannot be independently determined.

In this review, we will first discuss the original model-free approach, and then examine methods descended from it, including discussion of our own detector analysis, a relatively new approach that also provides a model-free analysis in the spirit of the original Lipari-Szabo approach, but can extract the full information content of relaxation data sets in instances where the model-free approach cannot. We discuss analysis of microsecond motions using R_1ρ relaxation, and finally consider how other methods, in particular molecular dynamics (MD) simulation, may be used to supply the information that NMR lacks, thus improving the interpretation of NMR parameters.

Model-Free

While dynamics analysis methods have existed for application to solid-state NMR for some years now (Chevelkov et al., 2009b; Schanda et al., 2010; Zinkevich et al., 2013; Lamley et al., 2015a; Smith et al., 2016; Lakomek et al., 2017; Kurauskas et al., 2017), most of the approaches applied have evolved from methodology first developed for solution-state NMR. Probably the most important advance in solution-state analysis was the development of the model-free approach (Lipari and Szabo, 1982a; Lipari and Szabo, 1982b), and related two-step techniques (Wennerström et al., 1979; Halle and Wennerström, 1981; Brown, 1982). Then, we begin by reviewing some of the existing methodology, to understand advantages and disadvantages to various approaches.

Model-Free Theory

Typical solution-state NMR data sets consist of relaxation rate constants for R₁ (1/T₁), R₂ (1/T₂), and nuclear Overhauser effect (NOE, ${\sigma }_{\text{IS}}$), acquired at one or more magnetic fields. The rate constants describe the signal decay ($I\left(t\right)=I_0{e}^{-R_{\zeta }t}$) or recovery ($I\left(t\right)=I_{\text{eq}}+\left(I_0-I_{\text{eq}}\right)e^{-R_{\zeta }t}$). In solid-state NMR, this behavior can be multi-exponential, whereas we use the rate constant that describes the powder-averaged value (Krushelnitsky et al., 2018). Relaxation is often driven by reorientation of a few anisotropic interactions, for example, for backbone ¹⁵N relaxation, a one-bond H–N dipole coupling and CSA are responsible for relaxation. For these experiments, the relaxation rate constants may be calculated from the spectral density, $J\left(\omega \right)$:

$\begin{array}{l}{R}_1^{\text{I}}=\underset{\text{dipolar}\text{relaxation}}{\underbrace{{\left(\frac{{\delta }^{\text{IS}}}{4}\right)}^2\left(J\left({\omega }_{\text{I}}-{\omega }_{\text{S}}\right)+3J\left({\omega }_{\text{I}}\right)+6J\left({\omega }_{\text{I}}+{\omega }_{\text{S}}\right)\right)}}+\underset{\text{CSA}\text{relaxation}}{\underbrace{\frac{1}{3}{\left({\omega }_{\text{I}}\Delta {\sigma }_{\text{I}}\right)}^{2}J\left({\omega }_{\text{I}}\right)}}\\ R_2^{\text{I}}=\frac{1}{2}{R}_1^{\text{I}}+\underset{\text{dipolar}\text{relaxation}}{\underbrace{{\left(\frac{{\delta }^{\text{IS}}}{4}\right)}^2\left(3J\left({\omega }_{\text{S}}\right)+2J\left(0\right)\right)}}+\underset{\text{CSA}\text{relaxation}}{\underbrace{\frac{2}{9}{\left({\omega }_{\text{I}}\Delta {\sigma }_{\text{I}}\right)}^{2}J\left(0\right)}}\\ {\sigma }_{\text{IS}}=\underset{\text{dipolar}\text{relaxation}}{\underbrace{{\left(\frac{{\delta }^{\text{IS}}}{4}\right)}^2\left(-J\left({\omega }_{\text{I}}-{\omega }_{\text{S}}\right)+6J\left({\omega }_{\text{I}}-{\omega }_{\text{S}}\right)\right)}}\end{array}(1)$

Here, ${\omega }_{\text{I}}$ is the Larmor frequency (in radians/s) of the nucleus being relaxed, ${\omega }_{\text{S}}$ the Larmor frequency of the coupled spin (usually ¹H), and ${\delta }^{\text{IS}}$ and $\Delta {\sigma }_{\text{I}}{\omega }_{\text{I}}$ are the anisotropies of the dipolar coupling and CSA, respectively (${\delta }^{\text{IS}}=-2\frac{{\mu }_0}{4\pi }\frac{h{\gamma }_{\text{I}}{\gamma }_{\text{S}}}{r_{\text{IS}}^2}$, with ${\mu }_0$ the vacuum permeability in T²m³/J, ${\gamma }_{\text{I}},{\gamma }_{\text{S}}$ the gyromagnetic ratios of the two spins in radians/s, h is Planck’s constant in J·s, and $r_{\text{IS}}$ the distance between the spins in meters, resulting in ${\delta }^{\text{IS}}$, which is the full breadth of the dipolar powder pattern in radians/s. $\Delta {\sigma }_{\text{I}}{\omega }_{\text{I}}$ is similarly the full breadth $(\Delta {\sigma }_{\text{I}}=\frac{3}{2}\left({\sigma }_{\text{zz}}-{\sigma }_{\text{iso}}\right))$ of the CSA powder pattern in radians/s when the Larmor frequency of spin I is given by ${\omega }_{\text{I}}$, also in radians/s (Schanda and Ernst, 2016)). The spectral density may be obtained from the Fourier transform of the correlation function of motion. The correlation function itself is the rank-2 tensor correlation function, and describes the reorientational behavior of an NMR interaction tensor in time. If we assume the correlation function is symmetric in time, we may replace $e^{i\omega t}$ with $\cos \left(\omega t\right)$ in the Fourier transform. We can also change the integration bounds from $\left(-\infty ,\infty \right)$ to $\left(0,\infty \right)$, and must multiply the integral by two in order to compensate for only integrating over half the space.

$\begin{array}{l}J\left(\omega \right)={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}C\left(t\right)e^{i\omega t}dt}\\ ={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\underset{\text{symmetric}\text{in}\text{time}}{\underbrace{C\left(t\right)\cos \left(\omega t\right)}}dt}+i{\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\underset{\text{antisymmetric}\text{in}\text{time}\to 0}{\underbrace{C\left(t\right)\sin \left(\omega t\right)}}dt}\\ J\left(\omega \right)=2{\displaystyle \underset{0}{\overset{\infty }{\int }}C\left(t\right)\cos \left(\omega t\right)dt}\end{array}(2)$

Then, model-free analysis makes a few assumptions about the correlation function:

1) The total motion of a given bond is the result of overall tumbling of the molecule in solution and internal motion of the bond within the molecule, and these two motions are statistically independent.

2) Decay of the correlation function due to internal motion is fast compared to all $\omega$ sampled by the set of experimental relaxation rate constants (i.e., the extreme narrowing limit).

The decay of the correlation due to internal motion does not need to be mono-exponential (or even multi-exponential, although we will later apply this assumption). Instead of the second assumption, we may assume that the correlation function due to internal motion is mono-exponential, in which case we do not require its decay to be fast (we will visit this case only briefly, as it is less likely to occur in practice). We also assume tumbling is isotropic, although this is not necessarily required. Note that separate methods exist in case overall tumbling and internal motion are coupled (Tugarinov et al., 2001), although we will not consider these here. As a set of equations, this yields

$\begin{array}{l}C\left(t\right)=C^{\text{intern}\text{.}}\left(t\right)\cdot C^{\text{rot}\text{.}}\left(t\right)\\ C^{\text{rot}\text{.}}\left(t\right)=\frac{1}{5}{e}^{-t/{\tau }_{\text{M}}}\\ C^{\text{intern}\text{.}}\left(t\right)=S^2+\left(1-S^2\right)G\left(t\right)\\ G\left(0\right)=1,\underset{t\to \infty }{\lim }G\left(t\right)=0\end{array}(3)$

The first equation is the result of statistical independence of internal and overall motion, such that we may write the total correlation function, $C\left(t\right)$, as a product of a correlation function resulting from the internal motion ($C^{\text{intern}\text{.}}\left(t\right)$), and a correlation function resulting from the overall rotational tumbling ($C^{\text{rot}\text{.}}\left(t\right)$). The overall motion may be described by a single decaying exponential, with correlation time ${\tau }_{\text{M}}$ if that overall motion is isotropic (occurring if the molecule is approximately spherical). For internal motion, $C^{\text{intern}\text{.}}\left(t\right)$ has an initial value of 1, and equilibrates at $S^2$. $S^2$ is referred to as the generalized order parameter, and is related to, but not always equal to order parameters that may be extracted from measurement of residual couplings, as will be discussed in Determining S². $G\left(t\right)$ is simply the decaying part of $C^{\text{intern}\text{.}}\left(t\right)$, normalized such that its initial value is 1, and final value is 0. If the second assumption, fast decay of the correlation function due to internal motion is fulfilled, we may calculate $J\left(\omega \right)$ using the parameters ${\tau }_{\text{M}}$, $S^2$, and $\langle {\tau }_e\rangle$, where

$\langle {\tau }_e\rangle ={\displaystyle \underset{0}{\overset{\infty }{\int }}{e}^{-t/{\tau }_{\text{M}}}G\left(t\right)dt}(4)$

We calculate $J\left(\omega \right)$ in order to see how it is a function of the parameters ${\tau }_{\text{M}}$, $S^2$, and $\langle {\tau }_e\rangle$.

$\begin{array}{l}J\left(\omega \right)=\frac{2}{5}{\displaystyle \underset{0}{\overset{\infty }{\int }}\left[S^2{e}^{-t/{\tau }_{\text{M}}}+\left(1-S^2\right)e^{-t/{\tau }_{\text{M}}}G\left(t\right)\right]\cos \left(\omega t\right)dt}\\ =\frac{2}{5}\left[\frac{S^2{\tau }_{\text{M}}}{1+{\left(\omega {\tau }_{\text{M}}\right)}^2}+\left(1-S^2\right){\displaystyle \underset{0}{\overset{\infty }{\int }}{e}^{-t/{\tau }_{\text{M}}}G\left(t\right)\underset{\approx 1}{\underbrace{\cos \left(\omega t\right)}}dt}\right]\\ =\frac{2}{5}\left[\frac{S^2{\tau }_{\text{M}}}{1+{\left(\omega {\tau }_{\text{M}}\right)}^2}+\left(1-S^2\right)\langle {\tau }_e\rangle \right]\end{array}(5)$

We see that if $e^{-t/{\tau }_{\text{M}}}G\left(t\right)$ decays quickly compared to $\omega$, then we may replace $\cos \left(\omega t\right)$ with 1, since the exponential approaches zero more quickly than the cosine term can evolve away from 1. Then, regardless of the precise form of $G\left(t\right)$, $J\left(\omega \right)$ may always be calculated from the parameters $S^2$, $\langle {\tau }_e\rangle$, and ${\tau }_{\text{M}}$. Furthermore, if ${\tau }_{\text{M}}$ is known (usually from the analysis of R₁ and R₂ throughout a molecule (Kay et al., 1989)), $J\left(\omega \right)$ becomes a linear function of the parameters $S^2$ and $\left(1-S^2\right)\langle {\tau }_e\rangle$.

Instead of assuming fast decay of $G\left(t\right)$, one may alternatively assume that it is mono-exponential ($G\left(t\right)=e^{-t/\tau }$), yielding

$\begin{array}{l}J\left(\omega \right)=\frac{2}{5}{\displaystyle \underset{0}{\overset{\infty }{\int }}\left[S^2{e}^{-t/{\tau }_{\text{M}}}+\left(1-S^2\right)e^{-t/{\tau }_{\text{M}}}{e}^{-t/\tau }\right]\cos \left(\omega t\right)dt}\\ {\langle {\tau }_e\rangle }^{-1}={\tau }_{\text{M}}^{-1}+{\tau }^{-1}\\ J\left(\omega \right)=\frac{2}{5}\left[\frac{S^2{\tau }_{\text{M}}}{1+{\left(\omega {\tau }_{\text{M}}\right)}^2}+\left(1-S^2\right)\frac{\langle {\tau }_e\rangle }{1+{\left(\omega \langle {\tau }_e\rangle \right)}^2}\right]\end{array}(6)$

In the extreme narrowing limit, where decay of the correlation function is fast, we have $\omega \langle {\tau }_e\rangle \ll 1$ such that this result equals the result in Eq. 5. The expression in Eq. 6 is equivalent to Eq. 1 in (Lipari and Szabo, 1982a), and is valid either in the case of mono-exponential decay or fast decay of the internal correlation function. However, we find the case of fast, multi-exponential decay the more likely scenario, and so focus on this assumption.

The notation $\langle {\tau }_e\rangle$ is used to indicate the average of the effective correlation time. To understand how the integral of $e^{-t/{\tau }_{\text{M}}}G\left(t\right)$ is related to this average, we must assume that $G\left(t\right)$ is the sum of decaying exponentials. This may be achieved with a sum over a discrete number of correlation times, weighted with $A_i$, or a continuous distribution, defined by the function $\theta \left(z\right)$.

$\begin{array}{l}G\left(t\right)={\displaystyle \sum _i{A}_i{e}^{-t/{\tau }_i}}\\ \text{where}{\Sigma }_i{A}_i=1\\ \text{–or–}\\ G\left(t\right)={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)e^{-t/\left({10}^z\cdot 1\text{ s}\right)}dz}\\ \text{where}{\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)dz}=1\end{array}(7)$

Since $G\left(0\right)=1$, it is clear that the sum of amplitudes ($A_i$) must be 1. For the former equation, we take a simple sum, and for the latter form, we use a distribution of correlation times, $\theta \left(z\right)$, given on a logarithmic scale, such that $z={\log }_{10}\left({\tau }_c/\text{s}\right)$. The distribution must similarly integrate to 1. The two forms can be treated equivalently. We have recently re-introduced the latter form (Smith et al., 2018), which was previously used to describe a variety of continuous correlation time distributions, e.g., see Beckmann (1988). We may insert this expression for $G\left(t\right)$ (Eq. 7) into Eq. 4 in order to obtain the relationship between $\theta \left(z\right)$ and $\langle {\tau }_e\rangle$.

$\begin{array}{l}\langle {\tau }_e\rangle ={\displaystyle \underset{0}{\overset{\infty }{\int }}{e}^{-t/{\tau }_{\text{M}}}\underset{G\left(t\right)}{\underbrace{{\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)e^{-t/\left({10}^z\cdot 1\text{ s}\right)}dz}}}}dt\\ ={\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)e^{-t\left({\tau }_{\text{M}}^{-1}+{\left({10}^z\cdot 1\text{s}\right)}^{-1}\right)}dz}dt}\\ {\left({\tau }_e\left(z\right)\right)}^{-1}={\tau }_{\text{M}}^{-1}+{\left({10}^z\cdot 1\text{s}\right)}^{-1}\\ \langle {\tau }_e\rangle ={\displaystyle \underset{0}{\overset{\infty }{\int }}{\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)e^{-t/{\tau }_e\left(z\right)}dz}dt}={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right){\left(-{\tau }_e\left(z\right)e^{-t/{\tau }_e\left(z\right)}\right)|}_{t=0}^{\infty }dz}\\ \langle {\tau }_e\rangle ={\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right){\tau }_e\left(z\right)dz}\\ \text{equivalently}:\langle {\tau }_e\rangle ={\displaystyle \sum _i{A}_i{\tau }_e^i},\text{for}{\left({\tau }_e^i\right)}^{-1}={\tau }_{\text{M}}^{-1}+{\tau }_i^{-1}\end{array}(8)$

${\left({\tau }_e^i\right)}^{-1}={\tau }_{\text{M}}^{-1}+{\left({10}^z\cdot 1\text{s}\right)}^{-1}$ is the effective correlation time, resulting from decay of both the correlation function due to the internal correlation time, $z={\log }_{10}\left({\tau }_c/s\right)$ and correlation time of the overall motion, ${\tau }_{\text{M}}$. Since $\theta \left(z\right)$ integrates to 1, ${\displaystyle {\int }_{-\infty }^{\infty }\theta \left(z\right){\tau }_e\left(z\right)dz}$ yields the weighted average of the effective correlation time, $\langle {\tau }_e\rangle$. Then, one fits experimental data to a correlation function having the following model:

$C^{}\left(t\right)=\frac{1}{5}(S^2{e}^{-t/{\tau }_M}+\left(1-S^2\right)e^{-t/\langle {\tau }_e\rangle })(9)$

Applying this model does not require that the true correlation function has exactly this form, but rather, the model correlation function simply must have the same values of $S^2$ and $\langle {\tau }_e\rangle$ as the true correlation function. In this sense, the analysis itself remains model-free, although equating $\langle {\tau }_e\rangle$ with the averaged effective correlation time requires the true correlation function to be a sum of decaying exponentials, as in Eq. 7.

A Few Notes on Linearity

We will later note that many of the methods used for analyzing relaxation rate constants result in parameters that are linear functions of the distribution of correlation times, $\left(1-S^2\right)\theta \left(z\right)$. Specifically, we mean that any parameter, $P_m$, is linear to $\left(1-S^2\right)\theta \left(z\right)$ if it can be written as

$P_m=\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)p_m\left(z\right)dz}(10)$

That is, for every correlation time, $z={\log }_{10}\left({\tau }_c/\text{s}\right)$, P increases proportionally to $\left(1-S^2\right)\theta \left(z\right)$ at that correlation time, where the proportionality is defined by $p_m\left(z\right)$. Furthermore, any linear combination of parameters, $P_m$, is then also linear to $\left(1-S^2\right)\theta \left(z\right)$, as we can see by integrating a sum of parameters, $P_m$, and swapping the order of the integration and the summation.

$\begin{array}{l}{\displaystyle \sum _m{a}_m{P}_m}={\displaystyle \sum _m{a}_m\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}{p}_m\left(z\right)\theta \left(z\right)dz}}\\ =\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\underset{=\Sigma \left(z\right)}{\underbrace{\left[{\displaystyle \sum _m{a}_m{p}_m\left(z\right)}\right]}}\theta \left(z\right)dz}\\ =\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\Sigma \left(z\right)\theta \left(z\right)dz}\end{array}(11)$

We define the function $\Sigma \left(z\right)$ to be the weighted sum of the sensitivities, $p_m\left(z\right)$, which then defines the linear relationship of the sum of the $P_m$ to $\left(1-S^2\right)\theta \left(z\right)$. This principle is one of the basic tenants of linear algebra. What can be less obvious is that a linear fit of parameters, $P_m$, defined by a matrix, M, to a new set of parameters, $Q_n$ is also linear to $\left(1-S^2\right)\theta \left(z\right)$. This is only the case if restrictions on the fit parameters, $Q_n$, are not applied (no priors are used). In this case, the parameters $Q_n$ should minimize the following equation.

$\begin{array}{l}\min \left[{\displaystyle \sum _m{\left|P_m-{\left[\mathbf{M}\right]}_{m,n}{Q}_n\right|}^2}\right]\\ Q_n={\displaystyle \sum _m{\left[{\mathbf{M}}^{-1}\right]}_{n,m}{P}_m}\end{array}(12)$

One may determine the $Q_n$ by computing the pseudoinverse of M (denoted ${\mathbf{M}}^{-1}$) and multiplying by the $P_m$. Linearity of the $Q_n$ to $\left(1-S^2\right)\theta \left(z\right)$ results from the fact that linear combinations defined by ${\mathbf{M}}^{-1}$ remain unchanged regardless of the value of the parameters being fit, $P_m$. However, if the allowed values of the $Q_n$ are restricted with priors, then it can be that some values of $P_m$ will result in the latter formula in Eq. 12 yielding $Q_n$ outside of the allowed range. In this case, a linear least squares algorithm will search for a different solution than that given by Eq. 12, such that the $Q_n$ are no longer defined by ${\mathbf{M}}^{-1}$, and no longer have a consistent linear relationship to $\left(1-S^2\right)\theta \left(z\right)$. Note that if priors are used, but Eq. 12 does not yield $Q_n$ outside of the bounds defined by the priors, then Eq. 12 still remains the best solution and linearity is maintained. In general, we will find analysis methods that rely on linear combination of data have more predictable behavior than those that do not.

Then, the model-free parameters $S^2$ and $\left(1-S^2\right)\langle {\tau }_e\rangle$ are linear to $\left(1-S^2\right)\theta \left(z\right)$, because one can fit experimental relaxation rate constants with $S^2$ and $\left(1-S^2\right)\langle {\tau }_e\rangle$ (see Eq. 5), where the relaxation rate constants themselves are linear to the spectral density (Eq. 1), the spectral density is linear to the correlation function (via Fourier transform, Eq. 2), and the correlation function is linear to the distribution of correlation times, $\left(1-S^2\right)\theta \left(z\right)$ (Eqs. 3, 7)). Assuming the correlation function decays quickly, this linear relationship is given by the following, where ${\tau }_e\left(z\right)$ is defined in Eq. 8.

$\begin{array}{l}{S}^2=1-\left[\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)dz}\right]\\ \left(1-S^2\right)\langle {\tau }_e\rangle =\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}{\tau }_e\left(z\right)\theta \left(z\right)dz}\end{array}(13)$

Note that $\langle {\tau }_e\rangle$ is not itself linear to $\left(1-S^2\right)\theta \left(z\right)$, but is easily obtained from the above parameters.

Fitting With Model-Free

In Figure 3, we test the performance of model-free fitting under a number of conditions. In Figure 3A, we calculate a number of relaxation rate constants from motion having a single internal correlation time and overall tumbling with ${\tau }_{\text{M}}$ = 4 ns, and then fit the results, assuming the model-free correlation function (Eq. 9). We may calculate the spectral density exactly, or we may assume that the correlation function decays quickly, by using the spectral density given in Eq. 5, resulting in a linear fit. The former method is shown as a blue, solid line, where the input parameters always exactly match the fit parameters, whereas using a linear fit (red, dashed line) results in disagreement of input and fit parameters when the correlation function does not decay quickly compared to the frequencies sampled ($\omega {\tau }_e\ll 1$); in this case, Eq. 5 is no longer a good estimate of the spectral density whereas Eq. 6 has the correct form.

FIGURE 3

FIGURE 3. Model-free fit parameters as a function of input parameters. For each plot, a data set is calculated, using the experiments found in from Table I of Lipari and Szabo (1982b), and the resulting rate constants are fit using the model-free approach, with the resulting ${\langle {\tau }_e\rangle }_{\text{fit}}$ and $\left(1-S^2\right)$ shown on the left and right, respectively. For all plots, the tumbling correlation time is ${\tau }_{\text{M}}$ = 4 ns and $\left(1-S^2\right)$ = 0.3. One correlation time of the internal motion is varied, and we plot ${\langle {\tau }_e\rangle }_{\text{in}}$ on the x-axis. In each plot, we fit using the full spectral density (blue, solid, see Eq. 6) and using a linear approximation (red, dashed, see Eq. 5). In (A), the input correlation function only has a single correlation time. In (B), one correlation time is fixed to 10 ps, and the second correlation time is swept. In (C), a log-Gaussian distribution (μ = 10 ps, σ = 0.75 order of magnitude) is combined with a correlation time that is varied (with total amplitude equal). On the left plots, black dotted lines indicate where the input value, ${\langle {\tau }_e\rangle }_{\text{in}}$, matches the fit, ${\langle {\tau }_e\rangle }_{\text{fit}}$. In all plots, vertical black dotted lines indicate where $\omega {\langle {\tau }_e\rangle }_{\text{in}}=0.5$ for $\omega /2\pi$ = 90 MHz, where this frequency corresponds to the highest field used for the data set.

In Figure 3B, we include two correlation times in the input, each with equal amplitude, where one correlation time is fixed (10 ps), and a second correlation time is swept. We calculate the mean effective correlation time directly on the x-axis (${\langle {\tau }_e\rangle }_{\text{in}}$), and compare this to the fitted parameters on the y-axis (${\langle {\tau }_e\rangle }_{\text{fit}}$, left plot, $1-S^2$, right plot). As expected, if the assumption that $\omega {\tau }_e\ll 1$ holds for all frequencies sampled and all correlation times present, the fit parameters are in good agreement with their input values, but when $\omega {\tau }_e\ll 1$, ${\langle {\tau }_e\rangle }_{\text{fit}}$ and $S^2$ no longer reproduce the correct values. Note that performing this fit with the full spectral density (blue, solid line) and using just a linear fit (red, dashed line) produces very similar results. In Figure 3C, we perform the same tests, but instead of fixing a correlation time to 10 ps, we have a log-Gaussian distribution of correlation times, centered at 10 ps, with a standard deviation of 0.75 orders of magnitude. Results are similar to those found in Figure 3B.

Determining S²

For model-free analysis, $\langle {\tau }_e\rangle$ is the average effective correlation time, and can be calculated from the distribution of correlation times. $S^2$, on the other hand, is determined from the distribution of orientations sampled by internal motion. By definition, it is equal to the correlation function of internal motion, taken as the limit of t goes to infinity. We may obtain $S^2$ by first considering the formula for the correlation function.

$C^{\text{intern}\text{.}}\left(t\right)={\langle P_2\left(\overrightarrow{\mu }\left(\tau \right)\cdot \overrightarrow{\mu }\left(t+\tau \right)\right)\rangle }_{\tau }(14)$

$P_2\left(x\right)$ is the second Legendre polynomial ($P_2\left(x\right)=\left(3x^2-1\right)/2$), and $\overrightarrow{\mu }\left(\tau \right)$ is a normalized vector that gives the direction of the principal component of an NMR interaction as a function of time, due to internal motion only (without tumbling). The dot product ($\overrightarrow{\mu }\left(\tau \right)\cdot \overrightarrow{\mu }\left(t+\tau \right)$) yields the cosine of the angle between the two vectors. The correlation function itself may take on a variety of complex forms, depending on the correlation times present, but $S^2$, its value as $t\to \infty$, depends only on the distribution of orientations sampled by the internal motion. This may be obtained by taking a weighted average over all possible starting orientations (p) and all possible final orientations (q), and calculating $P_2\left({\overrightarrow{\mu }}_p\cdot {\overrightarrow{\mu }}_q\right)$ for each pair. Defining $p_{\text{eq}}\left({\overrightarrow{\mu }}_p\right)$ to be the fraction of orientation ${\overrightarrow{\mu }}_p$ at thermal equilibrium, we obtain

$S^2={\displaystyle \sum _p{\displaystyle \sum _q{p}_{\text{eq}}\left({\overrightarrow{\mu }}_p\right)p_{\text{eq}}\left({\overrightarrow{\mu }}_q\right)P_2\left({\overrightarrow{\mu }}_p\cdot {\overrightarrow{\mu }}_q\right)}}(15)$

Then, if we have a precise description of the internal dynamics, we may calculate parameters $\langle {\tau }_e\rangle$ and $S^2$ using Eqs. 8, 15. We may not easily go backwards, to obtain a precise description of the dynamics from only these parameters. However, this is not a limitation of the method of analysis, but rather of the information content of the data.

In solid-state NMR, we no longer have overall tumbling motion, so the term $e^{-t/{\tau }_{\text{M}}}$ vanishes from the correlation function and Eq. 5 becomes simply

$J\left(\omega \right)=\frac{2}{5}\left(1-S^2\right)\langle \tau \rangle (16)$

This prevents us from separating $S^2$ and $\langle \tau \rangle$ via relaxation data alone (we drop the subscript e from $\tau$, since it is no longer an effective correlation time); however, one may measure the size of residual couplings in NMR (Chevelkov et al., 2009a; Schanda et al., 2011), often via DIPSHIFT (Munowitz et al., 1981) or REDOR (Gullion and Schaefer, 1989). In this case, the ratio of the anisotropies of the rigid interaction (${\delta }_{\text{rigid}\text{.}}$) to the motionally averaged interaction (${\delta }_{\text{resid}\text{.}}$) defines $S_{\text{resid}\text{.}}$.

$S_{\text{resid}\text{.}}={\delta }_{\text{resid}\text{.}}/{\delta }_{\text{rigid}}(17)$

One usually equates $S^2$ and $S_{\text{resid.}}^2$, although for motion that does not have at least a three-fold symmetry axis, these terms are not necessarily equal (Supplementary Section S3). Examples are found in Figures 2C-E, although we see the deviation is actually quite small (e.g., $S_{\text{resid}\text{.}}^2$ = 0.69, vs. $S^2$ = 0.7), so that this approach may be used to obtain good separation of $S^2$ and $\langle \tau \rangle$.

Alternative Methods

In the case that all internal motion is fast, such that the correlation function decays quickly, model-free analysis is an ideal approach for extracting dynamics information from relaxation data: the full information content of the relaxation data is captured in the parameters $S^2$ and $\langle {\tau }_e\rangle$, where these parameters have simple relationships to the distribution of correlation times, $\left(1-S^2\right)\theta \left(z\right)$ (parameters $S^2$ and $\left(1-S^2\right)\langle {\tau }_e\rangle$ are furthermore linearly related to $\left(1-S^2\right)\theta \left(z\right)$). In case the correlation function does not decay quickly compared to the sampled frequencies, our formula for the spectral density becomes significantly more complex. To obtain it, we begin from Eq. 5 (first expression), and insert the assumed form of $G\left(t\right)$, found in Eq. 7, yielding the equation for the solution-state spectral density.

$\begin{array}{l}J\left(\omega \right)=\frac{2}{5}{\displaystyle \underset{0}{\overset{\infty }{\int }}\left[S^2{e}^{-t/{\tau }_{\text{M}}}+\left(1-S^2\right)e^{-t/{\tau }_{\text{M}}}{\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)e^{-t/\left({10}^z\cdot 1\text{s}\right)}dz}\right]\cos \left(\omega t\right)dt}\\ {\left({\tau }_e\left(z\right)\right)}^{-1}={\tau }_{\text{M}}^{-1}+{\left({10}^z\cdot 1\text{s}\right)}^{-1}\\ z_e\left(z\right)={\log }_{10}\left({\tau }_e\left(z\right)/\text{s}\right)\\ J\left(\omega \right)=\frac{2}{5}{\displaystyle \underset{0}{\overset{\infty }{\int }}\left[S^2{e}^{-t/{\tau }_{\text{M}}}+\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)e^{-t/\left({10}^{z_e(z)}\cdot 1\text{s}\right)}dz}\right]\cos \left(\omega t\right)dt}\\ =\frac{2}{5}\left[\frac{S^2{\tau }_{\text{M}}}{1+{\left(\omega {\tau }_{\text{M}}\right)}^2}+\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)\frac{{10}^{z_e\left(z\right)}\cdot 1\text{s}}{1+{\left(\omega \cdot {10}^{z_e\left(z\right)}\cdot 1\text{s}\right)}^2}dz}\right]\end{array}(18)$

The first step is to combine the two exponential terms, where we define the log-effective correlation time, $z_e\left(z\right)$, as a function of the log-internal correlation time, z, and also the rotational correlation time, ${\tau }_{\text{M}}$. Subsequently, each exponential term is Fourier transformed to yield the familiar Lorentzian function. The spectral density for solid-state NMR can be similarly calculated, where the overall motion is omitted.

$\begin{array}{l}J\left(\omega \right)=\frac{2}{5}{\displaystyle \underset{0}{\overset{\infty }{\int }}\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)e^{-t/\left({10}^z\cdot 1\text{s}\right)}dz}\cos \left(\omega t\right)dt}\\ =\frac{2}{5}\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)\frac{{10}^z\cdot 1\text{s}}{1+{\left(\omega \cdot {10}^z\cdot 1\text{s}\right)}^2}dz}\end{array}(19)$

The integral has a complex dependence on $\omega$, and depends on the specific form of $\left(1-S^2\right)\theta \left(z\right)$, so that by using multiple relaxation experiments, we can extract more than two parameters describing the internal motion. However, we require a different approach to extract that information. We discuss four approaches developed for treating this case: the extended model-free approach (EMF), spectral density mapping (SDM), LeMaster’s approach, and IMPACT. Another approach that bears mentioning is the slowly relaxing local structure model (SRLS), which accounts for coupling of local motional modes to overall motion of a molecule in solution (Polimeno and Freed, 1992; Tugarinov et al., 2001; Mendelman and Meirovitch, 2021; Shapiro and Meirovitch, 2012). SRLS reduces to the model-free approach as coupling between local and overall motion vanishes. However, we do not include further comparison to the analytically simpler methods discussed here.

Extended Model-Free

Clore and coworkers found that when measuring relaxation data at higher fields (up to 600 MHz) that not all backbone motion could be well fit using the model-free approach for staphylococcal nuclease and interleukin-1β (Clore et al., 1990). They found that the simplest correlation function that could fit the data was obtained by adding another decaying exponential term, yielding the EMF correlation function.

$C^{\text{intern}\text{.}}\left(t\right)=\left(1-S_f^2\right)e^{-t/{\tau }_f}+S_f^2\left(1-S_s^2\right)e^{-t/{\tau }_s}+S_f^2{S}_s^2(20)$

In this correlation function, the total internal motion is separated into fast and slow components, with order parameters $S_f^2$ and $S_s^2$, and effective correlation times, ${\tau }_f$ and ${\tau }_s$, respectively. The product $S_f^2{S}_s^2$ should yield the total order parameter, $S^2$. Also note that the faster motion’s order parameter scales the influence of the slower motion, as seen in the term $S_f^2\left(1-S_s^2\right)e^{-t/{\tau }_s}$. Data analysis with EMF in solid- and solution-state NMR involves simply varying the parameters, $S_f^2$, $S_s^2$, ${\tau }_f$, and ${\tau }_s$, to find an optimal fit to experimental data. Often, one also performs a model selection step, where one may determine how many parameters should be included in the fit (Mandel et al., 1995; d’Auvergne and Gooley, 2003; Zinkevich et al., 2013; Gill et al., 2016). In Figure 4, the behavior of EMF parameters is shown for several correlation functions. In each subplot, all terms except one correlation time are fixed, and we observe the model behavior as we sweep through the variable correlation time. In Figure 4A, two correlation times are used, so that the input correlation function has the same form as the correlation function used for fitting; as expected, the fitted parameters perfectly match the input parameters, since the input and fit models match. In Figure 4B, three correlation times are input, where the fast and slow correlation times are fixed at 10 ps and 1 ns, and the intermediate correlation time is swept. In this case, when the intermediate correlation time is fast, the fitted ${\tau }_f$ falls in between the fast and intermediate correlation times, and the fitted amplitude for the fast motion is the sum of the input amplitudes for the fast and intermediate motions. However, for longer correlation times, the fitted ${\tau }_f$ again gets shorter, eventually equaling 10 ps, so that the fitted ${\tau }_s$ takes over the role of fitting the intermediate correlation time. This is especially well illustrated in Figure 4(B, right), where the amplitude corresponding to the slow motion increases from 0.1 to 0.2, indicating that the slow motion in the model fits both the input intermediate and slow motions. Similar behavior is observed in Figure 4C, where a distribution of correlation times is combined with a single correlation time that is swept.

FIGURE 4

FIGURE 4. EMF parameters as a function of input correlation time (solution-state). For each plot, a data set is calculated, using the set of experiments from Clore et al. (1990), and the resulting rate constants are fitted using the EMF approach. For all plots, ${\tau }_{\text{M}}$ = 8.3 ns, and the input $\left(1-S^2\right)$ = 0.3. In each subplot, the fitted correlation times (left) and amplitudes (right) are shown, as a function of an input correlation time (x-axis). In (A), the input correlation function has two correlation times (with equal amplitudes), with one fixed at 10 ps, and the other swept. In (B), the input correlation function has three correlation times, two fixed at 10 ps and 1 ns, and the third is swept. In (C), a log-Gaussian distribution of correlation times is used (μ = 100 ps, σ = 0.75 orders of magnitude), and a single correlation time is swept. Black dotted lines show the input correlation times (left plots).

To the best of our knowledge, the behavior of the fit parameters has no well-defined relationship to the distribution of correlation times, $\left(1-S^2\right)\theta \left(z\right)$: if we know $\left(1-S^2\right)\theta \left(z\right)$ precisely, our only way to obtain the EMF parameters from it would be to explicitly calculate a set of relaxation rate constants, and then fit the results to Eq. 20. This is in sharp contrast to the original model-free parameters. Similar limitations arise for the EMF approach in solid-state NMR, as seen in Figure 5. Note that typical solution-state data sets are fairly continuous in their sensitivity to motion as a function of correlation time (Smith et al., 2019a), whereas solid-state NMR has a “blind-spot” in sensitivity centered around ∼100 ns (Schanda, 2019), which results in some of the more unusual behavior for EMF in solids (see Case 1: Extended Model-Free for a detailed discussion of the behavior of typical model-free parameters in solid-state NMR).

FIGURE 5

FIGURE 5. EMF parameters as a function of input correlation time (solid-state). For each plot, a data set is calculated, including direct measurement of S_resid. via residual couplings (Eq. 17), ¹⁵N T₁ at 400, 500, and 850 MHz, and T₂ with MAS of 60 kHz. The resulting rate constants are fitted using the EMF approach. For all plots, $\left(1-S^2\right)$ = 0.3. In each subplot, the fitted correlation times (left) and amplitudes (right) are shown, as a function of an input correlation time (x-axis). In (A), the input correlation function has two correlation times (with equal amplitudes), with one fixed at 3.2 ps, and the other swept. In (B), the input correlation function has three correlation times, two fixed at 3.2 ps and 32 ns, and the third is swept. In (C), a log-Gaussian distribution of correlation times is used (μ = 100 ps, σ = 0.75 orders of magnitude), and a single correlation time is swept. Black dotted lines show the input correlation times (left plots).

Spectral Density Mapping

In contrast to EMF, SDM is achieved by simple linear combination of sets of relaxation data at a single magnetic field (Peng and Wagner, 1992; Ishima et al., 1999). From a set of R₁, R₂, and NOE relaxation rate constants, one calculates

$\begin{array}{l}J\left(0\right)=\frac{{\mathit{R}}_2-{\mathit{R}}_1/2-0.454{\sigma }_{\text{IS}}}{{\delta }_{\text{IS}}^2/2+2{\left(\Delta {\sigma }_{\text{I}}{\omega }_{\text{I}}\right)}^2}\\ J\left({\omega }_{\text{I}}\right)=\frac{{\mathit{R}}_1-1.249{\sigma }_{\text{IS}}}{3{\left({\delta }_{\text{IS}}/4\right)}^2+{\left(\Delta {\sigma }_{\text{I}}{\omega }_{\text{I}}\right)}^2/3}\\ J\left(0.870{\omega }_{\text{S}}\right)=16{\sigma }_{\text{IS}}/\left(5{\delta }_{\text{IS}}^2\right)\end{array}(21)$

The above expressions yield very close approximations of the spectral density at specific frequencies: 0, ${\omega }_{\text{I}}$, and 0.870${\omega }_{\text{S}}$, where ${\omega }_{\text{I}}$ is the nuclear Larmor frequency of the spin being relaxed, and ${\omega }_{\text{S}}$ is a spin which is dipole coupled to that spin (usually a directly bonded ¹H). Differences in the representations of the anisotropies (${\delta }_{\text{IS}}$, $\Delta {\sigma }_{\text{I}}{\omega }_{\text{I}}$) result in the different appearances of the normalization factors (denominators). These terms may be interpreted as being proportional to the amount of motion near the given frequency (which corresponds to the correlation time $\tau =1/\omega$), but otherwise they do not provide a more physical interpretation of the motion. One may subsequently fit the spectral densities to model-free parameters for better interpretation (Gill et al., 2016). If we have a precise description of the motion (e.g., $\left(1-S^2\right)\theta \left(z\right)$), the terms $J\left(\omega \right)$ are easily obtained:

$J\left(\omega \right)=\frac{2}{5}\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)\frac{{10}^z\cdot 1\text{s}}{1+{\left(\omega \cdot {10}^z\cdot 1\text{s}\right)}^2}dz}(22)$

The parameters resulting from SDM always behave the same way in response to a given correlation time, regardless of other correlation times present, and is the consequence of properties of linearity discussed in A Few Notes on Linearity. This is seen in Figure 6A, where we calculate relaxation rate constants resulting from a single correlation time and analyze with SDM. In Figure 6B, we split motion over two correlation times, and observe how the terms respond to sweeping one of them, and in Figure 6C, we split motion into a distribution and a single, swept correlation time and determine how the terms respond to the swept correlation time. The result is always identical (scaling by 0.5 results from dividing the total amplitude into two parts), a very useful property occurring when data is analyzed strictly by linear combination of data. Unlike EMF analysis, behavior of SDM is independent of the form of the distribution of correlation times.

FIGURE 6

FIGURE 6. Behavior of SDM as a function of correlation time. In each subplot, we calculate ¹⁵N T₁, T₂, and σ_NH at 600 MHz, and analyze the results using Eq. 21. In (A), the input total correlation function consists of a single decaying exponential term (with amplitude 1), where the terms $J\left(\omega \right)$ are plotted as the correlation time is varied (results are normalized). Black dotted lines show the spectral densities, $J\left(0\right)$, $J\left({\omega }_{\text{I}}\right)$, $J\left(0.870{\omega }_{\text{S}}\right)$, calculated with Eq. 22, and colored lines show the results of the data analysis, yielding an almost exact correspondence. In (B), the total correlation function now uses two correlation times (equal amplitudes), with one fixed at 10 ps, and the second swept (x-axis). On the y-axis, we plot contribution to the terms, $\Delta J\left(\omega \right)$, from the correlation time being varied. The resulting behavior is identical to that in (A), except that the amplitude is half as large, since we have split the total amplitude between the fixed and variable correlation time (dashed line marks 0.5). In (C), the same information is plotted, but the total correlation function includes a log-Gaussian distribution (μ = 630 ps, σ = 1 order of magnitude), and a single, variable correlation time.

Note that this approach describes the total motion, and does not separate out tumbling from internal motion in the case of solution-state NMR, which has an especially strong influence on $J\left(0\right)$. The original approach only incorporates data from one field, whereas later work has extended the method to include data from more than one field, although one still requires specific sets of experiments (Skrynnikov et al., 2002; Hsu et al., 2018).

LeMaster’s Approach

LeMaster proposed an alternative to SDM analysis of R₁, R₂, and NOE data from one field, in order to separate overall tumbling from internal motion (LeMaster, 1995). In this case, LeMaster proposed fitting data to the following correlation function:

$\begin{array}{l}C\left(t\right)=S_f^2{S}_{\text{H}}^2{S}_{\text{N}}^2{e}^{-t/{\tau }_{\text{M}}}+S_f^2\left(1-S_{\text{H}}^2\right)e^{-t/{\tau }_{\text{H}}}+S_f^2{S}_{\text{H}}^2\left(1-S_{\text{N}}^2\right)e^{-t/{\tau }_{\text{N}}}+{\left(1-S_f\right)}^2{e}^{-t/{\tau }_f}\\ {\tau }_{\text{H}}={\left({\omega }_{\text{H}}+{\omega }_{\text{N}}\right)}^{-1},{\tau }_{\text{N}}=|{\omega }_{\text{N}}{|}^{-1}\end{array}(23)$

It is assumed that ${\tau }_f$ is very short so that the term $\left(1-S_f^2\right)e^{-t/{\tau }_f}$ makes only negligible contributions to the spectral density, resulting in the following formula:

$\begin{array}{l}J\left(\omega \right)=\frac{2}{5}{S}_f^2\left[S_{\text{H}}^2{S}_{\text{N}}^2\frac{{\tau }_{\text{M}}}{1+{\left(\omega {\tau }_{\text{M}}\right)}^2}+\left(1-S_{\text{H}}^2\right)\frac{{\tau }_{\text{H}}}{1+\left(\omega {\tau }_{\text{H}}^2\right)}+S_{\text{H}}^2\left(1-S_{\text{N}}^2\right)\frac{{\tau }_{\text{N}}}{1+\left(\omega {\tau }_{\text{N}}^2\right)}\right]\\ =\frac{2}{5}[\frac{{\tau }_{\text{M}}}{1+{\left(\omega {\tau }_{\text{M}}\right)}^2}+\left(1-S_f^2\right)\left(-\frac{{\tau }_{\text{M}}}{1+{\left(\omega {\tau }_{\text{M}}\right)}^2}\right)+\\ S_f^2\left(1-S_{\text{H}}^2\right)\left(\frac{{\tau }_{\text{H}}}{1+\left(\omega {\tau }_{\text{H}}^2\right)}-\frac{{\tau }_{\text{M}}}{1+{\left(\omega {\tau }_{\text{M}}\right)}^2}\right)+S_f^2{S}_{\text{H}}^2\left(1-S_{\text{N}}^2\right)\left(\frac{{\tau }_{\text{N}}}{1+\left(\omega {\tau }_{\text{N}}^2\right)}-\frac{{\tau }_{\text{M}}}{1+{\left(\omega {\tau }_{\text{M}}\right)}^2}\right)]\end{array}(24)$

In the latter formulation, we find that the spectral density becomes a linear combination of terms, weighted by $\left(1-S_f^2\right)$, $S_f^2\left(1-S_{\text{H}}^2\right)$, and $S_f^2{S}_{\text{H}}^2\left(1-S_{\text{N}}^2\right)$. Then, one must fit these terms to the experimental relaxation rate constants. We do so in Figure 7 for calculated relaxation rate constants. Like SDM, responses as a function of correlation time are always identical (again, excepting a scaling factor of 0.5 resulting from splitting the total motion into components), although the functions themselves are different: this results from the fact that LeMaster’s approach characterizes the internal motion, and not the total motion, so that we obtain one amplitude, $\left(1-S_f^2\right)$, which captures information about the fastest correlation times (<30 ps), one amplitude, $S_f^2\left(1-S_{\text{H}}^2\right)$, which captures information for correlation times near to ${\tau }_{\text{H}}$, and one amplitude, $S_f^2{S}_{\text{H}}^2\left(1-S_{\text{N}}^2\right)$, which captures information for correlation times near to ${\tau }_{\text{N}}$.

FIGURE 7

FIGURE 7. Behavior of LeMaster’s approach as a function of correlation time. In each subplot, we calculate ¹⁵N T₁, T₂, and σ_NH at 600 MHz for motion with $\left(1-S^2\right)$ = 0.3 and tumbling correlation time of ${\tau }_{\text{M}}$ = 4 ns, and analyze the results using Eq. 24. In (A) the internal correlation function consists of a single decaying exponential term (with amplitude 0.3), where the fitted amplitudes are plotted as the correlation time is varied. In (B) the internal correlation function uses two correlation times (both amplitudes are 0.15), with one correlation time fixed at 10 ps, and the second swept (x-axis). On the y-axis, we plot contributions to the terms from the correlation time being varied. The resulting behavior is identical to that in (A), except that the amplitude is half, since we have split the total amplitude between the fixed and variable correlation time (dashed line marks 0.15). In (C), the same information is plotted, but the total correlation function includes a log-Gaussian distribution (μ = 630 ps, σ = 1 order of magnitude), and a single, variable correlation time.

LeMaster’s approach is a linear fit, without priors; as discussed in A Few Notes on Linearity, this means that the fitted parameters may also be obtained by a linear combination of the experimental relaxation rate constants. Therefore, the parameters $\left(1-S_f^2\right)$, $S_f^2\left(1-S_{\text{H}}^2\right)$, and $S_f^2{S}_{\text{H}}^S\left(1-S_{\text{N}}^2\right)$ are linear to $\left(1-S^2\right)\theta \left(z\right)$. The parameters $S_{\text{H}}^2$ and $S_{\text{N}}^2$ themselves are not linear to $\left(1-S^2\right)\theta \left(z\right)$, but may be obtained by simple arithmetic from the linear parameters. Like SDM, LeMaster’s approach is limited to data acquired at a single field.

Interpretation of Motions by a Projection onto an Array of Correlation Times Approach

Limitations of the approaches above have led Ferrage and coworkers to develop the interpretation of motions by a projection onto an array of correlation times (IMPACT) approach (Khan et al., 2015), which was applied to a protein with intrinsically disordered regions (IDR). A challenge of IDRs is that the lack of structure potentially yields a large number of distinct motions and therefore many correlation times, so that EMF approach is not appropriate for data analysis, but the limited number of parameters obtained with SDM fails to provide a complete description of the dynamics. Then, the IMPACT approach allows analysis of large, multi-field data sets, by taking the total correlation function to be a sum of several fixed correlation times, ${\tau }_k$, such that

$C\left(t\right)={\displaystyle \sum _k{A}_k{e}^{-t/{\tau }_k}}(25)$

Because $C\left(0\right)=1$ and decays to 0, the $A_k$ must sum to 1. For the Engrailed 2 protein, ¹⁵N T₁, NOE (${\sigma }_{\text{NH}}$), and transverse and longitudinal cross-relaxation rate constants at five fields (400, 500, 600, 800, 1,000 MHz) could be fit to an array of six correlation times, log-spaced between 21 ps and 21 ns. When fitting to Eq. 25, one restricts the amplitudes to remain between zero and one, and the sum of amplitudes must be set to one.

Following our procedures for SDM and LeMaster’s approach, we also examine the behavior of the IMPACT approach in Figure 8. When fitting a correlation function having a single correlation time in Figure 8A, we obtain ideal behavior from the IMPACT approach. When the input correlation time matches one of the correlation times in the IMPACT array, the corresponding amplitude is one, and all other amplitudes are zero. When the input correlation time is in between correlation times in the IMPACT array, then only the two nearest correlation times to the input value have non-zero amplitudes, and those two amplitudes sum to one (a minor deviation from this behavior occurs at 10 ns). However, if we input two correlation times in Figure 8B, or one correlation time and one distribution in Figure 8C, with motion split equally between the two correlation times or correlation time and distribution, the fit parameters’ response to the swept correlation time is not an exact reproduction of the behavior in Figure 8A, in contrast to SDM and LeMaster’s approach. While SDM and LeMaster’s approach are both linear combinations of relaxation rate constants, IMPACT is a linear fit for which its behavior depends heavily on restricting the values of the fit parameters (priors), which as discussed in A Few Notes on Linearity, means that the fit parameters are no longer linear to $\left(1-S^2\right)\theta \left(z\right)$. The result is that the response of the parameters $A_k$ to a given correlation time do depend weakly on other motions present, thus not fully obtaining the ideal, linear behavior of SDM and LeMaster’s approach. However, IMPACT provides a good approximation to this behavior, and is more generally applicable than SDM and LeMaster’s approach.

FIGURE 8

FIGURE 8. Behavior of the IMPACT approach as a function of correlation time. In each plot, we fit calculated relaxation rate constants, and fit the amplitudes in Eq. 25 according to the IMPACT procedure, using the set of experiments from Khan et al. (2015). In (A), the input total correlation function consists of a single decaying exponential term (with amplitude 1), where the amplitudes are plotted as the correlation time is varied. In (B), the total correlation function uses two correlation times (equal amplitudes), with one fixed at 1 ns, and the second swept (x-axis). On the y-axis, we plot contributions to the $A_k$ from the correlation time being varied. In (C), the same information is plotted, but the total correlation function includes a log-Gaussian distribution (μ = 630 ps, σ = 1 order of magnitude), and a single, variable correlation time.

IMPACT has not been developed for application to solid-state NMR, but it is worth investigating how such a method could work. In Figure 9A, we use an IMPACT-type approach to fitting R₁ at three fields and S², using an array of three correlation times. We restrict the fitted amplitudes ($A_k$) to fall between zero and one, but it does not make sense to require the $A_k$ to sum to one, since the correlation function in solid-state NMR does not usually decay to zero. Here, we assume a motion with just one correlation time, and $\left(1-S^2\right)$ = 0.3. Then, we find that IMPACT in solids is similar to its solution-state behavior. Note that the amplitudes corresponding to 1.4 and 5 ns capture motion near those correlation times, whereas the amplitude corresponding to 1 ps captures all motion not in proximity to 1.4 and 5 ns, including very slow motions. As with solution-state NMR, if we split the motion over two correlation times, and determine the response to one of the two correlation times Figure 9B, the response changes compared to fitting just the single correlation time. However, as discussed in A Few Notes on Linearity, and demonstrated with SDM and LeMaster’s approach, this dependence on other motions present vanishes if we eliminate restrictions on the fit parameters. Then, in Figure 9C, we repeat the fit from Figure 9A, without restrictions on the fit parameters, yielding reasonable behavior, excepting some negative amplitudes in the $A_k$. Figure 9D shows similar results, when fitting S² and R_1ρ, although the fitted correlation times must be in the sensitive range of the R_1ρ rate constants. Unfortunately, when we attempt to fit R₁ and R_1ρ simultaneously in Figure 9E, using the same correlation times as in Figures 9C,D, we find extremely unstable behavior. Apparently, we cannot simultaneously fit data on both sides of the solid-state NMR blind spot.

FIGURE 9

FIGURE 9. IMPACT behavior in solids. In each plot, we test the behavior of the amplitudes, $A_k$, using calculated solid-state NMR data (S², ¹⁵N R₁, ¹⁵N R_1ρ, with experimental conditions taken from Smith et al. (2016)). (A) plots the behavior of fitting S² and three R₁ rate constants to three correlation times (1 ps, 1.4 ns, 5 ns), where the input correlation function has a single correlation time ($\left(1-S^2\right)$=0.3), while restricting the$A_k$ to fall between 0 and 1. (B) shows fits under the same conditions, but includes two correlation times, with one fixed at 1 ns, and the other swept (x-axis). The y-axis plots the change in the $A_k$ due to the swept correlation time. (C) shows fits under the same conditions as (A), without restricting the values of the $A_k$. (D) also removes restrictions on the $A_k$, but fits S² and R_1ρ data, using correlation times of (1 ps, 2.5 μs, and 17.8 μs). (E) fits all data (S², R₁, R_1ρ) simultaneously without restrictions on the $A_k$, with correlation times of 1 ps, 1.4 ns, 5 ns, 2.5 μs, and 17.8 μs (F) fits R₁ data, but uses one very short correlation time (32 ps), and one very long correlation time (100 ns).

In Figures 9C,D, we have fairly good performance, excepting that some of the amplitudes become slightly negative. Interestingly, these negative amplitudes may be eliminated by placing two correlation times further away from each other. Then, in Figure 9F, we fit only R₁ data, using correlation times of 32 ps and 100 ns. The fitted correlation times no longer correspond to the center of the sensitive range of the $A_k$ (750 ps, 6.2 ns), and the amplitudes also far exceed the input value for $\left(1-S^2\right)$. Fitting while also including $S^2$ data allows using an additional correlation time (1 ps), but the corresponding $A_k$ becomes large and negative (not shown). From this final result, we could simply renormalize the amplitudes to have a maximum of one, and report the center of the sensitive range instead of the correlation times to which we actually fitted. The result would still be a linear combination of the experimental data, and therefore linear to $\left(1-S^2\right)\theta \left(z\right)$, but the result would have very little to do with the correlation times chosen to obtain that linear combination.

A New Approach for Solid-State Nuclear Magnetic Resonance

In the previous section, we investigated the behavior of a number of approaches to processing relaxation data. Of those approaches, model-free, SDM, and LeMaster’s approach provide parameters which are linear to the distribution of correlation times, $\left(1-S^2\right)\theta \left(z\right)$ (in some cases, some additional arithmetic operations are required to obtain the reported parameters, e.g., $\langle {\tau }_e\rangle$ is calculated from $S^2$ and $\left(1-S^2\right)\langle {\tau }_e\rangle$). IMPACT approximates this behavior, although heavy reliance on priors prevents perfect linearity. However, each approach is limited in its application to solid-state NMR data. Therefore, we have developed the detector analysis (Smith et al., 2018), which is a general method for processing relaxation data that maintains a linear relationship between fit parameters and the distribution of correlation times.

Linear Combination of Data

As we have emphasized for the above examples, one may obtain parameters that have a well-defined (linear) relationship to the distribution of correlation times by taking linear combinations of relaxation rate constants. Thus far, we have limited ourselves to very specific linear combinations: combinations that yield the spectral density, or combinations that are related to specific correlation times. However, why shouldn’t we use any linear combination that is optimized to give an ideal linear relationship to the distribution of correlation times, $\left(1-S^2\right)\theta \left(z\right)$? We first recall that the correlation function has been defined here as being a linear combination of decaying exponentials, defined by $\left(1-S^2\right)\theta \left(z\right)$, and its Fourier transform (also a series of linear combinations) must then also be linear to $\left(1-S^2\right)\theta \left(z\right)$.

$\begin{array}{l}C\left(t\right)=S^2+\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)e^{-t/\left({10}^z\cdot 1\text{s}\right)}dz}\\ J^{\left(\theta ,S\right)}\left(\omega \right)=\frac{2}{5}\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)\frac{{10}^z\cdot 1\text{s}}{1+{\left(\omega \cdot {10}^z\cdot 1\text{s}\right)}^2}dz}\end{array}(26)$

Here, we take $J^{\left(\theta ,S\right)}\left(\omega \right)$ to be the spectral density resulting from $\left(1-S^2\right)\theta \left(z\right)$. Then, any relaxation rate constant is a weighted sum of terms from the spectral density.

$\begin{array}{l}{R}_{\zeta }^{\left(\theta ,S\right)}={\displaystyle \sum _p{a}_p^{\zeta }{J}^{\left(\theta ,S\right)}\left({\omega }_p\right)}\\ ={\displaystyle \sum _p{a}_p^{\zeta }\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)}}\frac{{10}^z\cdot 1\text{s}}{1+{\left(\omega \cdot {10}^z\cdot 1\text{s}\right)}^2}dz\\ =\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)\underset{=R_{\zeta }\left(z\right)}{\underbrace{{\displaystyle \sum _p{a}_p\frac{{10}^z\cdot 1\text{s}}{1+{\left(\omega \cdot {10}^z\cdot 1\text{s}\right)}^2}}}}dz}\\ R_{\zeta }^{\left(\theta ,S\right)}=\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)R_{\zeta }\left(z\right)dz}\end{array}(27)$

$R_{\zeta }^{\left(\theta ,S\right)}$ is the relaxation rate constant for an experiment, indexed ζ, resulting from the distribution of correlation times, $\left(1-S^2\right)\theta \left(z\right)$. Coefficients $a_p^{\zeta }$ indicate the weightings of the spectral density for experiment ζ, sampled at frequencies ${\omega }_p$. Insertion of $J^{\left(\theta ,S\right)}\left(\omega \right)$ into this linear combination allows us to express $R_{\zeta }^{\left(\theta ,S\right)}$ as a linear function of $\left(1-S^2\right)\theta \left(z\right)$, where $R_{\zeta }\left(z\right)$ defines the linear relationship (we refer to this as the sensitivity).

Then, as is the case for model-free, SDM, and LeMaster’s approach, any sum of relaxation constants maintains linearity. Following our previous convention (Smith et al., 2018), we denote the sum as ${\rho }_n^{\left(\theta ,S\right)}$.

$\begin{array}{l}{\rho }_n^{\left(\theta ,S\right)}={\displaystyle \sum _{\zeta }{b}_{\zeta }{R}_{\zeta }^{\left(\theta ,S\right)}}\\ ={\displaystyle \sum _{\zeta }{b}_{\zeta }\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)R_{\zeta }\left(z\right)dz}}\\ =\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)\underset{={\rho }_n\left(z\right)}{\underbrace{{\displaystyle \sum _{\zeta }{b}_{\zeta }{R}_{\zeta }\left(z\right)}}}dz}\\ {\rho }_n^{\left(\theta ,S\right)}=\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right){\rho }_n\left(z\right)dz}\end{array}(28)$

Then, ${\rho }_n\left(z\right)$ defines the linear relationship between $\left(1-S^2\right)\theta \left(z\right)$ and ${\rho }_n^{\left(\theta ,S\right)}$. The subsequent question is, how do we find the best linear combinations of the experimental relaxation rate constants for analyzing our relaxation data?

Optimizing Detectors: The Relaxation-Rate Space Approach

While any linear combination of experimental relaxation rate constants yields a linear relationship between $\left(1-S^2\right)\theta \left(z\right)$ and the resulting ${\rho }_n^{\left(\theta ,S\right)}$, not all combinations are equally good choices. A few guidelines are, first, non-negativity of ${\rho }_n\left(z\right)$; we would like ${\rho }_n^{\left(\theta ,S\right)}$ to always increase when amplitude of motion increases, whereas negative regions of ${\rho }_n\left(z\right)$ could cause ${\rho }_n^{\left(\theta ,S\right)}$ to decrease with increasing amplitudes. Second, narrowness: we would like each ${\rho }_n^{\left(\theta ,S\right)}$ to report on a specific range of correlation times. Third, when the full set of relaxation data is analyzed, one should be able to back-calculate the experimental data (within some tolerance) from the parameters ${\rho }_n^{\left(\theta ,S\right)}$. This ensures that one captures all information in the experimental data (clearly, if the ${\rho }_n^{\left(\theta ,S\right)}$ can reproduce the experimental data, then the ${\rho }_n^{\left(\theta ,S\right)}$ must have retained the information in the experiments).

The question, then, is how to obtain optimized linear combinations satisfying the above requirements. Our initial answer to this question is the result of identifying a similar problem in a completely different field: When one sees the color of an object, its appearance depends on the distribution of wavelengths reflected (or emitted) by the object. The distribution of wavelengths is given by the spectral power distribution, $S\left(\lambda \right)$. Whereas $S\left(\lambda \right)$ is an infinite-dimensional description of the spectral power vs. wavelength, what is “seen” is a projection of that distribution onto a three dimensional space, corresponding to the three cones that detect color in the eye. This 3D space is often described using the CIE (Commission internationale de l'Eclairage) XYZ color space (Smith and Guild, 1931; Judd, 1951; Vos, 1978).

$\begin{array}{l}X={\displaystyle \underset{0}{\overset{\infty }{\int }}S\left(\lambda \right)\overline{x}\left(\lambda \right)d\lambda }\\ Y={\displaystyle \underset{0}{\overset{\infty }{\int }}S\left(\lambda \right)\overline{y}\left(\lambda \right)d\lambda }\\ Z={\displaystyle \underset{0}{\overset{\infty }{\int }}S\left(\lambda \right)\overline{z}\left(\lambda \right)d\lambda }\end{array}(29)$

The functions $\overline{x}\left(\lambda \right)$, $\overline{y}\left(\lambda \right)$, and $\overline{z}\left(\lambda \right)$ are plotted in Figure 10B. Based on the color one sees, one cannot define $S\left(\lambda \right)$ precisely, but certainly we learn something about the distribution of wavelengths. In the same way, based on a set of relaxation rate constants, we cannot fully define $\left(1-S^2\right)\theta \left(z\right)$, but certainly we can learn something about the dynamics. The matching forms of Eqs. 27, 29 further highlight the relationship between these problems.

FIGURE 10

FIGURE 10. Similarity between the CIE XYZ colorspace and the relaxation rate constant space. (A) plots the XYZ colorspace, black lines indicate where single wavelengths fall in the colorspace (z not shown, space is normalized such that x + y + z = 1). Points connected by a triangle indicate the definition of red, green, and blue colors as defined by the sRGB standard (Anderson et al., 1996). (B) plots the sensitivity of the $\overline{x}\left(\lambda \right)$, $\overline{y}\left(\lambda \right)$, and $\overline{z}\left(\lambda \right)$ color matching functions as a function of wavelength (λ). (C) plots sRGB sensitivities resulting from transformation from the XYZ to sRGB spaces. Points connected by triangles correspond to definitions of ${\overrightarrow{r}}_1$, ${\overrightarrow{r}}_2$, and ${\overrightarrow{r}}_3$ that define the detector space. (D) shows the normalized relaxation rate space for ¹³C R₁ at 300 and 800 MHz and H–C NOE at 800 MHz. (E) shows the sensitivities of each of these experiments a function of correlation time. (F) shows detector sensitivities resulting from transformation from the relaxation rate constant space to detector space (defined by the points in (D)).

The XYZ color space can be represented as a 2D space, shown in Figure 10A. Only x and y are shown, and z is selected so that $x+y+z=1$ (then, a third dimension would vary this sum, corresponding to brightness). By marking points in the color space, one can indicate how the color space may be represented in another basis. Here, we have marked points corresponding to red, green, and blue of the sRGB standard (Anderson et al., 1996). Colors within the resulting triangle may be obtained with positive linear combinations of the red, green, and blue of sRGB, so that this triangle is a good estimate of colors that may be obtained with a color monitor (which creates color by combining red, green and blue pixels–this means that in Figure 10A, colors outside the triangle are not correctly represented on your screen). These points also define a transformation from the XYZ color matching functions (Figure 10B) to the sRGB functions (Figure 10C). Note that any color may be represented in the sRGB space, but only those where $S\left(\lambda \right)$ results in positive R, G, and B values can actually be reproduced by a typical monitor.

Realizing that the mathematics of relaxation rate constants was essentially equivalent to color spaces, we created analogous relaxation rate constant spaces, replacing the X, Y, and Z values with normalized rate constants. However, instead of placing points within the relaxation rate space, we surrounded the space in Figure 10D, since we wanted to describe all points in the space with positive parameters. Interestingly, by surrounding the space as closely as possible, without crossing into the space, we obtained a transformation to functions with well-separated and non-negative sensitivities, see Figure 10F. In the example here, we use three points to transform the three experimental sensitivities into detector sensitivities, resulting in three detectors. However, redundancy in the information of larger data sets often results in the space becoming narrow in a given dimension, so that the full space may also be approximately described using fewer points, resulting in fewer detectors than experimental data points, but better signal-to-noise in the resulting parameters. Full details of this approach are described in Smith et al. (2018).

Optimizing Detector Sensitivities: Automated Approach

Investigating the relaxation rate space is a powerful way to grasp the information content of a relaxation data set, however, detector optimization using this method requires manual selection of points in the space. This quickly became excessively tedious for large data sets, as is the case for analysis of relaxometry data (Smith A. A. et al., 2021), so that we have also automated the optimization of linear combination (Smith et al., 2019a).

For automation, one still has the requirements that we capture the information in the experiments (that is, we can fit the data), while minimizing the number of parameters to describe that data, and second, that we obtain detector sensitivities that are narrow and non-negative. The first requirement may be met using singular value decomposition (Golub and Kahan, 1965). Suppose we have a matrix, M, for which each row is a sensitivity of one of our experiments ($R_{\zeta }\left(z\right)$), where we perform a normalization to prioritize fitting of higher quality data (procedure: first, we normalize all sensitivities to a maximum of one, second we multiply the sensitivity by the median of the experimental rate constants, and third we divide by the median standard deviation of those rate constants). Each column then corresponds to a correlation time. For N experiments, we obtain the best approximation of M that can be achieved with a linear combination of t vectors, defined by

$\begin{array}{l}\mathbf{M}\approx \tilde{\mathbf{M}}={\mathbf{U}}_t\cdot {\Sigma }_t\cdot {\mathbf{V}}_t^{\prime }\\ {\mathbf{V}}_t^{\prime }={\Sigma }_t^{-1}\cdot {\mathbf{U}}_t^{\prime }\cdot \mathbf{M}\end{array}(30)$

The t rows of ${\mathbf{V}}_t^{\prime }$ are linear combinations of the rows of M, with recombination defined by the product ${\Sigma }_t^{-1}\cdot {\mathbf{U}}_t^{\prime }\cdot \mathbf{M}$ (${\mathbf{U}}_t$, ${\mathbf{V}}_t^{\prime }$ are unitary matrices, and ${\Sigma }_t$ is diagonal, with the largest n singular values along the diagonal). Linear combination of the rows of ${\mathbf{V}}_t^{\prime }$ to yield the rows of M is an approximate relationship, but the inverse, recombination of the rows of M to yield ${\mathbf{V}}_t^{\prime }$, is exact. Then, the closer $\tilde{M}$ is to M, the better the data can be fit, but this requires t to be larger, and thus more noise is also present in the final analysis. In principle, this linear recombination could be directly applied to the experimental data, to obtain detectors with sensitivities given by the rows of the ${\mathbf{V}}_t^{\prime }$. The result would capture (approximately) the maximum amount of information possible from the experiment with t parameters. However, the sensitivities found in the rows of ${\mathbf{V}}_t^{\prime }$ are not narrow, and usually have large negative regions. On the other hand, a linear recombination of the vectors in ${\mathbf{V}}_t^{\prime }$ would maintain the information content and fit quality, but allows one to optimize the detector sensitivities to be separated and non-negative.

$\left[\begin{array}{c}{\rho }_1\left(z\right)\\ {\rho }_2\left(z\right)\\ ⋮\end{array}\right]=\mathbf{T}\cdot {\mathbf{V}}_t^{\prime }=\mathbf{T}\cdot {\Sigma }_t^{-1}\cdot {\mathbf{U}}_t^{\prime }\cdot \mathbf{M}(31)$

Then, T defines the linear recombination of the ${\mathbf{V}}_t^{\prime }$ to yield the ${\rho }_n\left(z\right)$, where T is a square matrix. The product of a row of T with ${\mathbf{V}}_t^{\prime }$ defines one of the detectors sensitivities, ${\rho }_n\left(z\right)$. A row of T is determined in order to optimize a detector sensitivity, first by choosing a single correlation time, $z_{\text{max}}={\log }_{10}\left({\tau }_c/\text{s}\right)$, for which we optimize a linear combination of the rows of ${\mathbf{V}}_t^{\prime }$ such that ${\rho }_n\left(z_{\text{max}}\right)$ = 1, while simultaneously minimizing ${\rho }_n\left(z\right)$ for all other correlation times, and requiring that all ${\rho }_n\left(z\right)$ remain non-negative. This can be quickly solved using a linear programming algorithm (Kantorovich, 1960; Dantzig, 1982; Virtanen, 2020). However, if we sweep through an array of correlation times, performing this optimization at each correlation time, we find that we are only successful at t correlation times (we consider the minimization as having failed if for some z, we find that ${\rho }_n\left(z\right)$ exceeds 1). Currently, we find the best t detectors by sweeping over a large array of correlation times (200), although this algorithm could be improved to reduce the number of optimizations required (spaces method and automated method both implemented in MATLAB, download from https://difrate.sourceforge.io).

In the detector analysis, once we have optimized the detectors, we apply the same linear combination to the experimental relaxation rate constants as were applied to the sensitivities in order to obtain optimized detector responses. Note in practice that this is implemented as a fit, allowing one to prioritize fitting relaxation rate constants with lower measurement error. Furthermore, we place bounds on the fitted detector responses, ${\rho }_n^{\left(\theta ,S\right)}$. In A Few Notes on Linearity, we noted that bounds (priors) on the fit parameters can cause the fit parameters to not be linear to $\left(1-S^2\right)\theta \left(z\right)$. This is only the case if the priors exclude the best fit. Detectors are constructed such that any allowed set of relaxation rate constants will not result in parameters that violate the priors. Allowed rate constants are any set that may occur for an arbitrary form of $\left(1-S^2\right)\theta \left(z\right)$. If, due to noise or measurement error, a dis-allowed set of relaxation rate constants is measured, then the priors will force the fitted relaxation rate constants to fall in the allowed space.

Model-Free, or Not?

We see that the original model-free approach, SDM, LeMaster’s approach, and detector analysis all belong to a family of methods that yield parameters with well-defined relationships to the distribution of correlation times, here defined by $\left(1-S^2\right)\theta \left(z\right)$. For SDM and detectors, the final parameters ($J\left(\omega \right)$, ${\rho }_n^{\left(\theta ,S\right)}$) are linearly related to $\left(1-S^2\right)\theta \left(z\right)$; for model-free, $S^2$ and $\left(1-S^2\right)\langle {\tau }_e\rangle$ are linear, and for LeMaster’s approach, $\left(1-S_f^2\right)$, $S_f^2\left(1-S_{\text{H}}^2\right)$, and $S_f^2{S}_{\text{H}}^2\left(1-S_{\text{N}}^2\right)$ are linear, whereas the final parameters ($S^2$, $\langle {\tau }_e\rangle$, $S_f^2$, $S_{\text{H}}^2$, and $S_{\text{N}}^2$) must be obtained via additional arithmetic operations. Response of EMF parameters, on the other hand, may react to changes in one motion differently, depending on other motions in the system. Still, its simplicity in analysis and interpretation—one to three pairs of correlation times and amplitudes—makes it an attractive choice for relaxation data analysis. Should we then compromise in some cases, and sacrifice well-defined parameters for more easily interpreted parameters?

Case 1: Extended Model-Free

Using detectors, we may better understand how EMF parameters in solid-state NMR depend on amplitudes of motion for particular windows of correlation times. We re-analyze relaxation data of HET-s (218–289) fibrils (Smith et al., 2016), by first performing a detector analysis on the data, shown in Figure 11A and then iteratively fitting detector responses to correlation times and amplitudes in Figure 11B, resulting in the EMF analysis in Figure 11C.

FIGURE 11

FIGURE 11. Model-free analysis from detectors. (A) shows a detector analysis of HET-s (218–289) fibrils (Smith et al., 2016), with sensitivities shown in (B) (amplitude scale not shown; sensitivities have a maximum of 1). (B) illustrates the procedure to convert 273Ser detector responses into model-free parameters. Bars give the detector responses (y-axis), plotted at the center of the corresponding detector’s sensitivity (x-axis, note that ρ₀, blue, does not have a well-defined center). At top, we find the ratio of ${\rho }_3^{\left(\theta ,S\right)}/{\rho }_2^{\left(\theta ,S\right)}$ is consistent with a correlation time of 34 ns, with corresponding amplitude of 0.075 (intermediate motion). After subtracting the contribution of this correlation time to ${\rho }_0^{\left(\theta ,S\right)}$ (middle), we find the ratio ${\rho }_1^{\left(\theta ,S\right)}/{\rho }_0^{\left(\theta ,S\right)}$ is consistent with a correlation time of 49 ps, and amplitude of 0.12 (fast motion). Using a fixed correlation time of 14.7 μs, we find an amplitude for the slow motion of 1.8 × 10⁻³ (bottom). (C) shows the results of EMF analysis for all residues using the procedure in (B).

Using the following procedure, we are able to reproduce our previous model-free results, illustrated in Figure 11B for residue 273Ser. The procedure is given below as a set of simple equations, where results are a good reproduction of our previous direct fit using the model-free approach.

$\begin{array}{l}\begin{array}{lll}\begin{array}{l}\mathbf{Step}\mathbf{1}:z_{\text{i}}\\ \frac{{\rho }_2^{\left(\theta ,S\right)}}{{\rho }_3^{\left(\theta ,S\right)}}=\frac{{\rho }_2\left(z_{\text{i}}\right)}{{\rho }_3\left(z_{\text{i}}\right)}\end{array}\hfill & \hfill & \begin{array}{l}\mathbf{Step}\mathbf{2}:A_{\text{i}}\\ A_{\text{i}}=\frac{{\rho }_2^{\left(\theta ,S\right)}}{{\rho }_2\left(z_{\text{i}}\right)}=\frac{{\rho }_3^{\left(\theta ,S\right)}}{{\rho }_3\left(z_{\text{i}}\right)}\end{array}\hfill \\ \begin{array}{l}\mathbf{Step}\mathbf{3}:z_{\text{i}}\\ \frac{{\rho }_0^{\left(\theta ,S\right)}-A_{\text{i}}{\rho }_0\left(z_{\text{i}}\right)}{{\rho }_1^{\left(\theta ,S\right)}-A_{\text{i}}{\rho }_1\left(z_{\text{i}}\right)}=\frac{{\rho }_0\left(z_{\text{f}}\right)}{{\rho }_1\left(z_{\text{f}}\right)}\end{array}\hfill & \hfill & \begin{array}{l}\mathbf{Step}\mathbf{4}:A_{\text{f}}\\ A_{\text{f}}=\frac{{\rho }_0^{\left(\theta ,S\right)}-A_{\text{i}}{\rho }_0\left(z_{\text{i}}\right)}{{\rho }_0\left(z_{\text{f}}\right)}=\frac{{\rho }_1^{\left(\theta ,S\right)}-A_{\text{i}}{\rho }_1\left(z_{\text{i}}\right)}{{\rho }_1\left(z_{\text{f}}\right)}\end{array}\hfill \end{array}\\ \mathbf{Step}\mathbf{5}:A_{\text{s}}\\ A_{\text{f}}=\frac{{\rho }_4^{\left(\theta ,S\right)}-A_{\text{i}}{\rho }_4\left(z_{\text{i}}\right)-A_{\text{f}}{\rho }_4\left(z_{\text{f}}\right)}{{\rho }_4\left(z_{\text{s}}\right)}\\ \begin{array}{ccc}\begin{array}{l}{S}_{\text{f}}^2=1-A_{\text{f}}\\ S_{\text{i}}^2=1-A_{\text{i}}/S_{\text{f}}^2\\ S_{\text{s}}^2=1-A_{\text{s}}/\left(S_{\text{f}}^2{S}_{\text{i}}^2\right)\end{array}& & \begin{array}{l}{\tau }_{\text{f}}={10}^{z_{\text{f}}}\cdot 1\text{ s}\\ {\tau }_{\text{i}}={10}^{z_{\text{i}}}\cdot 1\text{ s}\\ {\tau }_{\text{s}}={10}^{z_{\text{s}}}\cdot 1\text{ s}\end{array}\end{array}\end{array}(32)$

In the first and second steps, we find a correlation time for which the ratio of sensitivities of ρ₂ and ρ₃ matches the ratio of the detector responses, and then subsequently find the correct amplitude to reproduce these correlation times. With ${\rho }_2^{\left(\theta ,S\right)}$ = 2.0 × 10⁻², and ${\rho }_3^{\left(\theta ,S\right)}$ = 2.2 × 10⁻³, we find ${\tau }_{\text{i}}$ = 34 ns. Our first concern with this fit is that the intermediate correlation time, $z_{\text{i}}={\log }_{10}\left({\tau }_{\text{i}}/\text{s}\right)$, is a compromise between a detector sensitive to motions around 6 ns and a second sensitive around 2 μs. It seems unlikely that the same motion can really explain these two detector responses, which have sensitivities separated by three orders of magnitude. The second problem is because we use a compromise correlation time, both detector sensitivities are very low at this correlation time, which must be counterbalanced by using a large amplitude ($A_{\text{i}}$) in the model-free fit. Then, in our example, $A_{\text{i}}$ = 0.075 is significantly larger than the detector responses, ${\rho }_2^{\left(\theta ,S\right)}$ and ${\rho }_3^{\left(\theta ,S\right)}$, from which it results, so that we are very likely overestimating the amplitude of this motion.

In the third and fourth steps, we subtract the contributions from $z_{\text{i}}$ and $A_{\text{i}}$ from ${\rho }_0^{\left(\theta ,S\right)}$ and ${\rho }_1^{\left(\theta ,S\right)}$, and similarly use the ratios of the remainder of the detector responses to obtain $z_{\text{f}}$, and their amplitudes to obtain $A_{\text{f}}$. Again, it is not clear if these detectors should be treated as if they describe a single motion. In particular, the relatively uniform behavior of ${\rho }_0^{\left(\theta ,S\right)}$ likely is a result of primarily local librational motion, which will not be described by the same amplitudes and correlation times of motions leading to greater variation in ${\rho }_1^{\left(\theta ,S\right)}$. Interestingly, because the amplitudes do not vary in the same way, the variation in amplitude of ${\rho }_1^{\left(\theta ,S\right)}$ cannot be reproduced in the trends for $A_{\text{f}}$, but instead has to be fitted by variation in correlation time (${\tau }_{\text{f}}$). The result is that amplitude trends in $S_{\text{i}}^2=1-A_{\text{i}}/S_{\text{f}}^2$, shown in Figure 11(C, middle) correlate well with trends in ${\tau }_{\text{f}}$, especially near breaks between the β-sheets of HET-s (near 235Glu, 271Gly). However, this correlation is actually coming from similar amplitude trends observed for ${\rho }_1^{\left(\theta ,S\right)}$ and ${\rho }_2^{\left(\theta ,S\right)}$. The corresponding detector sensitivities are centered at 760 ps and 6.1 ns, and in fact overlap, suggesting that they may describe the same or at least related motions. EMF attributes these detector responses to different motions, having median correlation times of 22 ps and 42 ns (taken over all residues), thus being separated by three orders of magnitude.

In the final step, one fixes the slow correlation time to 14.7 μs (based on a fit optimization over the whole data set). In this case, the amplitude of ${\rho }_4^{\left(\theta ,S\right)}$ determines $A_{\text{s}}$ alone; the proximity of 14.7 μs to the center of ${\rho }_4$ (24 μs) results in fairly reasonable amplitudes (for 273Ser, ${\rho }_4^{\left(\theta ,S\right)}$ and $A_{\text{s}}$ fall within rounding error, yielding 1.8 × 10⁻³).

Then, the major problems with this EMF analysis are intermediate correlation times falling within the NMR blind spot (∼20–600 ns), along with correspondingly inflated amplitudes, as well as similar problems due to fitting fast correlation times to ${\rho }_0^{\left(\theta ,S\right)}$ and ${\rho }_1^{\left(\theta ,S\right)}$, which requires a compromise correlation time between librational motions (∼ps) with nanosecond motions. Furthermore, this behavior prevented comparison of EMF parameters for HET-s to MD results, whereas detectors yielded reasonable agreement (Smith et al., 2019b). As we have previously pointed out (Smith et al., 2017), the model-free parameters in HET-s fibrils are far from being atypical, in fact they are fairly consistent across multiple protein systems, likely due most studies utilizing similar data sets and analysis methodology (Chevelkov et al., 2009b; Schanda et al., 2010; Haller and Schanda, 2013; Zinkevich et al., 2013; Lamley et al., 2015a).

Case 2: Model-Free Analysis of μs-Motion

Microsecond motion is the result of processes having higher free-energy cost than nanosecond and picosecond dynamics. We suggest dividing these motions into local and collective motions, where the free energy cost of local motions comes from higher amplitude motions (∼10°) that require traversing a large energy barrier. In contrast, collective motions tend to be very low amplitude motion, where the high free-energy cost of the motion is not due to large amplitude dynamics or a significant energy barrier, but rather diffusive dynamics involving large numbers of atoms. Such dynamics are characterized by modes of motion, where a continuum of possible correlation lengths leads to a distribution of correlation times. In contrast, some local microsecond dynamics can be reasonably well approximated as a hopping motion between two orientations, and therefore described with a single correlation time (although effort should be made to determine whether relaxation might be due to multi-site exchange, and understand how this changes the interpretation of data analysis).

Local Dynamics

The availability of R_1ρ data, including formulas for its analysis (Trott and Palmer, 2002; Abergel and Palmer, 2003; Miloushev and Palmer, 2005; Kurbanov et al., 2011; Rovo and Linser, 2017) and improving methods for its collection (Kurauskas et al., 2017; Lakomek et al., 2017; Keeler et al., 2018; Krushelnitsky et al., 2018) has recently resulted in considerable improvement in the ability to characterize local micro- to millisecond motions (Rovó, 2020). We consider two categories of R_1ρ experiments: the first is Bloch-McConnell relaxation dispersion experiments (BMRD), for which R_1ρ relaxation is the result of motion modulating the isotropic chemical shift, and the NEar Rotary-resonance Relaxation Dispersion (NERRD, (Kurauskas et al., 2017)), for which orientational changes in anisotropic tensors leads to R_1ρ relaxation. For two-site exchange, BMRD R_1ρ relaxation rate constants depend on exchange rate ($k_{\text{ex}}=1/{\tau }_c$), the change in isotropic chemical shift due to exchange ($\Delta {\omega }_{12}$), and the population ($p_1,p_2=1-p_1$). Rate constants further depend on the effective field strengths corresponding to each of the two chemical shifts, ${\omega }_{\text{e}1}$ and ${\omega }_{\text{e}2}$, as well as the effective field for the average chemical shift. Palmer and coworkers provide us with the following expression (Trott and Palmer, 2002; Trott et al., 2003; Miloushev and Palmer, 2005), which is valid in the fast or intermediate exchange regimes:

$\begin{array}{l}{R}_{1\rho }=\frac{R_1}{2}\left(1+{\cos }^2{\beta }_{\text{e}}\right)+R_{1\rho }^{\text{DD,CSA}}+R_{1\rho }^{\text{ex}}\\ R_{1\rho }^{\text{ex}}=\frac{{\sin }^2{\beta }_{\text{e}}{p}_1{p}_2\Delta {\omega }_{12}^2{k}_{\text{ex}}}{\frac{{\omega }_{\text{e}1}^2{\omega }_{\text{e}2}^2}{{\omega }_{\text{e}}^2}+k_{\text{ex}}^2-{\sin }^2{\beta }_{\text{e}}{p}_1{p}_2\Delta {\omega }_{12}^2\left(1+\frac{2k_{\text{ex}}^2(p_1{\omega }_{\text{e}1}^2+p_2{\omega }_{\text{e}2}^2)}{{\omega }_{\text{e}1}^2{\omega }_{\text{e}2}^2+{\omega }_{\text{e}}^2{k}_{\text{ex}}^2}\right)}\\ \Omega =p_1{\Omega }_1+p_2{\Omega }_2\\ {\omega }_{\text{e}}^2={\omega }_1^2+{\Omega }^2\\ {\omega }_{\text{e}1}^2={\omega }_1^2+{\Omega }_1^2,{\omega }_{\text{e}2}^2={\omega }_2^2+{\Omega }_2^2\end{array}(33)$

The total $R_{1\rho }$ relaxation has contributions from longitudinal relaxation (R₁), transverse relaxation from dipole and CSA tensors ($R_{1\rho }^{\text{DD,CSA}}$), and from chemical exchange ($R_{1\rho }^{\text{DD,CSA}}$). (Kurbanov et al., 2011) give the formula for $R_{1\rho }^{\text{DD,CSA}}$.

$\begin{array}{l}{R}_{1\rho }^{\text{DD},\text{CSA}}={\sin }^2{\beta }_e\times [{\left(\frac{\delta }{4}\right)}^2\left(J\left({\omega }_{\text{S}}\right)+\frac{1}{3}J\left(2{\omega }_{\text{r}}-{\omega }_{\text{e}}\right)+\frac{2}{3}J\left({\omega }_{\text{r}}-{\omega }_{\text{e}}\right)+\frac{2}{3}J\left({\omega }_{\text{r}}+{\omega }_{\text{e}}\right)+\frac{1}{3}J\left(2{\omega }_{\text{r}}+{\omega }_{\text{e}}\right)\right)\\ +\frac{2}{27}{\left({\omega }_{\text{I}}\Delta {\sigma }_{\text{I}}\right)}^2\left(\frac{1}{2}J\left(2{\omega }_{\text{r}}-{\omega }_{\text{e}}\right)+J\left({\omega }_{\text{r}}-{\omega }_{\text{e}}\right)+J\left({\omega }_{\text{r}}+{\omega }_{\text{e}}\right)+\frac{1}{2}J\left(2{\omega }_{\text{r}}+{\omega }_{\text{e}}\right)\right)]\end{array}(34)$

${\omega }_{\text{r}}$ is the magic angle spinning frequency, and ${\omega }_{\text{e}}$ is the effective field as defined above. If one assumes the microsecond dynamics are dominated by two-site hoping, the spectral density is given simply by

$J\left(\omega \right)=\frac{2}{5}3p_1{p}_2\left(1-{\cos }^2\theta \right)\frac{k_{\text{ex}}}{k_{\text{ex}}^2+{\omega }^2}=\frac{2}{5}\left(1-S^2\right)\frac{{\tau }_c}{1+{\left(\omega {\tau }_c\right)}^2}(35)$

Then, the question is, how may we most efficiently extract the exchange rate ($k_{\text{ex}}=1/{\tau }_c$), populations ($p_1,p_2$), chemical shift changes ($\Delta {\omega }_{12}$), and angle changes ($\theta$). Fitting of $k_{\text{ex}}$ has been fairly well established using both NERDD or BMRD (Trott et al., 2003; Ma et al., 2014; Rovo and Linser, 2017; Marion et al., 2019), and combining both methods should improve the accuracy of the resulting $k_{\text{ex}}$. However, separation of populations from either $\theta$ (NERDD) or $\Delta {\omega }_{12}$ (BMRD) is non-trivial. Supposing we already know $k_{\text{ex}}$, a given experiment’s relaxation rate constant then depends on the populations and either $\theta$ or $\Delta {\omega }_{12}$ (at sufficiently fast MAS, a given effective field usually results in either $R_{1\rho }^{\text{DD,CSA}}$ or $R_{1\rho }^{\text{ex}}$ being dominant, although in principle both terms are active in the same experiments). Inspecting Eqs. 34, 35 we note that terms $p_1,p_2$, and $\theta$, only appear once as a product of terms, $3p_1{p}_1\left(1-{\cos }^2\theta \right)$. Then, based on NERRD data alone, these parameters are inseparable. This is seen in Figure 12, where we plot $R_{1\rho }$ as a function of $p_1$ and $\theta$. We also calculate $R_{1\rho }^{\text{DD,CSA}}$ specifically for $p_1$ = 0.25 and $\theta ={16}^{°}$, and then indicate all other positions resulting in the same value of $R_{1\rho }^{\text{DD,CSA}}$ Figures 12A,B as a black contour. In Figures 12E,F, we only show contours where $R_{1\rho }^{\text{DD,CSA}}$ matches the value obtained for $p_1$ = 0.25 and $\theta ={16}^{°}$, but show several different experimental conditions (varying the field strength, ${\nu }_1={\omega }_1/2\pi$). Because this results in identical contours, we are unable to disentangle these parameters based on NERRD experiments under different conditions.

FIGURE 12

FIGURE 12. Separating population from hop angle and change in chemical shift in NERRD and BMRD experiments. Relevant parameters are shown as insets (${\omega }_r/2\pi$=40 kHz for NERDD plots). In (A–D), contour plots are shown for NERDD and BMRD relaxation rate constants under various conditions, and in each plot, a contour shows all values of $p_1$ and $\theta$ or $\Delta {\omega }_1$ that yield $R_{1\rho }$ equal to the value obtained for $p_1$ = 0.25 and $\theta ={16}^{°}$ or $\Delta {\omega }_{\text{I}}$ = 500 Hz (marked as a cross on each plot). In (E–L), we only show the contour, but for a range of experimental conditions (five experiments, linearly spaced, with range indicated in the plot). In some cases, this yields nearly identical contours, such that we only see one of the five contours.

In contrast, $R_{1\rho }$ relaxation resulting from chemical exchange has a more complex dependence on the various parameters. In particular, effective fields for each of the two states in exchange, ${\omega }_{\text{e}1}$ and ${\omega }_{\text{e}2}$ depend on the different offsets, ${\Omega }_1$, ${\Omega }_2$, but do not depend on the populations, in principle making the terms separable. Indeed, several plots in Figure 12 show that different experimental conditions lead to different contours for $p_1$ vs. $\Delta {\omega }_{12}$ (contours correspond to $R_{1\rho }$ that is equal to $R_{1\rho }$ obtained for $p_1$ = 0.25 and $\Delta {\omega }_{12}$ = 500 Hz, where contour intersections yield the input values). We are then able to identify the critical conditions required for separating population from chemical shift change. First, we see that if $k_{\text{ex}}\gg \Delta {\omega }_{12}$, contours are fully overlapped so that we are not able to separate the terms, shown in Figure 12G. This is because, in Eq. 33, $k_{\text{ex}}^2$ must be much larger than the last term in the denominator. If it is also larger than ${\omega }_{\text{e}1}^2{\omega }_{\text{e}2}^2/{\omega }_{\text{e}}^2$, then the critical dependence of the $R_{1\rho }$ on ${\omega }_{\text{e}1}$ or ${\omega }_{\text{e}2}$ is lost. In case $k_{\text{ex}}^2$ is not larger than ${\omega }_{\text{e}1}^2{\omega }_{\text{e}2}^2/{\omega }_{\text{e}}^2$, then the effective field must be much larger than $\Delta {\omega }_{12}$, so that this term converges on ${\omega }_e^2$, again losing dependence on ${\omega }_{\text{e}1}$ and ${\omega }_{\text{e}2}$ (i.e. the denominator simplifies to ${\omega }_{\text{e}}^2+k_{\text{ex}}^2$ (Trott and Palmer, 2002)). In any case, if the effective fields become large, ${\omega }_{\text{e}1}\to {\omega }_e$, ${\omega }_{e2}\to {\omega }_e$, similarly preventing separation in terms. For example, see Figures 12J,K, where a large offset or large field strength on the spin-locking field results in overlapping contours. Finally, note that we require a frequency offset to be applied in order to obtain the sign of $\Delta {\omega }_{12}$. If no frequency offset is applied, then all contours are symmetric as in Figure 12H.

Separability occurs only when $k_{\text{ex}}$, $\Delta {\omega }_{12}$, and ${\omega }_{\text{e}}$ are of similar size. Restricting ${\omega }_{\text{e}}$ is particularly challenging in solid-state NMR, where coherent effects may contribute to relaxation when the spin-locking field becomes too small (Öster et al., 2019). One approach would be to use increasing spinning frequencies (Penzel et al., 2015; Lakomek et al., 2017), although we note that some of the most clear improvements in Figure 12 occur in Figure 12L, where the field strength is only a few times bigger than the H–N J-couplings, which cannot be averaged by spinning.

In case we are in the fast exchange limit for BMRD experiments, we are left only with the terms $p_1{p}_2\left(1-{\cos }^2\theta \right)$ from NERRD experiments and $p_1{p}_2\Delta {\omega }_{12}^2$ from BMRD experiments. In this case, there is little to be done to fully separate population from the other parameters. If values of $\Delta {\omega }_{12}$ may be bounded, it is then possible to also bound $p_1{p}_2$, and therefore one finds a restricted range for possible values of $\theta$ (the reverse approach also works). However, if we are in the range of intermediate exchange, then we may separate populations from $\Delta {\omega }_{12}$, and use the result to also obtain $\theta$ (note that inclusion of NERDD data should additionally improve the accuracy of $k_{\text{ex}}$, which in turn improves separation of $\Delta {\omega }_{12}$ from populations based on the BMRD data). Note that Marion et al. have recently presented similar arguments (Marion et al., 2019), although separation of terms was apparently achieved by combining NERRD and BMRD data for fairly fast exchange ($k_{\text{ex}}$=18,000 s⁻¹, $\Delta {\omega }_{12}/2\pi$ = 240 Hz). While we agree that using both data sets together is beneficial, the information to separate population must come from the BMRD data and this is only possible in the intermediate exchange regime (Marion et al. calculated R_1ρ for a set of conditions, and via a coarse grid search, were able to find the initial conditions, however, other solutions along contours as in our Figure 12 likely were overlooked in the grid search).

A final consideration when analyzing BMRD and NERRD data is whether a two-site exchange model is reasonable. In a true two-site exchange, all moving residues should have identical exchange rates and populations, but differing $\Delta {\omega }_{12}$ and $\theta$ values. Then, validation of the two-site model could be achieved by independently analyzing all residues and establishing that all fits have approximately the same $p_1$, $p_2$, and $k_{\text{ex}}$ (or just the same $k_{\text{ex}}$ if populations cannot be determined). In case the true behavior is, for example, three-site exchange, fitting to the two-site exchange model will yield exchange rates that are a weighted average of the two non-zero eigenvalues of the 3 × 3 exchange matrix, where weighting will depend on the chemical shifts of the three sites and/or the angles sampled. In this case, it may be appropriate to apply a three-site exchange model, while jointly fitting all residues using a common set of rate constants (four to six independent parameters, depending on the model chosen). Such an approach has been demonstrated using CPMG relaxation in solution-state NMR (Korzhnev et al., 2005; Neudecker et al., 2006), with the general equations solved for CPMG (Koss et al., 2018).

Collective Dynamics

NERDD relaxation also appears throughout the whole protein in the absence of BMRD relaxation, depending on sample conditions, and is attributed to low amplitude rocking of the whole protein. This is observed very weakly in GB1 crystals (Krushelnitsky et al., 2018), and strongly in GB1 complexed with IgG (Lamley et al., 2015b), HET-s (218–289) (Smith et al., 2016), ubiquitin crystals with amplitude depending heavily on crystal form (Ma et al., 2015; Kurauskas et al., 2017; Lakomek et al., 2017), and SH3 (Krushelnitsky et al., 2018). The apparent global nature of this motion led all of these studies, with the exception of Lakomek et al., to fit R_1ρ relaxation using a slow motion with a single correlation time for all residues. In our HET-s analysis, we proposed fitting R_1ρ data using a global, slow correlation time, where the corresponding order parameter could vary, and additionally an offset term that would account for faster motion that could not be fully parameterized from R_1ρ data alone. Kurauskas et al. also followed this procedure, whereas Krushelnitsky and coworkers included explicit fitting of an additional fast motion with a distribution of correlation times. By including an offset term, and using a single correlation time globally, we again have a linear fit.

$R_{1\rho }=\underset{\text{ns}\text{motions}}{\underbrace{R_{1\rho }^0}}+\left(1-S_{\text{s}}^2\right)\underset{\mu \text{s}\text{motion},\text{fixed}}{\underbrace{R_{1\rho }\left({\tau }_{\text{s}}\right)}}(36)$

Then, for each residue, $R_{1\rho }^0$ and $\left(1-S_{\text{s}}^2\right)$ are varied, where $R_{1\rho }^0$ in principle compensates for relaxation due to fast, nanosecond motion, and $\left(1-S_{\text{s}}^2\right)$ should determine the effective amplitude of the global motion, with correlation time ${\tau }_s$, on the given residue. Practically, what happens is that the R_1ρ rate constants measured for a given residue have certain ratios. If those ratios match the ratios calculated for ${\tau }_{\text{s}}$, then $R_{1\rho }^0=0$ and the relaxation rate constants are fitted only with ($1–S_{\text{s}}^2$). In contrast, if all rate constants are approximately equal, then $\left(1-S_{\text{s}}^2\right)=0$ and $R_{1\rho }^0$ accounts for the full relaxation. However, in most cases, the ratios are closer to one than predicted by ${\tau }_{\text{s}}$, but not exactly one and so by including contributions from $R_{1\rho }^0$ and $\left(1-S_{\text{s}}^2\right)R_{1\rho }\left({\tau }_{\text{s}}\right)$, the data may be fit. One may investigate in more detail how the two terms vary as a function of correlation time (as in Figures 6–8). We show the behavior for the ¹⁵N and ¹³Cα R_1ρ data sets found in Smith et al. (2016) in Figure 13.

FIGURE 13

FIGURE 13. Behavior of fitting R_1ρ data to an offset and a fixed correlation time. (A) shows the offset term, $R_{1\rho }^0$, divided by 2000, and the order parameter for the slow correlation time, $(1-S_{\text{s}}^2)$ resulting from fitting calculated relaxation rate constants as a function of correlation time to Eq. 36. Experiments are ¹⁵N R_1ρ acquired with MAS frequency of 60 kHz and spin-lock strengths of 11, 16, 25, 38, and 51 kHz, and ${\tau }_{\text{s}}$ is fixed at 18.5 μs. (B) shows detector sensitivities optimized using the same data set. (C) shows $R_{1\rho }^0$ divided by 2000 and $\left(1-S_{\text{s}}^2\right)$ for ¹³C R_1ρ acquired with MAS frequency of 60 kHz and spin-lock strengths of 9, 18, 35, and 48 kHz, as well as an additional experiment with MAS frequency of 40 kHz and spin-lock strength of 25 kHz ${\tau }_{\text{s}}$ is fixed at 7.0 μs. (D) shows detector sensitivities optimized using the same data set.

In Figure 13A, we calculate R_1ρ relaxation rate constants for ¹⁵N relaxation, and fit to Eq. 36. $\left(1-S_{\text{s}}^2\right)$ reaches a maximum of approximately one at 19 μs, so that this parameter describes motion at and around the fixed correlation time of ${\tau }_{\text{s}}$ = 18.5 μs. On the other hand, the offset term, $R_{1\rho }^0$, actually is most sensitive at 2.5 μs, far from fitting primarily fast, nanosecond motion. We see that the functional forms are similar to detector sensitivities optimized from the same data set, Figure 13B. In Figure 13C, the behavior is less ideal: $\left(1-S_{\text{s}}^2\right)$ reaches a maximum of 1.28 at 13 μs, somewhat removed from the fixed correlation time of 7.0 μs, and the offset term becomes negative for correlation times around 18 μs.

The sensitivity of the offset term in Figure 13A to motion near 2.5 μs as opposed to faster motions may be surprising, although perhaps it should not be. NERRD experiments are most sensitive in the μs-range of correlation times, and rate constants under different experimental conditions have nearly converged to the same value at 1.9 μs (all rate constants within 5% of each other)– only slightly faster than the 2.5 μs where we find the maximum. Then, we would expect the offset term to be sensitive both near where R_1ρ is most sensitive, but also near where it converges, which is roughly what we find.

It is then important to note that fitting $R_{1\rho }$ to contributions from an offset term and a fixed correlation time results in an offset term that is most sensitive not to fast (nanosecond) motions, but rather to slower (microsecond) motions. In some cases, $\left(1-S_{\text{s}}^2\right)$ may be overly sensitive to some correlation times, with sensitivity exceeding one at positions that are removed from ${\tau }_{\text{s}}$. Detectors are also a better choice for characterizing broad distributions of correlation times, if one does not know the form of the distribution. In fact, we suspect that global rocking motion is the result of collective dynamics over varying length scales, where increasing the correlation length also increases the correlation time, and therefore yields a broad distribution of correlation times. We demonstrated the relationship between correlation length and correlation time window for HET-s fibrils on the nanosecond timescale using a combination of NMR and MD simulation (Smith et al., 2019b), however, the question remains whether similar behavior can fully explain rocking motion of crystalline proteins; for example, Schanda and coworkers argue that a coupling between overall rocking motion and local loop motion may exist in crystalline ubiquitin (Kurauskas et al., 2017).

Outlook: Combining Methods

We have seen that relaxation data in NMR may be processed by a variety of different methods, however, only some of these methods can really be thought of as “model-free,” such that we can establish a well-defined (linear) behavior for each parameter as a function of correlation time, independent of the actual model of the correlation function. These methods are the original model-free analysis, under the assumption that $\omega {\tau }_i\ll 1$, spectral density mapping, LeMaster’s approach, and detector analysis. Of these, only detector analysis is generally applicable to solid-state NMR.

So, are detectors the last word in NMR dynamics analysis? We certainly hope not. Each detector response provides a “window” into the total reorientational motion of some NMR tensor, with the window width and center defined by ${\rho }_n\left(z\right)$. Still, such a description is not very precise: a moderate detector response could result from a low amplitude motion near where ${\rho }_n\left(z\right)$ reaches its maximum, it could result from a high amplitude motion where ${\rho }_n\left(z\right)$ is small, or (and we suspect this is often the case), it characterizes a distribution of correlation times that overlaps the range of sensitivity of that detector. A collection of detectors, and their behavior as a function of position in a molecule gives further hints at the dynamics of a molecule, but leaves much to be desired in terms of details of motion. What we would rather have is better models of motion. If we use a good model, based on knowledge of the dynamics obtained from other methods, the information added to our experimental data should improve our interpretation of the experiment.

Molecular dynamics simulation is particularly powerful as a complimentary method to NMR. One obtains positions of all atoms as a function of time, allowing first, the direct calculation of the NMR-relevant correlation functions, and second, in principle allowing one to connect those correlation functions to specific motion in the molecule. $C\left(t\right)$ is explicitly calculated as

$\begin{array}{l}C\left(t_n\right)=\frac{1}{N}{\displaystyle \sum _{i=0}^{N-n-1}{P}_2\left(\overrightarrow{\mu }\left({\tau }_i\right)\cdot \overrightarrow{\mu }\left({\tau }_{i+n}\right)\right)}\\ \approx S^2+\left(1-S^2\right){\displaystyle \underset{-\infty }{\overset{\infty }{\int }}\theta \left(z\right)\underset{R_{C\left(t_n\right)}\left(z\right)}{\underbrace{e^{-t_n/\left({10}^z\cdot 1\text{s}\right)}}}dz}\end{array}(37)$

This is the discrete form of Eq. 14, as would be applied to an MD trajectory. To obtain the nth time point in the correlation function, $C\left(t_n\right)$, we simply average over all pairs of frames separated by n frames. The latter equation is our assumed form for the correlation function, where we note that a given time point of the correlation function, $C\left(t_n\right)$, is related to the distribution of correlation times with the same functional form as the relaxation rate constants (excepting the offset, $S^2$, see Eq. 27). This allows one to calculate detectors from the collection of time points in MD-derived correlation functions using a procedure nearly identical to that described in Optimizing Detector Sensitivities: Automated Approach, where the sensitivity, $R_{\zeta }\left(z\right)$ (Eq. 27), is replaced by the term $R_{C\left(t_n\right)}\left(z\right)=\exp \left(-t_n/\left({10}^z\cdot 1\text{s}\right)\right)$. In fact, not only may detector analysis be easily modified to analyze MD-derived data, but it is a general approach to numerically solving the inverse Laplace transform, which avoids some of the pitfalls of more common regularization approaches (Tikhonov and Arsenin, 1977).

When analyzing MD with detectors, one has two options: find the optimal set of detectors for describing correlation time distributions found with MD (that is, as many as possible with good signal-to-noise, and as narrow/non-overlapping as possible), or optimize the detectors to match some or all of the NMR-derived detectors. The latter approach is shown in Figure 14, where sensitivities of seven NMR experiments in Figure 14A are optimized to yield five detectors in Figure 14C, and the linear combination used to yield ${\rho }_2\left(z\right)$ is explicitly illustrated in Figure 14B. From MD, time points in the correlation function may also be linearly combined (sensitivities for 11 time points shown in Figure 14D), to match the NMR-derived detectors Figures 14E,F. Note that in Figure 14F, the linear combination is a very good match for ρ₁ and ρ₂, with moderate success for ρ₀, but detector sensitivities in the microsecond range are badly reproduced. The detector optimization indicates (correctly) that a 1 μs trajectory cannot reasonably predict dynamics in the range of several microseconds (ρ₃, ρ₄ in red, violet), thus providing a means for determining what information can and cannot be compared across methods. Where sensitivities agree, quantitative comparison of dynamics in MD and NMR is possible. Note that in Figures 14D–F, we only show 11 time points for illustrative purposes, but this procedure is equally valid for ∼10⁶ time points (for such a long correlation function, calculating Eq. 37 takes much longer than evaluating its result with detectors).

FIGURE 14

FIGURE 14. Combining NMR and MD. (A) plots normalized NMR sensitivities for a selection of experiments (S², ¹⁵N R₁ at 400, 500, 850 MHz, ¹⁵N R_1ρ at 850 MHz, 60 kHz MAS, ν₁ at 10, 25, and 45 kHz). (B) shows a linear combination of the normalized sensitivities (x, y positions shifted to reduce plot overlap), which yields the sensitivity of ρ₂, shown in (C) (green, bold). In color are the weighted contributions from each rate constant, and grey shows the cumulative sum (summing all sensitivities at and below the grey line). (C) shows the five sensitivities optimized from NMR data. (D) plots sensitivities of time points from MD-derived correlation functions (0 s, 10 points log-spaced from 50 ps to 1 μs). (E) shows a linear combination of those sensitivities, optimized to match the sensitivity of ρ₂ (x, y positions shifted to reduce plot overlap). (F) shows detectors optimized to match the NMR-derived detectors in (C). (G) shows spatial correlation of motion in a helix as a function of correlation time (windows for <20, ∼20, ∼100, ∼800 ns). Color intensity and bond radii indicate the correlation coefficient between that residue’s H–N motion and the motion of the black residue. (H) illustrates frames used to separate transformation from the PAS to the lab frame into four steps: a peptide plane frame, a helix frame, and a molecule frame (illustration inspired by Brown (1996), molecule plots created with ChimeraX (Pettersen et al., 2021)).

The detector analysis then provides a very reliable means of comparing NMR results to MD simulation. The ability to easily compare results across multiple methods is one of the primary advantages of detector analysis. We should note that carefully executed fitting of MD-derived correlation functions, followed by calculation of relaxation rate constants should yield similarly reliable rate constants, if the trajectory is sufficiently long (Mollica et al., 2012). However, the rate constants themselves are sensitive to a broader range of correlation times than detectors, so that the comparison has lower timescale resolution than detectors.

With MD and NMR data sets, one may then use NMR data via detectors (or relaxation rate constants) as a means of validating the MD, and potentially refining it; methods include selecting sections of trajectories that best reproduce experiment (Salvi et al., 2016), selecting the best force fields for a system (Antila et al., 2021), or validating the refinement of a force-field itself (Hoffmann et al., 2018a; Hoffmann et al., 2018b). One may also use NMR data (specifically order parameters) as a means of directing the simulation, so that the simulation returns parameters matching the experiment (Hansen et al., 2014). One should note that a major challenge of combining NMR and MD data is that, while NMR is highly sensitive to microsecond motions, for example, via R_1ρ measurements, it is challenging to obtain accurate dynamics on the microsecond timescale from MD simulations. Although MD simulations now regularly extend for multiple microseconds, or longer via enhanced sampling (Bernardi et al., 2015), one still lacks sufficient statistics to obtain reliable dynamics behavior. Consider, if we investigate a 1 μs motion, using a 10 μs trajectory, we should observe 10 events, but the variance in number of events is also 10 (assuming Poisson statistics), so that large errors easily occur. Additionally, correct replication of slower motions requires all the faster motion leading up to the slow motion to occur at approximately the correct rates, so that the slower motions are more susceptible to influences like force field inaccuracies, starting structure of the system, etc. This remains a significant challenge for combining experiment in simulation, requiring creative solutions to take advantage of simulation where reproduction of experimental observables is poor.

With experimental validation or refinement of an MD simulation, one may analyze the simulation further, with improved confidence of the accuracy of the simulation. However, we want to use the simulation specifically to improve our interpretation of the experimental parameters. For example, we recently showed that it was possible to calculate the spatial correlation of motions within a given detector window between different residues in HET-s (218–289) fibrils (Smith et al., 2019b), using a modified iRED analysis (Prompers and Brüschweiler, 2001; Prompers and Brüschweiler, 2002). The result is that we could see that detector windows corresponding to longer correlation times tended to result in correlation over longer distances, providing at least some explanation for the presence of slow, low amplitude motion in fibrils. A similar correlation analysis is shown in Figure 14G, in this case for residues in an α-helix, where similarly, detector windows corresponding to longer correlation times yield longer correlation lengths. We suspect this behavior to be nearly universal: even in well-defined structures, there is always some residual flexibility. Then, both short- and long-range modes of motion should be thermally populated (in terms of modes, these are more accurately described as having short and long wavelengths). However, the longer-range modes usually have longer correlation times, resulting in the trends in Figure 14G. Note that this implies that there should almost always be distributions of correlation times due to varying correlation length, further complicating the interpretation of the two to three correlation times provided by the EMF approach.

For fairly rigid regions of a molecule, we expect detector-specific correlation analysis to help explain dynamic trends. However, what should we do for regions that are more mobile, with multiple types of motion contributing? Having all of the atom positions in an MD simulation should provide the detail that would allow us to separate different motions. Then, we could define the total motion of a bond as resulting from the product of these motions. For example, for an H–N dipole coupling in an α-helix, the total rotation of the dipole is the result of the reorientation of the principal axis system (PAS) of the dipole within the peptide plane (PP), the peptide plane reorienting within the helix, the helix reorienting with the molecule, and the molecule reorienting within the lab frame.

$\begin{array}{l}\overrightarrow{v}\left(t+\tau \right)=\mathbf{R}\left({\Omega }_{\tau ,t+\tau }\right)\cdot \overrightarrow{v}\left(\tau \right)\\ =\mathbf{R}\left({\Omega }_{\tau ,t+\tau }^{\text{mol}\text{.}\to \text{lab}}\right)\cdot \mathbf{R}\left({\Omega }_{\tau ,t+\tau }^{\text{helix}\to \text{mol}\text{.}}\right)\cdot \mathbf{R}\left({\Omega }_{\tau ,t+\tau }^{\text{PP}\to \text{helix}}\right)\cdot \mathbf{R}\left({\Omega }_{\tau ,t+\tau }^{\text{PAS}\to \text{PP}}\right)\cdot \overrightarrow{v}\left(\tau \right)\end{array}(38)$

This concept is illustrated in Figure 14H. In the case that it is possible to derive a correlation function from each rotation, one then may effectively achieve an in silico model-free type separation of the correlation functions motion. A similar approach for the specific separation of librations, φ/ψ reorientation, and peptide plane tumbling in intrinsically disorded proteins has been demonstrated by Salvi et al. (2017), however we find that it is possible to fully generalize this concept for separation of arbitrary definitions of independent motions (manuscript under revision, (Smith et al., 2021b)). Then, separated motions may also be analyzed with detectors, to determine how both experimental and simulated detector responses depend on both timescale and position in the molecule. Separation of motions could also be coupled with mode analyses such as iRED (Prompers and Brüschweiler, 2001, 2002) or principal component analysis (Amadei et al., 1993; Altis et al., 2007), providing a method to better characterize distributions of correlation times arising from different motions and complex mode-like dynamics. In each proposed case, comparison of the different MD analyses is possible via the detector analysis. Our eventual goal is that one may extract enough detail from the MD to build explicit models of motion for direct application to the NMR experimental results, so that the final characterizations are no longer model-free at all, but rather yield highly detailed models based on the combined information from experiment and simulation.

Conclusion

We show that the original model-free approach, SDM, LeMaster’s approach, and detectors all belong to a class of methods where fit parameters are resulting from a linear combination of experimental relaxation rate constants (potentially requiring an additional arithmetic step to yield the final parameters). IMPACT is a close approximation to this behavior, whereas EMF parameters exhibit significantly different behavior. Analysis methods belonging to this class are particularly useful because it is straightforward to estimate the resulting parameters if the distribution of correlation times, $\left(1-S^2\right)\theta \left(z\right)$, is known. This is particularly advantageous when determining if a model is consistent with experimentally determined parameters, and also allows easy comparison of multiple methods.

The detector analysis is the most general of these approaches, being applicable to any collection of NMR relaxation experiments probing reorientational motion, and can be generalized for other methods such as MD simulation, requiring very little modification of the analysis. Then, the resulting detector responses from NMR and MD are easily compared. With experimental validation of MD, one may then use the wealth of detail in MD simulation to better understand how experimentally derived parameters are related to specific motion, via correlation of motion, separation of motion, and other existing and yet-to-be developed techniques. This has the potential to lead to improved models of motion for NMR analysis, which in turn can help obtain a more fundamental understand of dynamics in biomolecular systems.

Author Contributions

AS has prepared the manuscript text and KZ has prepared the figures.

Funding

The authors acknowledge funding by the Deutsche Forschungsgemeinschaft (DFG) grant SM 576/1-1 and by European Social Funds (ESF) and the Free State of Saxony (Junior Research Group UniDyn, Project No. SAB 100382164). The authors also acknowledge support from the German Research Foundation (DFG) and Universität Leipzig within the program of Open Access Publishing.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmolb.2021.727553/full#supplementary-material

References

Abergel, D., and Palmer, A. G. (2003). On the Use of the Stochastic Liouville Equation in Nuclear Magnetic Resonance: Application to R1ρ Relaxation in the Presence of Exchange. Concepts Magn. Reson. 19A, 134–148. doi:10.1002/cmr.a.10091

CrossRef Full Text | Google Scholar

Altis, A., Nguyen, P. H., Hegger, R., and Stock, G. (2007). Dihedral Angle Principal Component Analysis of Molecular Dynamics Simulations. J. Chem. Phys. 126, 244111. doi:10.1063/1.2746330

CrossRef Full Text | Google Scholar

Amadei, A., Linssen, A. B. M., and Berendsen, H. J. C. (1993). Essential Dynamics of Proteins. Proteins 17, 412–425. doi:10.1002/prot.340170408

PubMed Abstract | CrossRef Full Text | Google Scholar

Anderson, M., Motta, R., Chandrasekar, S., and Stokes, M. (1996). “Proposal for a Standard Default Color Space for the Internet—sRGB,” in Color and Imaging Conference (Society for Imaging Science and Technology), 238–245.

Google Scholar

Antila, H. S., Ferreira, T. M., Ollila, O. H. S., and Miettinen, M. S. (2021). Using Open Data to Rapidly Benchmark Biomolecular Simulations: Phospholipid Conformational Dynamics. J. Chem. Inf. Model. 61, 938–949. doi:10.1021/acs.jcim.0c01299

CrossRef Full Text | Google Scholar

Beckmann, P. A. (1988). Spectral Densities and Nuclear Spin Relaxation in Solids. Phys. Rep. 171, 85–128. doi:10.1016/0370-1573(88)90073-7

CrossRef Full Text | Google Scholar

Bernardi, R. C., Melo, M. C. R., and Schulten, K. (2015). Enhanced Sampling Techniques in Molecular Dynamics Simulations of Biological Systems. Biochim. Biophys. Acta (Bba) - Gen. Subjects 1850, 872–877. doi:10.1016/j.bbagen.2014.10.019

CrossRef Full Text | Google Scholar

Brown, M. (1982). Theory of Spin-Lattice Relaxation in Lipid Bilayers and Biological Membranes. 2H and 14N Quadrupolar Relaxation. J. Chem. Phys. 77 (3), 1576–1599. doi:10.1063/1.443940

CrossRef Full Text | Google Scholar

Brown, M. F. (1996). “Biological Membranes,” in Biological Membranes: A Molecular Perspective From Computation And Experiment. Editors K. M. Merz, and B. Roux (Birkhäuser Basel). doi:10.1007/978-1-4684-8580-6

CrossRef Full Text | Google Scholar

Chevelkov, V., Fink, U., and Reif, B. (2009a). Accurate Determination of Order Parameters from 1H,15N Dipolar Couplings in MAS Solid-State NMR Experiments. J. Am. Chem. Soc. 131, 14018–14022. doi:10.1021/ja902649u

CrossRef Full Text | Google Scholar

Chevelkov, V., Fink, U., and Reif, B. (2009b). Quantitative Analysis of Backbone Motion in Proteins Using MAS Solid-State NMR Spectroscopy. J. Biomol. NMR 45, 197–206. doi:10.1007/S10858-009-9348-5

CrossRef Full Text | Google Scholar

Clore, G. M., Szabo, A., Bax, A., Kay, L. E., Driscoll, P. C., and Gronenborn, A. M. (1990). Deviations from the Simple Two-Parameter Model-free Approach to the Interpretation of Nitrogen-15 Nuclear Magnetic Relaxation of Proteins. J. Am. Chem. Soc. 112, 4989–4991. doi:10.1021/ja00168a070

CrossRef Full Text | Google Scholar

d'Auvergne, E. J., and Gooley, P. R. (2003). The Use of Model Selection in the Model-free Analysis of Protein Dynamics. J. Biomol. NMR 25, 25–39. doi:10.1023/A:1021902006114

CrossRef Full Text | Google Scholar

Dantzig, G. B. (1982). Reminiscences about the Origins of Linear Programming. Operations Res. Lett. 1, 43–48. doi:10.1016/0167-6377(82)90043-8

CrossRef Full Text | Google Scholar

Gill, M. L., Byrd, R. A., and Palmer, A. G. (2016). Dynamics of GCN4 Facilitate DNA Interaction: a Model-free Analysis of an Intrinsically Disordered Region. Phys. Chem. Chem. Phys. 18, 5839–5849. doi:10.1039/c5cp06197k

PubMed Abstract | CrossRef Full Text | Google Scholar

Golub, G., and Kahan, W. (1965). Calculating the Singular Values and Pseudo-inverse of a Matrix. J. Soc. Ind. Appl. Maths. Ser. B Numer. Anal. 2, 205–224. doi:10.1137/0702016

CrossRef Full Text | Google Scholar

Gullion, T., and Schaefer, J. (1989). Rotational-Echo Double-Resonance NMR. J. Magn. Reson. (1969) 81, 196–200. doi:10.1016/0022-2364(89)90280-1

CrossRef Full Text | Google Scholar

Halle, B., and Wennerström, H. (1981). Interpretation of Magnetic Resonance Data from Water Nuclei in Heterogeneous Systems. J. Chem. Phys. 75, 1928–1943. doi:10.1063/1.442218

CrossRef Full Text | Google Scholar

Haller, J. D., and Schanda, P. (2013). Amplitudes and Time Scales of Picosecond-To-Microsecond Motion in Proteins Studied by Solid-State NMR: a Critical Evaluation of Experimental Approaches and Application to Crystalline Ubiquitin. J. Biomol. NMR 57, 263–280. doi:10.1007/s10858-013-9787-x

CrossRef Full Text | Google Scholar

Hansen, N., Heller, F., Schmid, N., and van Gunsteren, W. F. (2014). Time-averaged Order Parameter Restraints in Molecular Dynamics Simulations. J. Biomol. NMR 60, 169–187. doi:10.1007/s10858-014-9866-7

CrossRef Full Text | Google Scholar

Hoffmann, F., Mulder, F. A. A., and Schäfer, L. V. (2018a). Accurate Methyl Group Dynamics in Protein Simulations with AMBER Force Fields. J. Phys. Chem. B. 122, 5038–5048. doi:10.1021/acs.jpcb.8b02769

PubMed Abstract | CrossRef Full Text | Google Scholar

Hoffmann, F., Xue, M., Schäfer, L. V., and Mulder, F. A. A. (2018b). Narrowing the gap between Experimental and Computational Determination of Methyl Group Dynamics in Proteins. Phys. Chem. Chem. Phys. 20, 24577–24590. doi:10.1039/c8cp03915a

PubMed Abstract | CrossRef Full Text | Google Scholar

Hsu, A., O'Brien, P. A., Bhattacharya, S., Rance, M., and Palmer, A. G. (2018). Enhanced Spectral Density Mapping through Combined Multiple-Field Deuterium 13CH2D Methyl Spin Relaxation NMR Spectroscopy. Methods 138-139, 76–84. doi:10.1016/j.ymeth.2017.12.020

PubMed Abstract | CrossRef Full Text | Google Scholar

Ishima, R., Louis, J. M., and Torchia, D. A. (1999). Transverse 13C Relaxation of CHD2 Methyl Isotopmers to Detect Slow Conformational Changes of Protein Side Chains. J. Am. Chem. Soc. 121, 11589–11590. doi:10.1021/ja992836b

CrossRef Full Text | Google Scholar

Judd, D. B. (1951). Report of U.S. Secretariat Committee on Colorimetry and Artificial Daylight, Proceedings of the Twelfth Session of the CIE. Stockholm: Paris: Bureau Central de la CIE, 11.

Google Scholar

Kantorovich, L. (1960). On the Calculation of Production Inputs. Probl. Econ. 3, 3–10. doi:10.2753/PET1061-199103013

CrossRef Full Text | Google Scholar

Kay, L. E., Torchia, D. A., and Bax, A. (1989). Backbone Dynamics of Proteins as Studied by Nitrogen-15 Inverse Detected Heteronuclear NMR Spectroscopy: Application to Staphylococcal Nuclease. Biochemistry 28, 8972–8979. doi:10.1021/bi00449a003

PubMed Abstract | CrossRef Full Text | Google Scholar

Keeler, E. G., Fritzsching, K. J., and McDermott, A. E. (2018). Refocusing CSA during Magic Angle Spinning Rotating-Frame Relaxation Experiments. J. Magn. Reson. 296, 130–137. doi:10.1016/j.jmr.2018.09.004

CrossRef Full Text | Google Scholar

Khan, S. N., Charlier, C., Augustyniak, R., Salvi, N., Déjean, V., Bodenhausen, G., et al. (2015). Distribution of Pico- and Nanosecond Motions in Disordered Proteins from Nuclear Spin Relaxation. Biophysical J. 109, 988–999. doi:10.1016/j.bpj.2015.06.069

PubMed Abstract | CrossRef Full Text | Google Scholar

Kinosita, K., Kawato, S., and Ikegami, A. (1977). A Theory of Fluorescence Polarization Decay in Membranes. Biophysical J. 20, 289–305. doi:10.1016/S0006-3495(77)85550-1

CrossRef Full Text | Google Scholar

Korzhnev, D. M., Neudecker, P., Mittermaier, A., Orekhov, V. Y., and Kay, L. E. (2005). Multiple-Site Exchange in Proteins Studied with a Suite of Six NMR Relaxation Dispersion Experiments: An Application to the Folding of a Fyn SH3 Domain Mutant. J. Am. Chem. Soc. 127, 15602–15611. doi:10.1021/ja054550e

CrossRef Full Text | Google Scholar

Koss, H., Rance, M., and Palmer, A. G. (2018). General Expressions for Carr-Purcell-Meiboom-Gill Relaxation Dispersion for N-Site Chemical Exchange. Biochemistry 57, 4753–4763. doi:10.1021/acs.biochem.8b00370

PubMed Abstract | CrossRef Full Text | Google Scholar

Krushelnitsky, A., Gauto, D., Rodriguez Camargo, D. C., Schanda, P., and Saalwächter, K. (2018). Microsecond Motions Probed by Near-Rotary-Resonance R1ρ 15N MAS NMR Experiments: the Model Case of Protein Overall-Rocking in Crystals. J. Biomol. NMR 71, 53–67. doi:10.1007/s10858-018-0191-4

CrossRef Full Text | Google Scholar

Kurauskas, V., Izmailov, S. A., Rogacheva, O. N., Hessel, A., Ayala, I., Woodhouse, J., et al. (2017). Slow Conformational Exchange and Overall Rocking Motion in Ubiquitin Protein Crystals. Nat. Commun. 8, 145. doi:10.1038/s41467-017-00165-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Kurbanov, R., Zinkevich, T., and Krushelnitsky, A. (2011). The Nuclear Magnetic Resonance Relaxation Data Analysis in Solids: General R1/R1ρ Equations and the Model-free Approach. J. Chem. Phys. 135 (1–9), 184104. doi:10.1063/1.3658383

CrossRef Full Text | Google Scholar

Lakomek, N.-A., Penzel, S., Lends, A., Cadalbert, R., Ernst, M., and Meier, B. H. (2017). Microsecond Dynamics in Ubiquitin Probed by Solid-State 15N NMR Spectroscopy R1ρ Relaxation Experiments under Fast Mas (60-110 Khz). Chem. Eur. J. 23, 9425–9433. doi:10.1002/chem.201701738

PubMed Abstract | CrossRef Full Text | Google Scholar

Lamley, J. M., Lougher, M. J., Sass, H. J., Rogowski, M., Grzesiek, S., and Lewandowski, J. R. (2015a). Unraveling the Complexity of Protein Backbone Dynamics with Combined 13C and 15N Solid-State NMR Relaxation Measurements. Phys. Chem. Chem. Phys. 17, 21997–22008. doi:10.1039/c5cp03484a

PubMed Abstract | CrossRef Full Text | Google Scholar

Lamley, J. M., Öster, C., Stevens, R. A., and Lewandowski, J. R. (2015b). Intermolecular Interactions and Protein Dynamics by Solid‐State NMR Spectroscopy. Angew. Chem. Int. Ed. 54, 15374–15378. doi:10.1002/anie.201509168

PubMed Abstract | CrossRef Full Text | Google Scholar

LeMaster, D. (1995). Larmor Frequency Selective Model Free Analysis of Protein NMR Relaxation. J. Biomol. NMR 6, 366–374. doi:10.1007/BF00197636

CrossRef Full Text | Google Scholar

Lipari, G., and Szabo, A. (1982a). Model-free Approach to the Interpretation of Nuclear Magnetic Resonance Relaxation in Macromolecules. 1. Theory and Range of Validity. J. Am. Chem. Soc. 104, 4546–4559. doi:10.1021/ja00381a009

CrossRef Full Text | Google Scholar

Lipari, G., and Szabo, A. (1982b). Model-free Approach to the Interpretation of Nuclear Magnetic Resonance Relaxation in Macromolecules. 2. Analysis of Experimental Results. J. Am. Chem. Soc. 104, 4559–4570. doi:10.1021/ja00381a010

CrossRef Full Text | Google Scholar

Ma, P., Haller, J. D., Zajakala, J., Macek, P., Sivertsen, A. C., Willbold, D., et al. (2014). Probing Transient Conformational States of Proteins by Solid-State R1ρ Relaxation-Dispersion NMR Spectroscopy. Angew. Chem. Int. Ed. 53, 4312–4317. doi:10.1002/anie.201311275

CrossRef Full Text | Google Scholar

Ma, P., Xue, Y., Coquelle, N., Haller, J. D., Yuwen, T., Ayala, I., et al. (2015). Observing the Overall Rocking Motion of a Protein in a Crystal. Nat. Commun. 6, 8361. doi:10.1038/ncomms9361

PubMed Abstract | CrossRef Full Text | Google Scholar

Mandel, A. M., Akke, M., and Palmer, A. G. (1995). Backbone Dynamics of Escherichia coli Ribonuclease HI: Correlations with Structure and Function in an Active Enzyme. J. Mol. Biol. 246, 144–163. doi:10.1006/jmbi.1994.0073

CrossRef Full Text | Google Scholar

Marion, D., Gauto, D. F., Ayala, I., Giandoreggio-Barranco, K., and Schanda, P. (2019). Microsecond Protein Dynamics from Combined Bloch-McConnell and Near-Rotary-Resonance R1p Relaxation-Dispersion MAS NMR. ChemPhysChem 20, 276–284. doi:10.1002/cphc.201800935

PubMed Abstract | CrossRef Full Text | Google Scholar

Mendelman, N., and Meirovitch, E. (2021). Structural Dynamics from NMR Relaxation by SRLS Analysis: Local Geometry, Potential Energy Landscapes, and Spectral Densities. J. Phys. Chem. B. 125, 6130–6143. doi:10.1021/acs.jpcb.1c02502

CrossRef Full Text | Google Scholar

Miloushev, V. Z., and Palmer, A. G. (2005). R1ρ Relaxation for Two-Site Chemical Exchange: General Approximations and Some Exact Solutions. J. Magn. Reson. 177, 221–227. doi:10.1016/j.jmr.2005.07.023

CrossRef Full Text | Google Scholar

Mollica, L., Baias, M., Lewandowski, J. R., Wylie, B. J., Sperling, L. J., Rienstra, C. M., et al. (2012). Atomic-Resolution Structural Dynamics in Crystalline Proteins from NMR and Molecular Simulation. J. Phys. Chem. Lett. 3, 3657–3662. doi:10.1021/Jz3016233

PubMed Abstract | CrossRef Full Text | Google Scholar

Munowitz, M. G., Griffin, R. G., Bodenhausen, G., and Huang, T. H. (1981). Two-dimensional Rotational Spin-echo Nuclear Magnetic Resonance in Solids: Correlation of Chemical Shift and Dipolar Interactions. J. Am. Chem. Soc. 103, 2529–2533. doi:10.1021/ja00400a007

CrossRef Full Text | Google Scholar

Neudecker, P., Korzhnev, D. M., and Kay, L. E. (2006). Assessment of the Effects of Increased Relaxation Dispersion Data on the Extraction of 3-site Exchange Parameters Characterizing the Unfolding of an SH3 Domain. J. Biomol. NMR 34, 129–135. doi:10.1007/s10858-006-0001-2

CrossRef Full Text | Google Scholar

Öster, C., Kosol, S., and Lewandowski, J. R. (2019). Quantifying Microsecond Exchange in Large Protein Complexes with Accelerated Relaxation Dispersion Experiments in the Solid State. Sci. Rep. 9, 11082. doi:10.1038/s41598-019-47507-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Palmer, A. G. (2004). NMR Characterization of the Dynamics of Biomacromolecules. Chem. Rev. 104, 3623–3640. doi:10.1021/cr030413t

PubMed Abstract | CrossRef Full Text | Google Scholar

Peng, J. W., and Wagner, G. (1992). Mapping of Spectral Density Functions Using Heteronuclear NMR Relaxation Measurements. J. Magn. Reson. (1969) 98, 308–332. doi:10.1016/0022-2364(92)90135-T

CrossRef Full Text | Google Scholar

Penzel, S., Smith, A. A., Agarwal, V., Hunkeler, A., Org, M.-L., Samoson, A., et al. (2015). Protein Resonance Assignment at MAS Frequencies Approaching 100 kHz: a Quantitative Comparison of J-Coupling and Dipolar-Coupling-Based Transfer Methods. J. Biomol. NMR 63, 165–186. doi:10.1007/s10858-015-9975-y

CrossRef Full Text | Google Scholar

Pettersen, E. F., Goddard, T. D., Huang, C. C., Meng, E. C., Couch, G. S., Croll, T. I., et al. (2021). UCSF ChimeraX : Structure Visualization for Researchers, Educators, and Developers. Protein Sci. 30, 70–82. doi:10.1002/pro.3943

PubMed Abstract | CrossRef Full Text | Google Scholar

Polimeno, A., and Freed, J. H. (1992). Advances In Chemical Physics. John Wiley & Sons, 89–206. doi:10.1002/9780470141410.ch3A Many-Body Stochastic Approach to Rotational Motions in Liquids.

CrossRef Full Text | Google Scholar

Prompers, J. J., and Brüschweiler, R. (2002). General Framework for Studying the Dynamics of Folded and Nonfolded Proteins by NMR Relaxation Spectroscopy and MD Simulation. J. Am. Chem. Soc. 124, 4522–4534. doi:10.1021/ja012750u

CrossRef Full Text | Google Scholar

Prompers, J. J., and Brüschweiler, R. (2001). Reorientational Eigenmode Dynamics: A Combined MD/NMR Relaxation Analysis Method for Flexible Parts in Globular Proteins. J. Am. Chem. Soc. 123, 7305–7313. doi:10.1021/ja0107226

CrossRef Full Text | Google Scholar

Rovó, P., and Linser, R. (2017). Microsecond Time Scale Proton Rotating-Frame Relaxation under Magic Angle Spinning. J. Phys. Chem. B. 121, 6117–6130. doi:10.1021/acs.jpcb.7b03333

CrossRef Full Text | Google Scholar

Rovó, P. (2020). Recent Advances in Solid-State Relaxation Dispersion Techniques. Solid State. Nucl. Magn. Reson. 108, 101665. doi:10.1016/j.ssnmr.2020.101665

PubMed Abstract | CrossRef Full Text | Google Scholar

Salvi, N., Abyzov, A., and Blackledge, M. (2017). Analytical Description of NMR Relaxation Highlights Correlated Dynamics in Intrinsically Disordered Proteins. Angew. Chem. Int. Ed. 56, 14020–14024. doi:10.1002/anie.201706740

PubMed Abstract | CrossRef Full Text | Google Scholar

Salvi, N., Abyzov, A., and Blackledge, M. (2016). Multi-Timescale Dynamics in Intrinsically Disordered Proteins from NMR Relaxation and Molecular Simulation. J. Phys. Chem. Lett. 7, 2483–2489. doi:10.1021/acs.jpclett.6b00885

PubMed Abstract | CrossRef Full Text | Google Scholar

Schanda, P., and Ernst, M. (2016). Studying Dynamics by Magic-Angle Spinning Solid-State NMR Spectroscopy: Principles and Applications to Biomolecules. Prog. Nucl. Magn. Reson. Spectrosc. 96, 1–46. doi:10.1016/j.pnmrs.2016.02.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Schanda, P., Meier, B. H., and Ernst, M. (2011). Accurate Measurement of One-Bond H-X Heteronuclear Dipolar Couplings in MAS Solid-State NMR. J. Magn. Reson. 210, 246–259. doi:10.1016/j.jmr.2011.03.015

CrossRef Full Text | Google Scholar

Schanda, P., Meier, B. H., and Ernst, M. (2010). Quantitative Analysis of Protein Backbone Dynamics in Microcrystalline Ubiquitin by Solid-State NMR Spectroscopy. J. Am. Chem. Soc. 132, 15957–15967. doi:10.1021/Ja100726a

CrossRef Full Text | Google Scholar

Schanda, P. (2019). Relaxing with Liquids and Solids - A Perspective on Biomolecular Dynamics. J. Magn. Reson. 306, 180–186. doi:10.1016/j.jmr.2019.07.025

CrossRef Full Text | Google Scholar

Shapiro, Y. E., and Meirovitch, E. (2012). Slowly Relaxing Local Structure (SRLS) Analysis of 15N-H Relaxation from the Prototypical Small Proteins GB1 and GB3. J. Phys. Chem. B. 116, 4056–4068. doi:10.1021/jp300245k

CrossRef Full Text | Google Scholar

Skrynnikov, N. R., Millet, O., and Kay, L. E. (2002). Deuterium Spin Probes of Side-Chain Dynamics in Proteins. 2. Spectral Density Mapping and Identification of Nanosecond Time-Scale Side-Chain Motions. J. Am. Chem. Soc. 124, 6449–6460. doi:10.1021/ja012498q

CrossRef Full Text | Google Scholar

Smith, A. A., Bolik-Coulon, N., Ernst, M., Meier, B. H., and Ferrage, F. (2021a). How Wide Is the Window Opened by High-Resolution Relaxometry on the Internal Dynamics of Proteins in Solution? J. Biomol. NMR 75, 119–131. doi:10.1007/s10858-021-00361-1

CrossRef Full Text | Google Scholar

Smith, A. A., Ernst, M., and Meier, B. H. (2017). Because the Light Is Better Here: Correlation-Time Analysis by NMR Spectroscopy. Angew. Chem. 129, 13778–13783. doi:10.1002/ange.201707316

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, A. A., Ernst, M., Meier, B. H., and Ferrage, F. (2019a). Reducing Bias in the Analysis of Solution-State NMR Data with Dynamics Detectors. J. Chem. Phys. 151, 034102. doi:10.1063/1.5111081

CrossRef Full Text | Google Scholar

Smith, A. A., Ernst, M., and Meier, B. H. (2018). Optimized "detectors" for Dynamics Analysis in Solid-State NMR. J. Chem. Phys. 148, 045104. doi:10.1063/1.5013316

CrossRef Full Text | Google Scholar

Smith, A. A., Ernst, M., Riniker, S., and Meier, B. H. (2019b). Localized and Collective Motions in HET‐s(218‐289) Fibrils from Combined NMR Relaxation and MD Simulation. Angew. Chem. 131, 9483–9488. doi:10.1002/ange.201901929

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, A. A., Testori, E., Cadalbert, R., Meier, B. H., and Ernst, M. (2016). Characterization of Fibril Dynamics on Three Timescales by Solid-State NMR. J. Biomol. NMR 65, 171–191. doi:10.1007/s10858-016-0047-8

CrossRef Full Text | Google Scholar

Smith, A. A., Vogel, A., Engberg, O., Hildebrand, P. W., and Huster, D. (2021b). A Method to Construct the Dynamic Landscape of a Bio-Membrane with experiment and Simulation. ResearchSquare, 1–28. doi:10.21203/rs.3.rs-645823/v1

CrossRef Full Text | Google Scholar

Smith, T., and Guild, J. (1931). The C.I.E. Colorimetric Standards and Their Use. Trans. Opt. Soc. 33, 73–134. doi:10.1088/1475-4878/33/3/301

CrossRef Full Text | Google Scholar

Tikhonov, A. N., and Arsenin, V. I. A. (1977). Solutions of Ill-Posed Problems. Washington, DC: Winston and Sons.

Google Scholar

Trott, O., Abergel, D., and Palmer, A. G. (2003). An Average-Magnetization Analysis of R1ρ Relaxation outside of the Fast Exchange Limit. Mol. Phys. 101, 753–763. doi:10.1080/0026897021000054826

CrossRef Full Text | Google Scholar

Trott, O., and Palmer, A. G. (2002). R1ρ Relaxation outside of the Fast-Exchange Limit. J. Magn. Reson. 154, 157–160. doi:10.1006/jmre.2001.2466

CrossRef Full Text | Google Scholar

Tugarinov, V., Liang, Z., Shapiro, Y. E., Freed, J. H., and Meirovitch, E. (2001). A Structural Mode-Coupling Approach to 15N NMR Relaxation in Proteins. J. Am. Chem. Soc. 123, 3055–3063. doi:10.1021/ja003803v

CrossRef Full Text | Google Scholar

Virtanen, P. (2020). SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python. Nat. Methods 17, 261–272. doi:10.1038/s41592-019-0686-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Vos, J. J. (1978). Colorimetric and Photometric Properties of a 2° Fundamental Observer. Color Res. Appl. 3, 125–128. doi:10.1002/col.5080030309

CrossRef Full Text | Google Scholar

Wennerström, H., Lindman, B., Soederman, O., Drakenberg, T., and Rosenholm, J. B. (1979). Carbon-13 Magnetic Relaxation in Micellar Solutions. Influence of Aggregate Motion on T1. J. Am. Chem. Soc. 101, 6860–6864. doi:10.1021/ja00517a012

CrossRef Full Text | Google Scholar

Zinkevich, T., Chevelkov, V., Reif, B., Saalwächter, K., and Krushelnitsky, A. (2013). Internal Protein Dynamics on ps to μs Timescales as Studied by Multi-Frequency 15N Solid-State NMR Relaxation. J. Biomol. NMR 57, 219–235. doi:10.1007/s10858-013-9782-2

CrossRef Full Text | Google Scholar