PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
J Magn Reson. Author manuscript; available in PMC 2010 November 1.
Published in final edited form as:
PMCID: PMC2763024
NIHMSID: NIHMS142117

Improvement and Analysis of Computational Methods for Prediction of Residual Dipolar Couplings

Abstract

We describe a new, computationally effcient method for computing the molecular alignment tensor based on the molecular shape. The increase in speed is achieved by re-expressing the problem as one of numerical integration, rather than a simple uniform sampling (as in the PALES method), and by using a convex hull rather than a detailed representation of the surface of a molecule. This method is applicable to bicelles, PEG/hexanol, and other alignment media that can be modeled by steric restrictions introduced by a planar barrier. This method is used to further explore and compare various representations of protein shape by an equivalent ellipsoid. We also examine the accuracy of the alignment tensor and residual dipolar couplings (RDC) prediction using various ab initio methods. We separately quantify the inaccuracy in RDC prediction caused by the inaccuracy in the orientation and in the magnitude of the alignment tensor, concluding that orientation accuracy is much more important in accurate prediction of RDCs.

Keywords: residual dipolar coupling, alignment tensor, ab initio prediction, PALES

1. Introduction

Knowledge of protein structure plays a critical role in our understanding of the molecular mechanisms underlying biological processes. One of the main methods for obtaining structural information at atomic-level resolution is the use of nuclear magnetic resonance (NMR) spectroscopy for determining structural constraints. The NMR-derived constraints, such as NOEs, hydrogen bonds, and torsion angles, are intrinsically local or short-range and could be insuffcient for accurate structure determination of biological macromolecules and their complexes due to the scarcity of long-distance structural information. Residual dipolar couplings (RDCs), resulting from partial alignment of solute molecules relative to the magnetic field, provide valuable structural information in terms of global, long-range orientational constraints [1]. A commonly used method for aligning molecules in solution takes advantage of the anisotropy of molecular shape by imposing steric restrictions on the allowed orientations of the molecule (e.g. by means of bicelles [2], stretched gels [3, 4], or PEG/hexanol-based media [5]). Such steric alignment can often be modeled as caused by planar obstacles, and we will refer to this simplified model of molecular alignment as the barrier model. The alignment of a rigid molecule can be described by the so-called molecular alignment tensor. Accurate prediction of the molecular alignment tensor, and with it of the RDCs, is important for NMR-based structure determination and validation as well as applications to dynamic and disordered systems (see e.g., [6, 7, 8, 9, 10, 11]). The sensitivity to molecular shape has the potential for improving structure characterization, especially in multidomain systems and macromolecular complexes (e.g., [12]), by fully integrating RDC prediction into structure refinement protocols to directly drive structure optimization. Future progress in this direction critically depends on the effciency and accuracy of the alignment tensor prediction.

Several methods for computing the molecular alignment tensor ab initio, i.e. based solely on the three-dimensional shape of the molecule, have recently been proposed. In the method by Zweckstetter and Bax [13, 14], implemented in a program called PALES, the alignment tensor is computed by uniformly sampling all orientations of a molecule (see e.g., [15]) at various distances away from a planar barrier, and averaging over only those orientations in which the molecule’s surface does not collide with the barrier. The computational effciency of this method is limited due to the fact that it must compute collisions between an arbitrary shape and a plane for every sample in the four-dimensional problem space.

Simpler methods based on the barrier model, but representing the shape of the molecule by an equivalent ellipsoid, have also been proposed. In Fernandes et al. [16], the alignment tensor was computed by approximating the molecule as an axially-symmetric prolate ellipsoid and analytically solving the barrier model for the alignment tensor. In Almond et al. [17] and Azurmendi et al. [18] the barrier method is also used, but the formulae are derived empirically.

Here we describe a new, computationally effcient method for computing the molecular alignment tensor based on the barrier model. The increase in speed is achieved by re-expressing the problem as one of numerical integration, rather than one of simple uniform sampling. This formulation allowed us to simplify the problem by reducing its dimensionality from four to two. In addition to the reduction in computational complexity, numerical integration has the advantage of (i) allowing control over the size of numerical error, and (ii) allowing a more effcient sampling of the problem space [19]. Computational geometry techniques are used to increase the computational speed further. We will refer to our method as PATI (Prediction of Alignment Tensor using Integration). PATI can also be used with an equivalent ellipsoid of the molecule instead of the full surface. We will refer to this simplification as PATI-E. This simplified method is used to explore and compare various representations of protein shape by a (fully anisotropic) equivalent ellipsoid: based on the gyration tensor [16], the actual molecular surface [20], or the minimum-volume ellipsoid.

Finally, we examine the accuracy of the proposed methods (PATI and PATI-E) and the existing ab initio methods for RDC prediction. This analysis separately quantifies the effect of inaccuracy in the predicted RDCs caused by the inaccuracy in the orientation or in the magnitude of the alignment tensor. The results obtained for several proteins show that (i) the predicted RDCs and their agreement with experimental data are very sensitive to errors in orientation of the alignment tensor, and (ii) all ab initio prediction methods tested here give a rather crude estimate of the RDCs.

2. Theory

For a rigid molecule, the molecular alignment tensor A with respect to the magnetic field B is described by a 3 × 3 symmetric traceless matrix [13], sometimes referred to as the Saupe matrix [21], with the following elements (i, j = 1, 2, 3):

Aij=12<Fij>,Fij=3cosθicosθjδij,
(1)

where θi is the angle between molecular axis i and the magnetic field B, < … > is the average over all possible orientations of the molecule in solution, and δ is the Kronecker delta.

The RDC value DPQ for a specific bond PQ is related to the alignment tensor and the bond’s orientation relative to the molecule’s coordinate frame by the following equation:

DPQ=CPQi,jAijcosφicosφj,CPQ=SLSμ0γPγQ4π2rPQ3,
(2)

where [var phi]i is the angle between the PQ bond and the molecular axis i, SLS is the Lipari-Szabo generalized order parameter, μ0 is the permeability of free space, γP and γQ are the gyromagnetic ratios of the corresponding nuclei, [variant Planck's over 2pi] is the reduced Planck’s constant, and rPQ is the length of the bond. PQ can represent bonds such as NH, CαHα, CαC′, and CN.

2.1. The Model for the Alignment Tensor

Given the 3D structure of an arbitrary molecule, we will focus on the computation of its alignment tensor A, defined in Equation (1). We model the planar barrier causing steric alignment of the molecule as a set of two infinite planes with the surface normals in the z direction, positioned at a distance 2h from each other. The molecule is centered around some point m (e.g., its center of mass), which lies somewhere inside the convex hull of the molecule’s surface. The direction of the magnetic field is given by a unit vector B, where

B=[b1b2b3],b12+b22+b32=1.
(3)

Figure 1 shows a schematic representation of the planar barrier model. Note that due to the symmetry of the system, the possible orientations of the molecule positioned between 0 and h along the z-axis are mirror images (over the xy plane) of the possible orientations when the molecule is between h and 2h. Thus we can simplify the model by considering only the bottom plane and positioning the molecule’s center at a height from 0 to h above this plane.

Figure 1
Planar barrier model for molecular alignment.

The orientation of the molecule’s coordinate frame relative to the Cartesian coordinate system in Figure 1 can be defined by three Euler angles α β, and γ, which determine the rotation matrix R(α, β, γ). (See Appendix A.)

For a specific molecule, we define S to be a finite set of sample points from its molecular surface (e.g. van der Waals surface or Richards molecular surface), and the center m of the molecule to be some point inside the convex hull of this molecular surface. Referring to Figure 1, to characterize the vertical extent of the molecule under the rotation R(α, β, γ) around its center m, we define η(α, β, γ), to be the difference between the z-coordinate of the center of the molecule and the minimum z-coordinate value of all the rotated points in S:

η(α,β,γ)=minskS{(R(α,β,γ)(skm))·[001]}.
(4)

Note that η(α, β, γ) sets the lower limit on the height of the center of the molecule at a given orientation.

We rewrite F, from Equation (1), in terms of the rotation matrix,

Fij(α,β,γ)=3(R1ib1+R2ib2+R3ib3)(R1jb1+R2jb2+R3jb3)δij,
(5)

and average the F values at height a with the mirror cases of 2ha into one equation

F¯ij(α,β,γ)=32(R1ib1+R2ib2+R3ib3)(R1jb1+R2jb2+R3jb3)+32(R1ib1+R2ib2R3ib3)(R1jb1+R2jb2R3jb3)δij,=3(R1ib1+R2ib2)(R1jb1+R2jb2)+3R3iR3jb32δij.
(6)

Due to the symmetry of the system, [F with macron] can be used instead of F to simplify our model to just one plane and a height from 0 to h.

For any rotation R(α, β, γ), the center of the molecule cannot be located at a height between 0 and η(α, β, γ). Therefore, the range of interest is from η(α, β, γ) to h. The alignment tensor A is then computed by summing [F with macron] weighted by the probability of the current height and orientation, for all allowed orientations and heights from η(α, β, γ) to h.

We assume equal a priori probabilities of all orientations at all heights, and write the analytical expression for Aij from Equation (1). To obtain a uniform distribution of the Euler angles we multiply our integrand by the Jacobian J = sin β/(8π2) [22] to obtain

Aij=12Nγ0γ1β0β1α0α1η(α,β,γ)hF¯ijJdzdαdβdγ,
(7)

where N is the normalization factor

N=γ0γ1β0β1α0α1η(α,β,γ)hJdzdαdβdγ,
(8)

and [α0, α1], [β0, β1], [γ0, γ1] are the ranges in which α,β, γ are defined.

2.2. Computation of the Alignment Tensor: General Case

We show in Appendix A, using the Euler z-y-z rotation, that the expression for the alignment tensor A, of an arbitrarily-shaped molecule can be simplified from a quadruple integral to a double integral. Specifically,

A11=Sc16Nπ02π11(3u2cos2α+3cos2α1)η(α,arccosu)dudα,A22=Sc16Nπ02π11(3u2cos2α3cos2α3u2+2)η(α,arccosu)dudα,A33=Sc16Nπ02π11(3u21)η(α,arccosu)dudα,A21=3Sc16Nπ02π11sinαcosα(1u2)η(α,arccosu)dudα,A31=3Sc16Nπ02π11u1u2cos(α)η(α,arccosu)dudα,A32=3Sc16Nπ02π11u1u2sin(α)η(α,arccosu)dudα,
(9)

where η is independent of the angle γ and

N=h14π02π11η(α,arccosu)dudα,Sc=13b32.
(10)

Because A is a traceless symmetric tensor [13], only A11 and A22, A21, A31, and A32 need to be computed, while A33 = −(A11 + A22), A12 = A21, A13 = A31, and A23 = A32. One can multiply the alignment tensor by −0.8 to account for the incomplete bicelle alignment the order parameter for the liquid crystal), and to match the sign returned by PALES. The height h can be determined by the formula d/(2Vf), where d is the barrier thickness (≈ 40Å for DMPC/DHPC bicelles) and Vf ([double less-than sign] 1) is the sample volume fraction occupied by the barriers (see [13, 14]).

Thus, all one needs to know in order to compute the alignment tensor is η(α, β), defined in Equation (4). Being an intrinsic geometric property of the molecule, η(α, β) can be computed separately, regardless of the barrier.

2.3. Computing η

In the PALES approach [13, 14], A is estimated based on forming a mesh of the molecular surface and then rotating all the mesh triangles of this surface to check if any part of the mesh is below the barrier. Observe that the complexity of each rotation is proportional to the number of triangles in the mesh. It is possible to simplify the mesh using mesh simplification (see [23, 24]); however even this is overly complex. An infinite planar barrier is not sensitive to cavities on the surface of the molecule; therefore, a convex hull of the molecule is a sufficient representation of the molecule’s surface. Additional mesh simplifications could be performed on the convex hull to further reduce the number of points.

To compute η for an arbitrary molecule under a rotation R, we simply compute the convex hull of the atom positions of the molecule and consider the vertices of the convex hull as the set S in Equation (4). We add the van der Waals radius of the atom associated with the minimum z-value to Equation (4) to form η for the rotation R. Figure 2A shows the convex hull around the Cyanovirin-N molecule. The number of points used to represent the molecule drops dramatically, from 40708 in the molecular surface representation (see [25, 26]), to just 57 in the convex hull representation. For any rotation R, the relative error in η between the two representations is less than 5% and the absolute error is less than 0.5Å. Also the alignment tensors and the RDCs predicted by PATI (our method) and PALES are almost identical, as shown below.

Figure 2
Convex hull and equivalent ellipsoids for Cyanovirin-N molecule drawn on top of its van der Waals surface. (A) The convex hull around the molecule. (B) The GE ellipsoid representation. (C) The MVE ellipsoid representation. (D) The PCAE ellipsoid representation ...

2.4. Special Case of an Ellipsoid

A potential simplification for computing the alignment tensor is to represent a molecule by an equivalent ellipsoid. In this section we examine several methods of deriving an equivalent ellipsoidal representation of an arbitrary molecule.

An ellipsoid ε in R3 is defined as

E(E,m)={p(pm)TE(pm)=1},
(11)

where E is a 3 × 3 symmetric positive definite matrix that defines the shape of the ellipsoid and m [set membership] R3 is its center.

The ellipsoid’s semi-principal axes can be derived by an eigendecomposition of E, such that

E=VΛVT=[V1V2V3][λ1000λ2000λ3][V1V2V3]T,
(12)

where the lengths of the semi-principal axes are a=1/λ1b=1/λ2c=1/λ3, and V1, V2, V3 are their associated directions (eigenvectors). The matrix V also represents the rotation matrix of the ellipsoid from the orientation where it is aligned to the principal axes of the coordinate frame to the actual orientation.

One method for finding an equivalent ellipsoid of a molecule is to find the minimum volume ellipsoid (MVE) that encloses all of the atoms in the molecule. Several methods exist for computing MVE [27, 28]. When we use the MVE in PATI, we denote the method as PATI-E.

A second method for constructing an equivalent ellipsoid of a molecule is to find an ellipsoid with the same gyration tensor as the molecule. The length of the ith axis of the equivalent ellipsoid based on the gyration tensor is computed as 5λi, where λi is the ith eigenvalue of the gyration tensor. We will refer to the ellipsoid with an equivalent moment of inertia as the Gyration Ellipsoid (GE) (see Fernandes et al. [16]).

A third method for finding an equivalent ellipsoid is to fit the molecular surface. This can be achieved via the Principal Component Analysis (PCA) of the coordinates of the points representing the surface of the molecule (see Ryabov et al. [20]). The length of the ith axis of the equivalent ellipsoid based on this method is 3λi, where λi is the ith eigenvalue of the corresponding covariance matrix. We refer to this method as PCAE.

Figure 2 illustrates the three ellipsoid models of the molecular surface for the Cyanovirin-N molecule.

As shown in Appendix B, for an arbitrary ellipsoid under an Euler rotation R(α, β, γ), η(α, β, γ) can be expressed as

η=R312a2+R322b2+R332c2.
(13)

3. Results

In this section we present a comprehensive comparison of several methods for computing RDCs ab initio. All the RDC data analyzed here are for the backbone NH bonds located in structurally well-defined regions of proteins, i.e. the α-helices and β-sheets. The RDC data were retrieved from the BMRB repository using the PDB code of the molecule. Only the RDC values measured using the neutral bicelle alignment medium (or, in the case of the B3 domain of protein G, the PEG/hexanol-based medium) are used. The 9 proteins and their codes in the Protein Data Bank are listed in Table 1.

Table 1
Quality Factors Q = Qs for the Experimental Data

We assess the quality of our results by computing the quality factor between the vector of experimental RDCs, Dexp, and our predicted RDCs for those same bonds, Dpred, as [29]:

Q=DexpDpred2Dexp2,
(14)

where ||·||2 is the l2-norm.

Note that the predicted magnitude of the RDC values depends on the experimental conditions (they determine the barrier height h) and selection of values for equation constants, e.g. CPQ. These factors affect all RDCs approximately uniformly, and hence can be repsesented by a scaling factor. Therefore, in order to make our analysis less sensitive to possible errors in experimental conditions and imperfect selection of values for constants, we also introduce the scaled quality factor to quantify the agreement between the experimental and predicted data with an unknown scaling factor. We define the scaled quality factor as

Qs=minρDexpρDpred2Dexp2,
(15)

where the scalar ρ can be computed by linear least squares. (Note that both PATI and PALES can predict the magnitude of RDC values with reasonable accuracy. See Supplementary Material for the values of ρ.)

First, we present the Q values for the experimental alignment tensor. We define the experimental alignment tensor, Ã, as the alignment tensor that optimally fits the data, i.e. gives the lowest Q value between the experimental and back-calculated RDC data. This quality factor allows us to examine whether the experimental data are well approximated by the theoretical equation for RDCs. We derive à by solving a linear least-squares problem of the form

[(v11)2(v31)2(v21)2(v31)22v11v2i2v11v312v21v31(v1i)2(v3i)2(v2i)2(v3i)22v1iv2i2v1iv3i2v2iv3i(v1n)2(v3n)2(v2n)2(v3n)22v1nv2n2v1nv3n2v2nv3n][A11A22A12A13A23][Dexp1/CNHDexpi/CNHDexpn/CNH],
(16)

where vi=[v1i,v2i,v3i] is the normalized vector representing the orientation of the ith bond relative to the molecular coordinate frame, n is the number of bonds, and CNH is the value of CPQ in Equation (2) for a NH bond. The linear least-squares problem can be solved by standard methods; see, e.g., [30].

Note that à can be decomposed into the experimental rotation (eigenvectors) V and experimental magnitudes (eigenvalues) Ã1, Ã2, Ã3, where

A=VΛVT=[V1V2V3][A1000A2000A3][V1V2V3]T.
(17)

The Q values for the experimental alignment tensor derived using Equation (16) are presented in Column 3 (‘LS’) in Table 1. The corresponding Q values for the best-fit alignment tensor derived from PALES are presented in Column 4, labeled PALES-LS. Naturally, Qs = Q for both methods. It is worth emphasizing here that this quality factor measures the actual quality of the experimental data (i.e. how well they fit the theoretical equation for RDCs) and therefore provides the baseline Q value for subsequent evaluation of the prediction methods. Note also that the values in the parentheses, the relative error in the alignment tensor, confirm that Equation (16) gives the same experimental alignment tensor à as the alignment tensor derived using PALES’s best-fit algorithm.

Table 1 shows that the experimental RDC data are of high quality and consistent with the theoretical formulation of the RDC (Equations (12)). This is not surprising given that these RDCs were used as constraints in the calculation/refinement of the corresponding protein structures. The quality of the agreement is illustrated in Figure 3A for Cyanovirin-N. (See also Supplementary Material.)

Figure 3
Comparison of the predicted vs. experimental 1H15N RDC values for the backbone amides in Cyanovirin-N, using various versions of the molecular alignment tensor derived from PATI. (A) The experimental alignment tensor was derived directly from the experimental ...

The results of our ab initio calculations are presented in Table 2, for PATI, PALES, and for the ellipsoidal approximation methods using the MVE model. The MVE ellipsoid data are used in this table as this model provides on average a slightly more accurate estimation of the alignment tensor compared to the other two equivalent ellipsoid models considered in this study. (See Supplementary Material.) Surprisingly, the scaled quality factor Qs was rather high for all prediction methods, indicating a generally marginal agreement with experimental data, as illustrated in Figure 3D and the Supplementary Material. PATI and PALES calculations gave on average a slightly better agreement with the data compared to the other methods. It should be emphasized here that PATI gives almost identical results to PALES, as evident from Figure 4. In order to understand the reasons for the observed inaccuracy in our predictions, we now break down the contributions to the errors into those due to the principal values of the alignment tensor and those due to inaccuracy in its orientation.

Figure 4
The agreement between RDC values predicted using PATI are those from PALES prediction. Shown are the 1H15N RDCs for all backbone amides for all molecules studied here. The (unscaled) quality factor Q between the two sets of RDC values is 0.05, the RMSD ...
Table 2
Quality factors Qs from RDC Prediction for ab initio Methods

In Table 3 we compare the Qs values for ‘synthetic’ alignment tensors that have the same orientation as the ab initio calculated tensors but the correct (experimental) principal values. We constructed these tensors by combining the rotation matrix V determined from our five models, GE, PCAE, MVE, PATI, and PALES, with the experimental magnitudes Ã1, Ã2, Ã3 of the alignment tensor derived from Ã. Such a comparison is expected to rank the methods based on the accuracy of prediction of the tensor’s orientation. Since the orientation of V for the equivalent-ellipsoid-based methods is derived directly from the orientation of the ellipsoid, this table also provides a direct comparison of the ellipsoid models. Note that there are six different combinations for V, since it is unknown a priori which Ãi is associated with which Vj. The smallest of the six Qs values is shown. Naturally, Qs = Q in this case.

Table 3
Quality of Prediction for the Orientation of Alignment Tensor

As evident from Table 3, correcting the principal values of the alignment tensor while keeping its predicted orientation did not improve the agreement with experimental data. (See also Figure 3B.) There are large variations among the various models in the accuracy of the predicted orientation of the alignment tensor. Of the three ellipsoid models tested here, MVE gave on average a somewhat better orientation (as documented in the Supplementary Material), while PATI and PALES yielded generally similar results.

We then constructed ‘synthetic’ alignment tensors that have the correct orientation (i.e. the V matrices derived from the experimental tensors Ã) but the same principal values (A1, A2, A3) as the ab initio calculated tensors. Table 4 displays the Qs values for five prediction methods, PATI, PALES, PATI-E, Almond, and PROLFIT. From this table, it is clear that using the correct orientation of the tensor dramatically improved the agreement with experimental data (cf. Table 2). This improvement is illustrated in Figure 3C for Cyanovirin-N and in the Supplementary Material for the other molecules.

Table 4
Quality of Prediction for the Magnitude of Alignment Tensor

Note that Column 3 (‘PATI-E’) and Column 5 (‘PROLFIT’) in Table 4 show that an additional degree of freedom provided by a fully anisotropic ellipsoid versus an axially-symmetric prolate ellipsoid approximation gives an improvement in the Qs.

Thus, the analysis presented above demonstrates that accurate prediction of the orientation of the alignment tensor is critical for the agreement with experimental data. Accurate prediction of the principal components of the tensor is important, too. However, when experimental RDCs are available, one can make an educated guess, based on the observed histogram/distribution of the data, about the magnitude of the tensor components (e.g., as described in [39]) and scale the predicted alignment tensor appropriately, whereas there is no obvious way to predict the orientation of the tensor.

4. Conclusions

We have reformulated the planar barrier model as a numerical integration problem and implemented it in a program called PATI. Our method has accuracy similar to PALES but is computationally more effcient and allows for finer control over numerical error. In addition, the convex hull provides a simpler representation of the surface, thus further increasing the computational effciency of the proposed method. This could allow PATI-based RDC prediction to be incorporated into the existing structure determination/refinement protocols. Because the molecular alignment tensor (and hence the RDC) is sensitive to the overall size and shape of the molecule, this would provide additional structural constraints that could potentially improve the accuracy of structure determination by NMR.

We compared several methods (old and new) for the computation of an equivalent ellipsoid of a molecule. We examined the accuracy of these equivalent ellipsoid models in predicting the alignment tensor and showed that the minimal volume ellipsoid gives on average a slightly better prediction of the alignment tensor orientation.

Finally, we compared all these methods against an extensive set of experimental RDC data. The analysis of the discrepancy between the experimental and predicted values emphasized the importance of the accurate prediction of the orientation of the alignment tensor. Possible sources of inaccuracy in ab initio alignment tensor prediction are the dynamic nature (structural flexibility) of protein molecules, not accounted for in the current prediction models, as well as the fact that the simple steric barrier model might not fully allow the correct alignment of all the molecules.

Supplementary Material

01

Acknowledgments

Supported by NIH grant GM065334 to D.F. The PATI program is available from the authors upon request.

Appendices

A. Derivation of Equation (9)

In this section we will simplify the expression for A for an arbitrary molecule from the quadruple integral of Equation (7) to a double integral.

Define the Euler z-y-z rotation matrix as R(α, β γ) = Rz(γ)Ry(β)Rz(α), where

Rz(α)=[cosαsinα0sinαcosα0001],Ry(β)=[cosβ0sinβ010sinβ0cosβ],Rz(γ)=[cosγsinγ0sinγcosγ0001].
(A1)

Multiplying the three matrices yields the full expression

R(α,β,γ)=[cosγcosβcosαsinγsinαcosγcosβsinαsinγcosαcosγsinβsinγcosβcosα+cosγsinαsinγcosβsinα+cosγcosαsinγsinβsinβcosαsinβsinαcosβ].
(A2)

We now write the equations for A11, A22, A33, A21, A31, A32, and N, recalling that the Jacobian is J = sin β/(8π2):

A11=116Nπ202π0π02πη(α,β)hF¯11sinβdzdαdβdγ,A22=116Nπ202π0π02πη(α,β)hF¯22sinβdzdαdβdγ,A33=116Nπ202π0π02πη(α,β)hF¯33sinβdzdαdβdγ,A21=116Nπ202π0π02πη(α,β)hF¯21sinβdzdαdβdγ,A31=116Nπ202π0π02πη(α,β)hF¯31sinβdzdαdβdγ,A32=116Nπ202π0π02πη(α,β)hF¯32sinβdzdαdβdγ,N=18π202π0π02πη(α,β,γ)hsin(β)dzdαdβdγ.
(A3)

We observe that γ does not contribute to the vertical size of the molecule, and redefine η(α, β, γ) as η(α, β). Integrating by γ and z first gives us

A11=Sc16Nπ02π0π(3cos2αcos2β3cosα+1)(hη(α,β))sinβdβdα,A22=Sc16Nπ02π0π(3cos2αcos2β3cos2α3cos2β+2)(hη(α,β))sinβdβdα,A33=Sc16Nπ02π0π(3cos2β1)(hη(α,β))sinβdβdα,A21=Sc16Nπ02π0π3cosαsinαsin2β(hη(α,β))sinβdβdα,A31=Sc16Nπ02π0π3cosαsinβcosβ(hη(α,β))sinβdβdα,A32=Sc16Nπ02π0π3sinαsinβcosα(hη(α,β))sinβdβdα,

where

N=14π02π0π(hη(α,β))sinβdβdα,Sc=13b32.
(A4)

We perform a change of variable, u = cosβ, obtaining

A11=Sc16Nπ02π11(3u2cos2α3cos2α+1)(hη(α,arccosu))dudα,A22=Sc16Nπ02π11(3u2cos2α+3cos2α2)(hη(α,arccosu))dudα,A33=Sc16Nπ02π11(13u2)(hη(α,arccosu))dudα,A21=3Sc16Nπ02π11sinαcosα(1u2)(hη(α,arccosu))dudα,A31=3Sc16Nπ02π11ucosα1u2(hη(α,arccosu))dudα,A32=3Sc16Nπ02π11usinα1u2(hη(α,arccosu))dudα,N=14π02π11hη(α,arccosu)dudα,Sc=13b32.
(A5)

Integrating the terms that do not involve η gives us Equation (9).

B. Derivation of η for an Ellipsoid

In this section we derive the analytical expression for η for an arbitrary ellipsoid, following the notation of Section 2.4. Due to the symmetry of the ellipsoid we consider only one octant in our analysis, expressing all points p on the ellipsoid in that octant as

p(x,y)=(x,y,z(x,y)),
(B1)

where

z(x,y)=c(1x2a2y2b2),
(B2)

for x [set membership] [0, a] and y[0,b(1x2/a2)].

A rotation of the ellipsoid by R(α, β, γ) transforms the coordinates of these points into

x=R11(α,β,γ)x+R12(α,β,γ)y+R13(α,β,γ)z(x,y),y=R21(α,β,γ)x+R22(α,β,γ)y+R23(α,β,γ)z(x,y),z=R31(α,β,γ)x+R32(α,β,γ)y+R33(α,β,γ)z(x,y).
(B3)

We observe that

η(α,β,γ)=z(x,y),
(B4)

where x*(α, β, γ) and y*(α, β, γ) minimize z′.

To find the minimum/maximum value of our rotated ellipsoid, we solve [nabla] z′(x, y) = 0:

z(x,y)x=R31R33cxa2(1x2a2y2b2)=0,
(B5)

z(x,y)y=R32R33cyb2(1x2a2y2b2)=0.
(B6)

Let

x=R31a2R312a2+R322b2+R332c2,
(B7)
y=R32b2R312a2+R322b2+R332c2.
(B8)

It is easy to verify that x* and y* solve Equations (B5) and (B6):

z(x,y)x=R31R33cxa2(1x2a2y2b2)=R31R33cR31a2R312a2+R322b2+R332c2a2R33cR312a2+R322b2+R332c2=0,
(B9)

z(x,y)y=R32R33cyb2(1x2a2y2b2)=R32R33cR32b2R312a2+R322b2+R332c2b2R33cR312a2+R322b2+R332c2=0.
(B10)

Therefore, x*, y* minimizes z′, and the optimal value of z is

z=R33c2R312a2+R322b2+R332c2,
(B11)

Since R312+R322+R332=1. Therefore, from Equations (B4) and (B3), our solution is

η(α,β,γ)=R31(α,β,γ)x+R32(α,β,γ)y+R33(α,β,γ)z=R312a2+R322b2+R332c2,
(B12)

where γ is arbitrary for the Euler rotation defined in Equation (A1).

C. Axis-Angle Representation

The axis-angle representation of a rotation parameterizes the rotation by two values: a unit vector indicating the orientation of the axis, u, and an angle, θ, describing the magnitude of the rotation about that axis. The direction of rotation around the axis u is determined by the right-hand rule. One advantage of the axis-angle representation over Euler angles representation is that one can easily quantify the magnitude of the rotation by the size of the rotation angle θ.

Given a rotation matrix R, θ is computed as

θ=arccos(12(R11+R22+R331)).
(C1)

Footnotes

Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

References

1. Bax A. Weak alignment offers new NMR opportunities to study protein structure and dynamics. Protein Science. 2003;12(1):1–16. [PubMed]
2. Tjandra N, Bax A. Direct measurement of distances and angles in biomolecules by NMR in a dilute liquid crystalline medium. Science. 1997;278(5340):1111–1114. [PubMed]
3. Tycko R, Blanco FJ, Ishii Y. Alignment of biopolymers in strained gels: A new way to create detectable dipole-dipole couplings in high-resolution biomolecular NMR. Journal of the American Chemical Society. 2000;122(38):9340–9341.
4. Sass H, Musco G, Stahl S, Wingfield P, Grzesiek S. Solution NMR of proteins within polyacrylamide gels: Diffusional properties and residual alignment by mechanical stress or embedding of oriented purple membranes. Journal of Biomolecular NMR. 2000;18(4):303–309. [PubMed]
5. Ruckert M, Otting G. Alignment of biological macromolecules in novel nonionic liquid crystalline media for NMR experiments. Journal of the American Chemical Society. 2000;122(32):7793–7797.
6. Bewley C, Clore G. Determination of the relative orientation of the two halves of the domain-swapped dimer of Cyanovirin-N in solution using dipolar couplings and rigid body minimization. Journal of the American Chemical Society. 2000;122(25):6009–6016.
7. Warren J, Moore P. Application of dipolar coupling data to the refinement of the solution structure of the Sarcin-Ricin loop RNA. Journal of Biomolecular NMR. 2001;20(4):311–323. [PubMed]
8. Bewley C. Rapid validation of the overall structure of an internal domain-swapped mutant of the anti-HIV protein Cyanovirin-N using residual dipolar couplings. Journal of the American Chemical Society. 2001;123(5):1014–1015. [PubMed]
9. Azurmendi H, Martin-Pastor M, Bush C. Conformational studies of Lewis X and Lewis A trisaccharides using NMR residual dipolar couplings. Biopolymers. 2002;63(2):89–98. [PubMed]
10. Bernadó P, Blanchard L, Timmins P, Marion D, Ruigrok R, Black-ledge M. A structural model for unfolded proteins from residual dipolar couplings and small-angle X-ray scattering. Proceedings of the National Academy of Sciences. 2005;102(47):17002–17007. [PubMed]
11. Mukrasch M, Markwick P, Biernat J, von Bergen M, Bernado P, Griesinger C, Mandelkow E, Zweckstetter M, Blackledge M. Highly populated turn conformations in natively unfolded tau protein identified from residual dipolar couplings and molecular simulation. Journal of the American Chemical Society. 2007;129(16):5235–5243. [PubMed]
12. Ryabov Y, Fushman D. Structural assembly of multidomain proteins and protein complexes guided by the overall rotational diffusion tensor. Journal of the American Chemical Society. 129(25):7894–7902. [PMC free article] [PubMed]
13. Zweckstetter M, Bax A. Prediction of sterically induced alignment in a dilute liquid crystalline phase: Aid to protein structure determination by NMR. Journal of the American Chemical Society. 2000;122(15):3791–3792.
14. Zweckstetter M. NMR: Prediction of molecular alignment from structure using the PALES software. Nature Protocols. 2008;3(4):679–690. [PubMed]
15. Eisenhaber F, Lijnzaad P, Argos P, Sander C, Scharf M. The double cubic lattice method: Effcient approaches to numerical integration of surface area and volume and to dot surface contouring of molecular assemblies. Journal of Computational Chemistry. 1995;16(3):273–284.
16. Fernandes MX, Bernado P, Pons M, Garcia de la Torre J. An analytical solution to the problem of the orientation of rigid particles by planar obstacles. application to membrane systems and to the calculation of dipolar couplings in protein NMR spectroscopy. Journal of the American Chemical Society. 2001;123(48):12037–12047. [PubMed]
17. Almond A, Axelsen J. Physical interpretation of residual dipolar couplings in neutral aligned media. Journal of the American Chemical Society. 2002;124(34):9986–9987. [PubMed]
18. Azurmendi H, Bush C. Conformational studies of blood group A and blood group B oligosaccharides using NMR residual dipolar couplings. Carbohydrate Research. 2002;337(10):905–915. [PubMed]
19. Van Loan CF. Introduction to Scientific Computing: A Matrix-Vector Approach Using MATLAB. Prentice-Hall, Inc; Upper Saddle River, NJ, USA: 1997.
20. Ryabov Y, Geraghty C, Varshney A, Fushman D. An effcient computational method for predicting rotational diffusion tensors of globular proteins using an ellipsoid representation. Journal of the American Chemical Society. 2006;128(48):15432–15444. [PubMed]
21. Saupe A, Englert G. High-resolution nuclear magnetic resonance spectra of orientated molecules. Physical Review Letters. 1963;11(10):462–464.
22. Murnaghan F. The element of volume of the rotation group. Proceedings of the National Academy of Sciences. 1950;36(11):670–672. [PubMed]
23. Lindstrom P, Turk G. Fast and memory effcient polygonal simplification. VIS ’98: Proceedings of the Conference on Visualization ’98; Los Alamitos, CA, USA: IEEE Computer Society Press; 1998. pp. 279–286.
24. Lindstrom P, Turk G. Evaluation of memoryless simplification. IEEE Transactions on Visualization and Computer Graphics. 1999;5(2):98–115.
25. Varshney A, Brooks F, Jr, Wright W. Linearly scalable computation of smooth molecular surfaces. IEEE Computer Graphics and Applications. 1994;14(5):19–25.
26. Varshney A, Brooks F., Jr Fast analytical computation of Richard’s smooth molecular surface. IEEE Visualization’93 Proceedings; 1993. pp. 300–307.
27. Welzl E. New Results and New Trends in Computer Science. Vol. 555. Springer-Verlag; 1991. Smallest enclosing disks (balls and ellipsoids) pp. 359–370.
28. Todd MJ, Yildirim EA. On Khachiyan’s algorithm for the computation of minimum volume enclosing ellipsoids. Discrete Applied Mathematics. 2007;155(13):1731–1744.
29. Cornilescu G, Marquardt JL, Ottiger M, Bax A. Validation of protein structure from anisotropic carbonyl chemical shifts in a dilute liquid crystalline phase. Journal of the American Chemical Society. 1998;120(27):6836–6837.
30. Losonczi JA, Andrec M, Fischer MWF, Prestegard JH. Order matrix analysis of residual dipolar couplings using singular value decomposition. Journal of Magnetic Resonance. 1999;138(2):334–342. [PubMed]
31. Cai M, Huang Y, Zheng R, Wei S, Ghirlando R, Lee M, Craigie R, Gronenborn A, Clore G. Solution structure of the cellular factor BAF responsible for protecting retroviral DNA from autointegration. Nature Structural Biology. 1998;5(10):903–909. [PubMed]
32. Kuszewski J, Gronenborn AM, Clore GM. Improving the packing and accuracy of NMR structures with a pseudopotential for the radius of gyration. Journal of the American Chemical Society. 1999;121(10):2337–2338.
33. Ulmer TS, Ramirez BE, Delaglio F, Bax A. Evaluation of backbone proton positions and dynamics in a small protein by liquid crystal NMR spectroscopy. Journal of the American Chemical Society. 2003;125(30):9179–9191. [PubMed]
34. Drohat A, Tjandra N, Baldisseri D, Weber D. The use of dipolar couplings for determining the solution structure of rat apo-S100B (ββ) Protein Science. 1999;8(4):800–809. [PubMed]
35. Bewley C, Gustafson K, Boyd M, Covell D, Bax A, Clore G, Gronenborn A. Solution structure of Cyanovirin-N, a potent HIV-inactivating protein. Nature Structural Biology. 1998;5(7):571–578. [PubMed]
36. de Alba E, De Vries L, Farquhar M, Tjandra N. Solution structure of human GAIP (Gα interacting protein): A regulator of G protein signaling. Journal of Molecular Biology. 1999;291(4):927–939. [PubMed]
37. Schwalbe H, Grimshaw S, Spencer A, Buck M, Boyd J, Dobson C, Redfield C, Smith L. A refined solution structure of hen lysozyme determined using residual dipolar coupling data. Protein Science. 2001;10(4):677–688. [PubMed]
38. Jain NU, Tjioe E, Savidor A, Boulie J. Redox-dependent structural differences in putidaredoxin derived from homologous structure refinement via residual dipolar couplings. Biochemistry. 2005;44(25):9067–9078. [PubMed]
39. Bax A, Kontaxis G, Tjandra N. Dipolar couplings in macromolecular structure determination. In: Thomas VD, James L, Schmitz U, editors. Part B: Nuclear Magnetic Resonance of Biological Macromolecules, Vol. 339 of Methods in Enzymology. Academic Press; 2001. pp. 127–174.