Computing the Rotational Diffusion of ... - ACS Publications

Article pubs.acs.org/JPCB

Computing the Rotational Diffusion of Biomolecules via Molecular Dynamics Simulation and Quaternion Orientations Po-chia Chen,*,† Maggy Hologne, and Olivier Walker* Université de Lyon, CNRS, Université Claude Bernard Lyon1, Ens de Lyon, Institut des Sciences Analytiques, UMR 5280, 5 rue de la Doua, F-69100 Villeurbanne, France S Supporting Information *

ABSTRACT: Rotational diffusion (Drot) is a fundamental property of biomolecules that contains information about molecular dimensions and solute−solvent interactions. While ab initio Drot prediction can be achieved by explicit all-atom molecular dynamics simulations, this is hindered by both computational expense and limitations in water models. We propose coarse-grained force fields as a complementary solution, and show that the MARTINI force field with elastic networks is sufficient to compute Drot in >10 proteins spanning 5−157 kDa. We also adopt a quaternion-based approach that computes Drot orientation directly from autocorrelations of best-fit rotations as used in, e.g., RMSD algorithms. Over 2 μs trajectories, isotropic MARTINI+EN tumbling replicates experimental values to within 10−20%, with convergence analyses suggesting a minimum sampling of >50 × τtheor to achieve sufficient precision. Transient fluctuations in anisotropic tumbling cause decreased precision in predictions of axisymmetric anisotropy and rhombicity, the latter of which cannot be precisely evaluated within 2000 × τtheor for GB3. Thus, we encourage reporting of axial decompositions Dx, Dy, Dz to ease comparability between experiment and simulation. Where protein disorder is absent, we observe close replication of MARTINI+EN Drot orientations versus CHARMM22*/TIP3p and experimental data. This work anticipates the ab initio prediction of NMR-relaxation by combining coarse-grained global motions with all-atom local motions.

■

fields.11,24−31 This has the advantage of explicitly incorporating conformational fluctuations and solvent dynamics that must otherwise be modeled in static approaches, and permits verification of the fidelity of model-free assumptions standard in NMR analysis.22 As a bonus, the NMR-relaxation parameters underlying experimental Drot methods can be directly calculated from trajectories.11,32−36 These relaxation parameters describe the dynamics and time scales of local disorder critical to understanding molecular motion. However, Drot from MD simulations is known to be overestimated due to the low viscosity of water models (TIP3p and SPC/E) associated with popular force fields.27,37 While SPC/E artifacts are a marginal 10−20%, TIP3p deviations are 2.5−3-fold, inspiring follow-up mitigation studies that modify either the water model itself,29 or the protein tumbling velocities.28 Meanwhile, simulated Drot in TIP3p simulations has been ignored for the purposes of NMR parameter prediction,35,36 in favor of a combination of simulated internal motions and externally supplied D rot information.

INTRODUCTION Rotational diffusion is a fundamental property describing the tumbling of biomolecules in solution, which affects processes such as enzyme kinetics, partner cognition, and macromolecular assembly.1−3 Its association with molecular shape also extends useful applications toward studies of molecular crowding,4,5 as well as structural/dynamics studies involving complexation6−8 and significant disorder,9−12 whose properties may hinder popular crystallization approaches. Thus, methods to determine molecular diffusive properties remain a staple in contemporary in vivo and structural investigations. Experimental and computational methods describe molecular tumbling via the rotational diffusion tensor Drot, or its corresponding tumbling time(s) τiso. For construction of a cohesive model, relaxation times are generally measured in NMR-heteronuclear relaxation13−15 and polarized fluorescence16−18 experiments, and then compared to predictions based on atomistic coordinates. The availability of such atomistic models allows both Drot and its orientation to be evaluated with respect to the coordinate axis.10,15,19−21 Where one is lacking, model-free formalisms can also be used,22,23 at the cost of losing orientation information. Explicit calculation of Drot can also be carried out via molecular dynamics simulation (MD) using atomistic force © XXXX American Chemical Society

Received: November 21, 2016 Revised: January 17, 2017 Published: February 3, 2017 A

DOI: 10.1021/acs.jpcb.6b11703 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 1. Quaternion-based definition of molecular orientation, illustrated using ubiquitin. An animated version is available as Supporting Information. Ubiquitin is shown a cartoon representation in four conformations A, B, C, and D, along with its moment of inertia axes X′, Y′, and Z′ (MoI). The lab/reference frame is aligned with the MoI frame of A. To the right of each conformer is a mapping of their respective orientation Q to the unit sphere, by plotting Qx, Qy, and Qz. The numerical values of Q are labeled above or below each figure. Thus, relative to A: Q(A) is the unit quaternion (1, 0, 0, 0), Q(B) is a 30° rotation around the negative-Y axis, Q(C) is a ∼ 60° rotation around a vector pointing out of the page, and Q(D) is an 120° rotation around (1/√3, 1/√3, 1/√3), resulting in a cyclic permutation of ubiqutin’s MoI.

expressed as a quaternion q with four components, defined with respect to a fixed coordinate frame:

Explicit calculations of Drot also require convergence of respective rotational autocorrelation functions, which generally require time scales of at least 1−2 orders of magnitude above the expected tumbling time.27,33,34 Given that proteins have τiso ranging from 3 ns up to μs for larger complexes, this requirement can render all-atom-based approaches infeasible. Thus, coarse-grained (CG) alternatives may be utilized to reach desired time scales.31 We propose the use of a general-purpose coarse-grained force field38−40 to derive Drot through explicit dynamics, which can tackle both convergence criteria and water model limitations (which are comparatively easy to adjust for CG force fields). The MARTINI protein force field41 with elastic networks is used as a proof-of-principle test, using anisotropic Drot predictions of 10 proteins between 5 and 157 kDa as a means of evaluation against literature data. The rotational diffusion tensor is calculated by describing molecular orientation as a quaternion Q, analogous to Delong et al.42 and Chevrot et al.,43 but using quaternion-rotations to directly derive diffusion properties. This is in contrast with Euler-angle-based descriptions from which diffusion cannot be as trivially obtained. Applicability of the quaternion method is also validated for all-atom simulations of a number of proteins, with a future expectation that ab initio prediction of NMRrelaxations can be carried out by MD simulations alone.

q = (qw , qx , qy , qz)

(1)

In NMR and atomistic simulation literature, this fixed coordinate frame is usually termed the “lab frame” and “reference frame”, respectively. We will adopt the latter convention in this paper. q has the property qw2 + qx2 + qy2 + qz2 = 1, which reduces its four parameters to the 3 degrees of freedom associated with rotation. We note that q can also be expressed in the more physically transparent axis-angle notation: q = (cos(θq/2), sin(θq/2)vx , sin(θq/2)vy , sin(θq/2)vz) = (cos(θq/2), sin(θq/2)vq) (2)

In this form, one interprets q to be a rotation of angle θq around a unit vector vq. Thus, q is the single-rotation equivalent to the popular three Euler-angle rotations α, β, and γ around coordinate axes (see Figure 1 for a visual aid). Since q and −q denote the same rotation in space, we restrict q to the region qw ≥ 0 without loss of generality. This theorical description of rotation is realized in atomistic simulations by associating atomic coordinates of the molecule with respective internal coordinate frames. The reference frame is first fixed by defining, e.g., a published PDB file, as the zerorotation. For an arbitrary position A of the same molecule, its relative translation and rotation to the reference can be found by minimizing the atomic root-mean-squared deviation (RMSD), using alignment procedures in established algorithms.47,49 As the optimal rotation qrot in these algorithms is generally well-defined and unique, one can adopt the optimal

■

THEORY We present here concisely the definitions of quaternion orientation and rotation and the derivation of the molecular rotational diffusion tensor from a simulation trajectory of the biomolecule. As the language is based on a structural biologist’s perspective, we encourage further reading for those interested in a comprehensive treatment of quaternions.44−48 Quaternions as Rotations and Orientations. As a molecule tumbles over time, its rotation in space can be B


Article


Figure 2. Autocorrelation functions of the quaternion components and PAF behavior, using ubiquitin as an illustration. The Δt period considered for exponential fit was cut off either at 128 ns in the main plots, or at 8 ns in the inset. (A) Isotropic rotational tumbling via θq (gray), along with the single-exponential fit (dotted black). (B) Corresponding anisotropic rotational tumbling via the components of vq, qx, qy, and qz (red, green, and blue respectively), af ter rotation into the PAF at 80 ps. Each single-exponential fit is plotted in dotted black. (C) Rotation qPB required to diagonalize ⟨qiqj⟩ at Δt (see eq 14). The components of qPB are shown in black (w), red (x), green (y), and blue (z). (D) MARTINI-model of ubiquitin in the MoI frame, with backbone colored according to fixed secondary structure elements and side-chain atoms colored according to residue type (green, polar; red, acidic; blue, basic; white, other). (E) Visualization of the MoI frame (gray axes) and the drifting PAF axes as a function of Δt. Lengths of axes are determined qualitatively by Dx (red), Dy (green), and Dz (blue). An experimental Dz is shown in transparent yellow for comparison.50

rotations between time t and t + Δt, sampled along the entire trajectory. The set of rotations from Q(t) to Q(t + Δt) is collected in the internal frame (BODY) of the molecule by a frame transform of eq 5 from REF to BODY:

coordinate-rotation from the reference to A as the definition of the orientation-quaternion QA: Q A ≔ qrot ⊗ Q REF

(3)

Q REF ≔ (1, 0, 0, 0)

(4)

BODY

Here, the symbol ⊗ denotes quaternion multiplication, and the use of upper case distinguishes orientations Q from rotations q between two orientations. Rotations from A to B in the reference frame (REF) can then be expressed as two rotations from A to REF then from REF to B REF

qBA = Q B ⊗

Q A−1

=Q (t )−1 ⊗ Q (t + Δt )

(7) (8)

q(t, Δt) contains information about the rotation diffusion tensor Drot, and its observed autocorrelation decays are plotted in Figure 2 using MARTINI ubiquitin as an example. The derivation of isotropic and anisotropic components of Drot are further described below. Isotropic Tumbling. In the isotropic approximation, rotational diffusion is a function of the total rotation angle θq alone, without regard to coordinate frame or rotation direction. According to rigid-body Brownian dynamics, the probability distribution of finding rotation θq after Δt after starting at zerorotation is (eq 51 from Furry44):

(5)

where the inverse rotation is q−1 = (qw , −qx , −qy , −qz)

q(t , Δt ) = Q (t )−1 ⊗REF q(t , Δt ) ⊗ Q (t )

(6)

We highlight the analogy of eq 5 to vector operations v = rB − rA, which derive distance v between two positions. The coordinate frame in which this rotation occurs is clarified where necessary by prepending a superscript to q. Extracting Rotational Diffusion Tensors from Quaternions. We begin by obtaining the distribution of molecular C


Article

The Journal of Physical Chemistry B ∞

As a corollary of eq 15, the trace of T in any coordinate frame is sufficient to extract isotropic rotational diffusion. To align BODYq with Drot, one requires a further frame rotation from BODY to PAF (termed qPB). Rather than following Favro45 to derive qPB, we instead diagonalize BODYT using eigenvalue decomposition. The resulting matrix of eigenvectors correspond to the axis vectors PAFei within BODY, from which qPB can be obtained. PAF T as constructed from PAFq now shares the same frame as Drot, enabling direct identification between the tensor elements. In this f rame, we argue by reason for symmetry that eq 11 and eq 15 together imply:

P(θq) = π −1 ∑ e−n(n + 1)DisoΔt (2n + 1)(cos nθq n=0

− cos(n + 1)θq)

(9)

Prand(θ ) = 2/π sin 2(θq/2)

(10)

Here, Diso is the isotropic diffusion constant, and Prand(θ) is the long-time limit of P(θq) where no correlation remains and the distribution becomes random.44,51 The ensemble average ⟨cos θq(Δt)⟩ removes all terms above n = 1, resulting in a single-exponential decay:

∫0

π

⟨cos θq(Δt )⟩ =

⟨qi 2⟩(Δt ) = 0.25(1 − e−2DiΔt )

PAF

dθq cos(θq)P(θq)

= −0.5 + 1.5e−2DisoΔt

⟨1 − 2Tii2(Δt )⟩ = 0.5e−2DiΔt + 0.5

PAF

(11)

We note that the above decomposition is analogous to the decomposition of eq 12 into eq 13 for unit vectors, and that eq 10 can also be utilized to demonstrate the correct convergence of ⟨qi2⟩(Δt) to 0.25. For the purposes of visual familiarity between isotropic and anisotropic decay, we will adopt eq 17 in this work. For flexible biomolecules in practice, it is nontrivial to define an appropriate PAF since qPB depends on Δt (Figure 2CE, further discussed below). This dependence is due to influences upon Drot by nonlinear diffusive behavior caused by conformational fluctuations. Finite sampling also contributes significantly to errors at larger time scales. For the sake of simplicity and consistency, we avoid here attempts to match the PAF as found by NMR approaches and choose instead to best preserve rigidbody assumptions by selecting qPB at small Δt (80 ps) and cut off the fitting of eq 17 at 8 ns.

Cl(Δt ) = ⟨Pl(n̂ (t ) ·n̂ (t + Δt ))⟩ (12)

C1(Δt ) = n̂ ·x̂e−(Dy + Dz)Δt + n̂ ·ŷ e−(Dx + Dz)Δt + n̂ ·zê −(Dx + Dy)Δt

(13)

Here, Pl(x) is the Legendre polynomial. It is worth noting that the l = 2 solution is used when applying the unit-vector method above to NMR and fluorescence, so as to describe the underlying physical process (i.e., τiso,NMR = (6Diso)−1). Since this work directly retrieves dynamics information from simulation trajectories, we adopt the simpler l = 1 solution where τiso = (2Diso)−1, and anisotropic diffusion is described by only three time-constants instead of five. Anisotropic Tumbling. The full diffusion tensor Drot is conventionally defined as a 3 × 3 tensor in the principal-axis frame (PAF), where its diagonal elements Dx ≤ Dy ≤ Dz describe the rate of tumbling along the respective axes, and offdiagonal elements are zero.45,53 The elements of Drot are also commonly reported in terms of symmetric expansions Diso = (Dx + Dy + Dz)/3, Δ = 2Dz/(Dx + Dy), and ρ = 3(Dy − Dx)/ (2Dz − (Dx + Dy)).21 The rotation of PAF relative to REF is also required to determine the correct orientation of Drot elements with respect to an atomic structure. The axis-angle form of q implies that the isotropic decay ⟨cos θq(Δt)⟩ can be decomposed into its anisotropic components. To do this, we construct the 3 × 3 tensor T from the vector elements of q Tij(Δt ) = ⟨qiqj⟩(Δt )

i, j = x, y, z

■

RESULTS AND DISCUSSIONS Isotropic Tumbling of MARTINI Proteins Are Comparable with Experiment. The predicted versus measured isotropic tumbling Diso for a large range of proteins is shown in Figure 3, standardized to 293 K in 0% D2O (see Methods section for details). We find variation of ∼10% between 2 μs simulations and NMR/fluorescence measurements for the majority of proteins studied, with exceptions either possessing significant disorder (UIM, VHS, and plexin) or being large, fluorescence-based measurements (ovalbumin and aldolase). These we discuss further below. We rationalize the satisfactory MARTINI performance by noting that its parametrization strategy matches fortuitously with the variables affecting Drot. The force field specifically targets partition and hydration free energies between hydrophobic membranes and aqueous solvent,39 which should also reproduce an appropriate balance of protein−solvent interactions across different amino acids. In addition, the shear viscosity of MARTINI water beads at 298 K is also close to experiment,58 indicating fidelity of solvent−solvent interactions. These together suggest that a compatibility of rotational diffusion should not be unexpected, although our results are in contrast to the ∼4-fold faster translational diffusion observed for MARTINI.41,59 The independence of translational versus rotational diffusion properties (also shown for SPC/E water24) implies that the two phenomena should not be conflated. Validation of Oligomeric States. To disambiguate the influence of alternate oligomerization states upon Drot measure-

(14)

noting that the diagonal elements of T correspond to the components of ⟨cos θq(Δt)⟩, by virtue of the identity in eq 2: ⟨cos θq⟩ = 1 − 2⟨sin 2(θq/2)⟩ x ,y,z

=1−2

∑ ⟨qi2⟩ = 1 − 2tr(T ) i

i = x, y, z (17)

Relation to the Unit-Vector Method. Equation 11 is analogous to the rotation diffusion of a unit vector n̂(t) in a moving body as used in NMR spin relaxation and fluoresence.27,52 We transcribe here eqs 4−6 showing the relationship between the autocorrelation of n̂ with Drot:27 = e−l(l + 1)DisoΔt

(16)

(15) D


Article


Figure 3. Isotropic tumbling Diso via quaternion components from 2 μs MARTINI-simulations, compared with measured values via NMR (circles) or fluorescence correlation (triangle) in literature. All values have been standardized to 293 K and 0% D2O, and dotted lines show area of agreement within 20%. Error bars show standard deviation of Diso calculated from 500 ns subchunks. Proteins examined are labeled on the plot and presented here in order of increasing molecular weight: Stam2 UIM-domain,8 GB3 domain,14 Ubq,54 Titin Ig27domain,9 plexin-B1 cytoplasmic domain,7 Stam2 VHS-domain plus its complex with Ubq,54 D-peroxiredoxin dimer,55 green fluorescent protein,56 ovalbumin,18 and aldolase.57

Figure 4. Top: Estimation of sampling sufficiency by comparing distributions of θQ at tsim = 2 μs vs Prand(θ) (solid black, see text), illustrated using GB3 (dotted indigo-blue), GFP (dotted light green), and aldolase (dashed−dotted red). Bottom: Unweighted χ-difference between θQ and Prand(θ) as a function of total simulation time tsim scaled by theoretical tumbling time τtheor (see Table 2) in the Methods section.

ments, we computed both monomeric and dimeric forms of GFP,60 ovalbumin,61 and D-Prx,55 as well as dimer/tetrameic aldolase. In the first three cases, tumbling values are consistent with the presence of a single dominant oligomeric state, rather than a mixture: D-Prx is dimeric, while GFP and ovalbumin exist as monomers. These observations are consistent with the original interpretations given in respective studies. In the case of aldolase, experimental tumbling suggests significant dissociation of the tetrameric into a dimeric variant, which was not considered by Pieper and Enderlein.57 This discrepancy warrants future exploration to exclude possible size-dependent artifacts in the MARTINI+EN simulations. Aldolase notwithstanding, the overall results support the potential usefulness of Drot in determining oligomerization states. Sufficiency of Finite-Length Trajectories for Diso Calculations. Given the relative ease of reaching longer time scales in CG-MD, we investigated possible refinements to existing advice on Drot convergence requirements, given as 20− 100 × τiso via all-atom MD,27 or >100 × τiso via random-walk modeling.62 To minimize the dependence of our evaluations upon the characteristic of individual proteins, we use the sampling of the unit sphere by molecular orientation θQ as a proxy for convergence, and scale all protein time scales according to its theoretical tumbling time τtheor calculated on the basis of molecular size alone (see Table 2 in the Methods section). Thus, we plot below the distribution of θQ for GB3, monomeric GFP, and tetrameric aldolase versus the expected random distribution Prand(θ), and the χ-squared difference as a function of sampling length tsim/τtheor (Figure 4).

The three examples in Figure 4A exhibit different levels of sampling extent due to the use of a fixed 2 μs trajectory length. GB3 is an overkill example where ∼400 × τtheor was simulated, whose corresponding θQ distribution clearly approaches Prand(θ). GFP was simulated for 96 × τtheor, and its resemblance to Prand(θ) is qualitative, whereas 2 μs corresponds to a mere 15 τtheor for aldolase, and is clearly unconverged. Curiously, extensions to 5 μs for these largest proteins reduce Diso uncertainty without appreciably altering the actual value. This suggests that Diso accuracy is significantly easier to achieve than its precision. The χ-deviation for all tested systems and Diso predictions (Figure 3) together suggest a criterion of χ < 0.1 and tsim > 50τtheor as a bare minimum for reliable calculation of isotropic tumbling. We further note that the smaller criterion here versus that of Lu and Bout62 is likely due to our adoption of an 8 ns cutoff with respect to fitting quaternion autocorrelations (eq 11). Anisotropic Rotational Diffusion. We now report anisotropic rotational diffusion parameters of studied proteins in Table 1 and Table S2, including both prolate and oblate axisymmetric expansions. Visual comparisons with available literature data are also provided in Figures S1 and S2. We find that both accuracy and precision of Δ are somewhat poorer compared to Diso; a number of Δ predictions deviate by >20% (Titin, VHS-Ubq, and both D-Prx forms), but the increased E


F

0.60 0.44 0.28 0.11 0.49 2.9 1.3 0.19 0.81 1.7 2.2 0.16 0.38 1.3 0.05 0.07 0.05 0.067 0.150 0.049 0.38 0.10 0.043 0.083 0.071 0.048 0.004 0.012 0.018 0.011 0.032

14.46 12.30 5.85 5.54 5.32 15.4 17.0 3.62 3.73 9.8 10.3 3.29 3.81 11.7 2.07 1.50 1.10 1.16 1.17 0.509 1.87 1.77 1.05 1.05 0.751 0.675 0.269 0.259 0.126 0.122 0.468

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

Dz 1.85 2.25 1.01 0.42 0.53 5.6 5.7 0.32 0.52 3.3 4.1 0.53 0.66 4.5 0.21 0.18 0.08 0.17 0.15 0.025 0.27 0.19 0.13 0.02 0.075 0.026 0.022 0.012 0.015 0.005 0.064

7.26 6.58 4.54 4.34 4.23 12.7 10.4 3.22 3.29 8.3 9.0 2.95 3.19 7.5 1.90 1.39 1.01 0.917 0.908 0.425 1.76 1.59 0.673 0.680 0.528 0.513 0.232 0.229 0.109 0.109 0.313

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± 0.53 0.78 0.34 0.19 0.17 1.6 1.8 0.12 0.12 1.5 2.1 0.15 0.24 1.2 0.09 0.11 0.11 0.091 0.078 0.023 0.09 0.06 0.052 0.035 0.021 0.020 0.014 0.002 0.007 0.003 0.036

Diso

ρ 0.01 ± 0.08 0.12 ± 0.12 0.25 ± 0.72 0.27 ± 0.27 0.31 ± 1.04 0.72 ± 3.96 0.20 ± 0.38 1.26 ± 3.81 1.28 ± 2.24 0.92 ± 4.24 0.21 ± 1.76 0.73 ± 0.37 0.17 ± 3.08 0.12 ± 0.41 0.40 ± 46.27 0.07 ± 63.11 1.32 ± 3.04 0.17 ± 0.55 0.47 ± 0.98 0.10 ± 1.51 0.69 ± 3.34 0.14 ± 1.34 0.02 ± 0.17 0.10 ± 0.19 0.09 ± 0.43 0.33 ± 0.31 0.37 ± 0.76 0.81 ± 1.33 1.13 ± 2.87 1.66 ± 1.74 0.11 ± 0.11

prolate 3.95 ± 0.78 3.31 ± 0.67 1.50 ± 0.26 1.48 ± 0.09 1.44 ± 0.18 1.36 ± 0.60 2.39 ± 0.88 1.20 ± 0.17 1.21 ± 0.30 1.29 ± 0.37 1.24 ± 0.38 1.18 ± 0.23 1.32 ± 0.27 2.20 ± 1.12 1.14 ± 0.10 1.13 ± 0.16 1.14 ± 0.09 1.46 ± 0.17 1.51 ± 0.25 1.33 ± 0.15 1.10 ± 0.23 1.18 ± 0.20 2.17 ± 0.38 2.13 ± 0.27 1.80 ± 0.34 1.56 ± 0.13 1.26 ± 0.08 1.21 ± 0.09 1.26 ± 0.14 1.18 ± 0.04 1.98 ± 0.19

Δ

0.87 ± 3.08 0.50 ± 1.35 2.60 ± 0.38

0.76 ± 0.14 0.79 ± 0.08 0.64 ± 0.06

ρ 2.97 ± 0.31 2.56 ± 0.39 2.19 ± 1.93 2.15 ± 0.69 2.06 ± 2.05 1.32 ± 1.41 2.33 ± 0.83 0.77 ± 0.91 0.76 ± 4.68 1.08 ± 2.55 2.31 ± 2.98 1.32 ± 0.71 2.41 ± 1.56 2.56 ± 2.73 1.86 ± 2.41 2.75 ± 4.40 0.72 ± 1.47 2.43 ± 1.51 1.73 ± 1.83 2.65 ± 2.00 1.36 ± 3.61 2.50 ± 60.85 2.91 ± 0.63 2.65 ± 0.64 2.68 ± 1.64 2.01 ± 0.71 1.93 ± 4.01

oblate 0.40 ± 0.03 0.41 ± 0.06 0.75 ± 0.07 0.76 ± 0.02 0.77 ± 0.08 0.74 ± 0.30 0.52 ± 0.18 0.80 ± 0.09 0.79 ± 0.14 0.76 ± 0.04 0.87 ± 0.21 0.86 ± 0.15 0.84 ± 0.22 0.59 ± 0.24 0.91 ± 0.05 0.94 ± 0.14 0.85 ± 0.19 0.78 ± 0.04 0.71 ± 0.04 0.85 ± 0.20 0.92 ± 0.21 0.91 ± 0.16 0.62 ± 0.07 0.61 ± 0.05 0.69 ± 0.09 0.72 ± 0.05 0.85 ± 0.07

Δ

0.233 ± 0.018

0.732 ± 0.084

0.751 ± 0.035

1.59 ± 0.15

1.26 ± 0.12 1.13 0.65 N/A

1.54

2.94

0.08 0.08 0.11 0.04

1.23

3.52

± ± ± ±

1.34 ± 0.01

4.61 ± 0.00

1.73 1.57 1.20 0.83

3.17

Δ

expt 8.01 ± 2.77

Diso

ρ

0.393 ± 0.059

N/A 0.71 0.63 N/A

N/A

0.44

0.26 ± 0.06

N/A

a Coarse-grained MD entries were simulated for 2 μs (marked A, B) with extensions explicitly noted, while all-atom MD entries using TIP3p were simulated for 500 ns. bD-values are standardized to 293 K and 0% D2O, presented in units of 107 s−1, and sorted by Dx ≤ Dy ≤ Dz. Anisotropy Δ and rhombicity ρ values are reported for both prolate and oblate axial-symmetric expansions, with ρ < 1 entries in bold. Uncertainties indicate standard deviation of calculations conducted over four equal sub-segments, using the same PAF-orientation qPB as for the full trajectory. Experimental sources as in Figure 3.

3.69 4.06 4.05 3.90 3.86 10.3 6.5 3.27 3.36 8.3 8.5 2.90 2.93 5.1 1.85 1.34 1.02 0.815 0.837 0.387 1.74 1.51 0.488 0.512 0.426 0.459 0.221 0.227 0.110 0.113 0.228

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

0.24 0.45 0.35 0.20 0.16 1.6 1.3 0.30 0.33 1.3 1.1 0.24 0.41 0.76 0.10 0.21 0.232 0.099 0.065 0.062 0.15 0.19 0.063 0.043 0.026 0.017 0.023 0.012 0.012 0.005 0.020

± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± ±

UIMA UIMB GB3A (10 μs) GB3B GB3b (+5 a.a. tail)b UbqA UbqB UbqAb UbqBb TitinA TitinB Titinb plexinA VHSA VHS+UbqA GFPA GFPB dim. GFPA mono D-PrxA mono D-PrxB dim D-PrxA dim D-PrxB mono ovaA mono ovaB dim ovaA (5 μs) tet aldA (5 μs) (dimer 5 μs)

3.63 3.36 3.72 3.57 3.52 12.3 7.8 2.76 2.80 6.9 8.2 2.65 2.82 5.60 1.78 1.33 0.904 0.775 0.714 0.379 1.66 1.49 0.480 0.476 0.406 0.406 0.207 0.203 0.090 0.093 0.244

Dy

Dx

protein

sim

Table 1. Anisotropic Tumbling of Proteins by Diagonalization of ⟨vq⟩ in a Chosen PAFa

The Journal of Physical Chemistry B Article


Article


secondary structure elements around the hydrophobic core. The 3-residue unstructured tail is also folded against the rest of Ubq in MARTINI. These effects perturb the moment of inertia and presumably also the resultant solution tumbling. Similar shifts also appear in GB3 (Figure 6), but its primitive structure

uncertainties imply that results remain consistent and are therefore qualitatively correct. This is contrasted with ρ predictions, which cannot be reliably computed according to our uncertainty-estimation method. Nevertheless, ρ qualitatively identifies the appropriate axisymmetric expansion except for ubiquitin (further discussed below). As an anecdotal exercise, we extended one GB3 simulation to 10 μs, or >2000 × τtheor in order to measure convergence improvements of all Drot components as a function of sampling extent (Figures S3 and S4). Despite both Δ and ρ reaching “agreement to within error bars” at 50 × τtheor, by which time qPB also converges, this GB3 trajectory required >500 × τtheor to breach 10% uncertainty in Δ, while uncertainty in ρ never reduces below 100%. Inspection of the precision in axial components of Drot reveals the limitations of Δ and ρ shown above. For nearspherical proteins, two Di components are often numerically close enough such that their uncertainties overlap. In such cases, the protein may transiently violate the axis order Dx < Dy < Dz calculated on the basis of the entire trajectory. This forces ρ to be poorly defined for such segments, and drastically alters observed anisotropies. Thus, convergence requirements of ρ (and to a lesser extent Δ) are challenging to fulfill via direct simulations, and we recommend the direct reporting of Dx, Dy, and Dz to best enable future comparisons, as is done by, e.g., Berlin et al.14 Accuracy Limitations of MARTINI+EN Models. The above observations on deteriorating consistency of anisotropic MARTINI+EN predictions might also be explained by systematic artifacts introduced via the coarse-grained force field. We further explore this avenue below by comparing the full Drot tensor including PAF orientations between MARTINI +EN, CHARMM22*, and experiment, for ubiquitin, GB3, and VHS. The oblate tumbling of ubiquitin about Dx instead of prolate tumbling about Dz (Figure 5) is revealed to be an under-

Figure 6. Comparison of GB3MoI (gray) and PAF (colored) axes, overlaid on top of cartoon backbone representations. Prepared similarly to Figure 5. (A) Coarse-grained 1IGD.pdb with experimental Drot.14 (B) Last frame of 2-μs-MARTINI simulations. (C) Last frame of 500 ns TIP3p/CHARMM simulations.

appears to limit the impact of structural fidelity upon Drot. On the basis of these observations, we suggest that MARTINI+EN broadly preserves the tumbling characteristics of folded proteins, with a caveat on the potential for coarse-graining artifacts that are likely repairable through stronger constraints. For three proteins with significant disorder, however, the artificial packing of unstructured elements strongly affects anisotropic tumbling behavior. We found that these elements either bind themselves stably against folded counterparts (plexin and VHS) or assemble into molten globules (in the case of UIM). The impact of this ordering upon Diso is surprising limited (Figure 3), which we suspect can be attributed to the greater dependence of isotropic tumbling upon the solute mass than to details such as mass distribution and solute−solvent interactions. Intuitively, one may argue that the artificial ordering reduces MARTINI-based tumbling to an equivalent rigid body; thus, by comparison with experiment one can deduce the net effect of disorder. A comparison in VHS (Figure 7) indicates that the slowest tumbling axis in experiment connects both N- and Ctermini, whose corresponding axis is fastest in MARTINI. This

Figure 5. Comparison of ubiquitin MoI (gray) and PAF (colored) axes, overlaid on top of cartoon backbone representations. Axis lengths have been normalized to 4Rgyr and scaled according to elements of the MoI and Drot tensors, respectively. Colors of axes: Dx (red), Dy (green), and Dz (blue). (A) Coarse-grained representation of 1D3Z.pdb, first frame, with experimental Drot.54 (B) Last frame of 2μs-MARTINI simulations. (C) Last frame of 500 ns TIP3p/ CHARMM simulations.

estimation of Dx and Dz relative to Dy. The MARTINI+EN Drot is rotated approximately −20° around Dx such that Dz is no longer pointed along the C-terminal tail. In contrast, the CHARMM22*/TIP3p ubiquitin tumbling replicates experimental Drot closely in both orientation and anisotropy. Closer inspection of MARTINI ubiquitin structure reveals a slight compression along the experimental Dz axis and minor shifts of

Figure 7. Comparison of VHS MoI (gray) and PAF (colored) axes, overlaid on top of cartoon backbone representations. Prepared similarly to Figure 5. MoI axes exclude unstructured tails. (A) Coarse-grained 1X5B.pdb with modeled tails and experimental Drot.54 (B) CG-MD simulations. The wrapping of both N- and C-terminal tails against the protein leads to incorrect estimation of Drot dimensions. G


Article

The Journal of Physical Chemistry B correlates with a previous analysis in Titin,9 where the authors investigated two variants with termini of different lengths. The presence of these termini proximal to the protein’s long-axis increased Δ by exerting frictional drag upon axes perpendicular to its relative location. We also confirmed this explicitly by adding a five-residue N-terminal tail to GB3 in all-atom simulations, which significantly increases the anisotropy in the direction of the attached tail (Table 1). While the above deductive reasoning can be leveraged to discuss the role of disorder despite force field limitations, it is no replacement for improved modeling. A dedicated effort to develop such force fields is desired, particularly in coarse-grained approaches so as to examine the potential impacts of disorder upon Drot convergence.

illustrate a potential extension of Drot descriptions toward fully disordered systems. Finally, the utilization of quaternions to calculate Drot provides a number of improvements over popular unit-vectorbased approaches. Global-tumbling derived from unit vectors is negatively affected by uneven sampling of vector orientations,21,70 and requires a computationally nontrivial calculation of autocorrelation functions for each selected unit vector. These autocorrelations are superfluous when a comparison with NMR-relaxation or order parameters are not required. In contrast, the quaternion approach removes internal motions visa-vis standard RMSD procedures, and requires only three or four autocorrelations from the components of q(Δt). Moreover, by prior elimination of vector-dependent internal motions, nonlinear tumbling phenomena can, in principle, be directly observed in q-autocorrelations, although this is in practice hidden by statistical errors from incomplete convergence. We note here that 2000 × τtheor appears to be sufficient to confirm that MARTINI+EN GB3 largely conserves rigid-body dynamics (Figure S4). In terms of Drot orientation, the PAF definition adopted in this paper is simplistic but pragmatic, as the lowest Δt data converge fastest and eliminate almost all internal dynamics that perturb rigid-body approximations. It is worth investigating more refined definitions that take into account motions on experimentally sensitive time scales for an improved comparison with, e.g., NMR-based approaches. In summary, an accurate computation of rotational diffusion advances our ability to interpret experimental measurements using structural and dynamical simulations, in particular for NMR-relaxation techniques where nanosecond dynamics can be equally probed by both approaches. Since ab initio predictions based on all-atom MD alone remain at present a computationally expensive task, we suggest the adoption of coarse-grained approaches to converge slower global motions. Using this, one can achieve full prediction of relaxation parameters by recombining global coarse-grained dynamics and local all-atom dynamics. This accelerates integrated NMR/MD investigations into the dynamics of large biomolecules, by also leveraging sparse labeling in NMR. Future improvements in force field and water models for disordered proteins will also unlock avenues toward the clarification of disordered−ordered interactions. We hope this work will inspire broader exploration of the dynamics−function relationship of biomolecules.

■

DISCUSSIONS AND SUMMARY The ongoing development of computational power and molecular force fields is now enabling routine simulations that stretch into microsecond time scales, capable of converging rotational tumbling calculations directly from MD trajectories. As a response to known Drot overestimation by popular TIP3p and SPC/E water models, we tested analogous calculations in a coarse-grained force field in order to examine its relative performance and improvements from more robust sampling alone. The MARTINI force field adopted for this task was found to predict experimental Diso to within ∼10% for a range of proteins up to 42 kDa (ovalbumin). This accuracy is achievable with a minimum sampling of >50 × τtheor; we note that an adequate prediction requires only that one achieves comprehensive sampling of rotations within the range of Δt considered. The anisotropy Δ values are reproduced with decreased precision, which can be improved with extensive sampling. Due to the strong influence of protein shape and disorder, however, we do not offer a single guidance value here. The rhombicity ρ cannot be precisely calculated due to overlap of Dx, Dy, and Dz uncertainties, and even a 10 μs simulation of GB3 produces an underwhelming ρ of 0.27 ± 0.27. We thus advise reporting of Dx, Dy, and Dz directly so as to expedite simulation−experiment comparisons and preserve comparable precision across elements of Drot. The particular sensitivity of Drot to alterations in molecular shape as well as disordered regions suggests its value as a constraint in both structure determination and force field refinement. In the former application, Drot functions similar to small-angle scattering, which measures orientationally averaged molecular shape plus the solvation layer. By comparing reconstructed scattering shapes with measured Drot, it may be possible to deduce the disordered contributions directly. In terms of force field refinement, Drot serves as independent validation of the bottom-up parameters describing atomic interactions within the force field. We reiterate that MARTINI performs surprisingly well in Drot calculations, given its usual application toward membranes and protein−lipid interactions. Development of parameters specific for disordered residues, and proposed improvements to MARTINI side-chain sampling,63 may establish a future application toward long-time scale soluble protein dynamics. It would be also of interest to test other solution-based coarse-grained force fields40,64−67 and approaches31 in this manner, as well as all-atom formulations specific to the study of disorder.68,69 Here, we note the experimental measurement of local tumbling times along different segments of α-synuclein by Parigi et al.,12 which

■

METHODS Force Fields and MD Software. The CHARMM22* force field71,72 and the MARTINI force field39,41 have been used for simulations in the GROMACS 5.1 suite.73,74 Additional constraints of MARTINI proteins with elastic networks were added, during which we took care to remove constraints between unstructured loops/termini and the folded domain. In coarse-grained simulations, Lennard-Jones and electrostatic potentials have been shifted between 0.9−1.2 nm and 0.0−1.2 nm, respectively, so as to vanish at the cutoff, via the Group cutoff scheme. The relative dielectric constant ϵrel was set to 15. For all-atom simulations, we adopt CHARMM-types standards of force-switching and PME electrostatics75 with dual-range 10−12 Å cutoffs. The Verlet cutoff scheme was adopted to take advantage of improved performance and GPU acceleration. Calculations of quaternion orientations from trajectories have been ported based on Fiorin et al.76 into a fork of PLUMED H


Article

The Journal of Physical Chemistry B 2,77 and currently hosted on GitHub (https://github.com/ zharmad/plumed2). The PLUMED implementation utilizes Kearsley’s algorithm78 to align protein frames, whose rotation matrix can be converted into the equivalent quaternion. A separate TcL-implementation for use in VMD79 is available on request. Simulation Protocol. Source coordinates for all proteins have been listed in Table 2, each of which was solvated and

4 ps were adopted, utilizing GROMACS’s virtual site construction. Ensembles were maintained by Nose−Hoover93 thermostats and Parrinello−Rahman barostats,92 setting τt = 2.5 ps, T = 298.15 K, τp = 5.0 ps, and β = 4.5 × 10−5. Separate thermobaths for proteins and solvent were also utilized in allatom simulations. Computational Resources. All simulations in this work were conducted on workstations equipped with desktop-level CPUs and GPUs. MARTINI simulations of all proteins span 500−2000 ns per day depending on system size, conducted without GPU acceleration. CHARMM simulations of proteins span 20−80 ns per day, conducted with GPU acceleration as provided in GROMACS 5.x. Theoretical Rotational Diffusion. τtheor values tabulated in Table 2 are computed by published equations given in the caption text. Reference/lab frames are chosen by aligning allatom energy-minimized proteins (below) according to their calculated moments of inertia, such that the longest axis points toward z. All orientations reported in this work follow this convention. For the reporting of rotational diffusion constants, values have been standardized to 293 K according to the relationship found in refs 20 and 27. This is equivalent to 17% decrease of D from 300 K. Comparison with Unit-Vector-Based Drot Calculations. To confirm that the quaternion-based method here yields consistent values with popular NMR-based approaches, we also calculated Drot for a majority of the systems investigated according to eqs 12 and 13, taking a fixed set of unit vectors attached to the body frame. Broad agreement was achieved except for several ρ values.

Table 2. Basic Statistics of Proteins Used in This Studya name (abbrev) Stam2 UIM-domain (UIM) GB3 domain (GB3) GB3 domain (GB3) + N-terminal tail ubiquitin (Ubq) Titin Ig27-domain (Titin) Plexin-B1 cytoplasm. domain (Plx) Stam2 VHS-domain (VHS) VHS-Ubq complex (VHS +Ubq) D-peroxiredoxin (Prx) monomer D-peroxiredoxin (Prx) dimer green fluorescent protein (GFP) GFP dimer ovalbumin (Ova) ovalbumin dimer aldolase dimer (Ald) aldolase tetramer (Ald) a

mol mass (kDa)

source (PDB)

Naa

ab initio 1P7F82 1IGD83

42 56 61

4.96 6.21 6.65

3.95 4.94 5.30

1XQQ84 1TIT85 2JPH86

76 89 122

8.56 9.78 13.4

6.83 7.81 10.8

1X5B 2L0T54

163 239

17.7 26.2

14.2 21.2

author55

161

16.9

13.6

55

322 229

33.8 25.8

27.4 20.8

458 387 757 682 1452

51.6 42.9 84.5 74.0 157

42.1 35.0 70.1 64.7 133

author

1GFL60 1OVA87 3BV488 1ZAH89

τtheor (ns)

■

ASSOCIATED CONTENT

S Supporting Information *

80

Molecular masses calculated via the ExPASY server using the amino acid sequence present in the simulation. The theoretical tumbling time τtheor is calculated via the Stokes−Einstein−Debye relation 1/2τ = Diso = kBT/8πηRhyd3, using a predicted viscosity of water at 293 K (η = 1.00156 × 10−3Pa s),20 and a hydrodynamic radius of the equivalent sphere inferred from the molecular mass m and a predicted density curve ρ, density = 1.4106 + 0.14528e−m/134 g cm−3 proposed by Fischer et al.81

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.6b11703. Derivation of Drot, graphical depiction of Δ and ρ, Drot analysis, PAF orientation, qPB, for all proteins, and computed Drot (PDF) Reference coordinates for all proteins (ZIP) Animated version of Figure 1 (MPG)

■

energy-minimized in CHARMM22* before coarse-graining. The coarse-graining procedure has been previously described,41 and utilizes the martinize.py script publicly available, which is invoked with elastic networks to limit internal motions. Secondary structure for coarse-graining were determined via the DSSP 2.0 server.90 Coarse-grained coordinates were briefly minimized in vacuum, solvated in dodecahedron box with 1.5 nm margins and periodic boundary conditions, ionized to 150 mM, and then equilibrated in 200 ps NVT followed by 200 ps NPT simulations. Time steps of 20 fs were adopted, while ensembles were maintained through velocity-rescaling thermostats91 and Parrinello−Rahman barostats,92 setting τt = 1.0 ps, T = 300 K, τp = 5.0 ps, and β = 4.5 × 10−5. Proteins and solvents have been coupled to separate baths to prevent spontaneous temperature gradient formation. Output trajectories were written once every 10 ps. All-atom coordinates were first solvated and ionized to equivalent parameters as coarse-grained simulations above, in standard TIP3p water. This was followed by energy minimization, and gradual relaxation under side-chain and backbone restraints over 2 ns in NPT conditions. Time steps of

AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected]. *E-mail: [email protected]. ORCID

Po-chia Chen: 0000-0003-2272-2082 Present Address †

EMBL Heidelberg, Meyerhofstrasse 1, 69126 Heidelberg, Germany. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work was supported by a grant overseen by the French National Research Agency (ANR) as part of the “OH risque” programme (Metadyn, ANR-14-OHRI-0006-01).

■

REFERENCES

(1) Zhang, N.; Wang, Q.; Ehlinger, A.; Randles, L.; Lary, J. W.; Kang, Y.; Haririnia, A.; Storaska, A. J.; Cole, J. L.; Fushman, D.; et al.

I


Article


(20) Garciá de la Torre, J.; Huertas, M. L.; Carrasco, B. HYDRONMR: Prediction of NMR Relaxation of Globular Proteins from Atomic-Level Structures and Hydrodynamic Calculations. J. Magn. Reson. 2000, 147, 138−146. (21) Walker, O.; Varadan, R.; Fushman, D. Efficient and Accurate Determination of the Overall Rotational Diffusion Tensor of a Molecule from 15N Relaxation Data Using Computer Program ROTDIF. J. Magn. Reson. 2004, 168, 336−345. (22) Lipari, G.; Szabo, A. Model-Free Approach to the Interpretation of Nuclear Magnetic Resonance Relaxation in Macromolecules. 1. Theory and Range of Validity. J. Am. Chem. Soc. 1982, 104, 4546− 4559. (23) d’Auvergne, E. J.; Gooley, P. R. Set Theory Formulation of the Model−Free Problem and the Diffusion Seeded Model−Free Paradigm. Mol. BioSyst. 2007, 3, 483−494. (24) Smith, P. E.; van Gunsteren, W. F. Translational and Rotational Diffusion of Proteins. J. Mol. Biol. 1994, 236, 629−636. (25) Peter, C.; Daura, X.; van Gunsteren, W. F. Calculation of NMRRelaxation Parameters for Flexible Molecules from Molecular Dynamics Simulations. J. Biomol. NMR 2001, 20, 297−310. (26) Nederveen, A. J.; Bonvin, A. M. J. J. NMR Relaxation and Internal Dynamics of Ubiquitin from a 0.2 Ms MD Simulation. J. Chem. Theory Comput. 2005, 1, 363−374. (27) Wong, V.; Case, D. A. Evaluating Rotational Diffusion from Protein MD Simulations. J. Phys. Chem. B 2008, 112, 6013−6024. (28) Anderson, J. S.; LeMaster, D. M. Rotational Velocity Rescaling of Molecular Dynamics Trajectories for Direct Prediction of Protein NMR Relaxation. Biophys. Chem. 2012, 168−169, 28−39. (29) Takemura, K.; Kitao, A. Water Model Tuning for Improved Reproduction of Rotational Diffusion and NMR Spectral Density. J. Phys. Chem. B 2012, 116, 6279−6287. (30) Zhang, L.; Bouguet-Bonnet, S.; Buck, M. In Allostery; Fenton, A. W., Ed.; Methods in Molecular Biology 796; Springer: New York, 2012; pp 235−259. (31) Copperman, J.; Guenza, M. G. Coarse-Grained Langevin Equation for Protein Dynamics: Global Anisotropy and a Mode Approach to Local Complexity. J. Phys. Chem. B 2015, 119, 9195− 9211. (32) Trbovic, N.; Kim, B.; Friesner, R. A.; Palmer, A. G. Structural Analysis of Protein Dynamics by MD Simulations and NMR SpinRelaxation. Proteins: Struct., Funct., Genet. 2008, 71, 684−694. (33) Maragakis, P.; Lindorff-Larsen, K.; Eastwood, M. P.; Dror, R. O.; Klepeis, J. L.; Arkin, I. T.; Jensen, M. O.; Xu, H.; Trbovic, N.; Friesner, R. A.; et al. Microsecond Molecular Dynamics Simulation Shows Effect of Slow Loop Dynamics on Backbone Amide Order Parameters of Proteins. J. Phys. Chem. B 2008, 112, 6155−6158. (34) Lindorff-Larsen, K.; Trbovic, N.; Maragakis, P.; Piana, S.; Shaw, D. E. Structure and Dynamics of an Unfolded Protein Examined by Molecular Dynamics Simulation. J. Am. Chem. Soc. 2012, 134, 3787− 3791. (35) Robustelli, P.; Trbovic, N.; Friesner, R. A.; Palmer, A. G. Conformational Dynamics of the Partially Disordered Yeast Transcription Factor GCN4. J. Chem. Theory Comput. 2013, 9, 5190−5200. (36) Gu, Y.; Li, D.-W.; Brüschweiler, R. NMR Order Parameter Determination from Long Molecular Dynamics Trajectories for Objective Comparison with Experiment. J. Chem. Theory Comput. 2014, 10, 2599−2607. (37) González, M. A.; Abascal, J. L. F. The Shear Viscosity of Rigid Water Models. J. Chem. Phys. 2010, 132, 096101. (38) Bond, P. J.; Holyoake, J.; Ivetac, A.; Khalid, S.; Sansom, M. S. P. Coarse-Grained Molecular Dynamics Simulations of Membrane Proteins and Peptides. J. Struct. Biol. 2007, 157, 593−605. (39) Marrink, S. J.; Risselada, H. J.; Yefimov, S.; Tieleman, D. P.; de Vries, A. H. The MARTINI Force Field: Coarse Grained Model for Biomolecular Simulations. J. Phys. Chem. B 2007, 111, 7812−7824. (40) Pasi, M.; Lavery, R.; Ceres, N. PaLaCe: A Coarse-Grain Protein Model for Studying Mechanical Properties. J. Chem. Theory Comput. 2013, 9, 785−793.

Structure of the S5a:K48-Linked Diubiquitin Complex and Its Interactions with Rpn13. Mol. Cell 2009, 35, 280−290. (2) Yang, S.; Blachowicz, L.; Makowski, L.; Roux, B. Multidomain Assembled States of Hck Tyrosine Kinase in Solution. Proc. Natl. Acad. Sci. U. S. A. 2010, 107, 15757−15762. (3) Munari, F.; Rezaei-Ghaleh, N.; Xiang, S.; Fischle, W.; Zweckstetter, M. Structural Plasticity in Human Heterochromatin Protein 1β. PLoS One 2013, 8, e60887. (4) Dayel, M. J.; Hom, E. F. Y.; Verkman, A. S. Diffusion of Green Fluorescent Protein in the Aqueous-Phase Lumen of Endoplasmic Reticulum. Biophys. J. 1999, 76, 2843−2851. (5) Balbo, J.; Mereghetti, P.; Herten, D.-P.; Wade, R. C. The Shape of Protein Crowders is a Major Determinant of Protein Diffusion. Biophys. J. 2013, 104, 1576−1584. (6) Göbl, C.; Madl, T.; Simon, B.; Sattler, M. NMR Approaches for Structural Analysis of Multidomain Proteins and Complexes in Solution. Prog. Nucl. Magn. Reson. Spectrosc. 2014, 80, 26−63. (7) Bouguet-Bonnet, S.; Buck, M. Compensatory and Long-Range Changes in Picosecond−Nanosecond Main-Chain Dynamics upon Complex Formation: 15N Relaxation Analysis of the Free and Bound States of the Ubiquitin-like Domain of Human Plexin-B1 and the Small GTPase Rac1. J. Mol. Biol. 2008, 377, 1474−1487. (8) Lange, A.; Castañeda, C.; Hoeller, D.; Lancelin, J.-M.; Fushman, D.; Walker, O. Evidence for Cooperative and Domain-Specific Binding of the Signal Transducing Adaptor Molecule 2 (STAM2) to Lys63Linked Diubiquitin. J. Biol. Chem. 2012, 287, 18687−18699. (9) Nicastro, G.; Margiocco, P.; Cardinali, B.; Stagnaro, P.; Cauglia, F.; Cuniberti, C.; Collini, M.; Thomas, D.; Pastore, A.; Rocco, M. The Role of Unstructured Extensions in the Rotational Diffusion Properties of a Globular Protein: The Example of the Titin I27 Module. Biophys. J. 2004, 87, 1227−1240. (10) Bae, S.-H.; Dyson, H. J.; Wright, P. E. Prediction of the Rotational Tumbling Time for Proteins with Disordered Segments. J. Am. Chem. Soc. 2009, 131, 6814−6821. (11) Zerbetto, M.; Buck, M.; Meirovitch, E.; Polimeno, A. Integrated Computational Approach to the Analysis of NMR Relaxation in Proteins: Application to Ps-ns Main Chain 15N-1H and Global Dynamics of the Rho GTPase Binding Domain of Plexin-B1. J. Phys. Chem. B 2011, 115, 376−388. (12) Parigi, G.; Rezaei-Ghaleh, N.; Giachetti, A.; Becker, S.; Fernandez, C.; Blackledge, M.; Griesinger, C.; Zweckstetter, M.; Luchinat, C. Long-Range Correlated Dynamics in Intrinsically Disordered Proteins. J. Am. Chem. Soc. 2014, 136, 16201−16209. (13) Dosset, P.; Hus, J.-C.; Blackledge, M.; Marion, D. Efficient Analysis of Macromolecular Rotational Diffusion from Heteronuclear Relaxation Data. J. Biomol. NMR 2000, 16, 23−28. (14) Berlin, K.; Longhini, A.; Dayie, T. K.; Fushman, D. Deriving Quantitative Dynamics Information for Proteins and RNAs Using ROTDIF with a Graphical User Interface. J. Biomol. NMR 2013, 57, 333−352. (15) Rezaei-Ghaleh, N.; Klama, F.; Munari, F.; Zweckstetter, M. Predicting the Rotational Tumbling of Dynamic Multidomain Proteins and Supramolecular Complexes. Angew. Chem., Int. Ed. 2013, 52, 11410−11414. (16) Jovin, T. M.; Bartholdi, M.; Vaz, W. L. C.; Austin, R. H. Rotational Diffusion of Biological Macromolecules by Time-Resolved Delayed Luminescence (Phosphorescence, Fluorescence) Anisotropy. Ann. N. Y. Acad. Sci. 1981, 366, 176−196. (17) Borst, J. W.; Laptenok, S. P.; Westphal, A. H.; Kühnemuth, R.; Hornen, H.; Visser, N. V.; Kalinin, S.; Aker, J.; van Hoek, A.; Seidel, C. A. M.; et al. Structural Changes of Yellow Cameleon Domains Observed by Quantitative FRET Analysis and Polarized Fluorescence Correlation Spectroscopy. Biophys. J. 2008, 95, 5399−5411. (18) Loman, A.; Gregor, I.; Stutz, C.; Mund, M.; Enderlein, J. Measuring Rotational Diffusion of Macromolecules by Fluorescence Correlation Spectroscopy. Photochem. Photobiol. Sci. 2010, 9, 627−636. (19) Lee, L. K.; Rance, M.; Chazin, W. J.; Palmer, A. G. Rotational Diffusion Anisotropy of Proteins from Simultaneous Analysis of 15N and 13Cα Nuclear Spin Relaxation. J. Biomol. NMR 1997, 9, 287−298. J


Article


(63) Herzog, F. A.; Braun, L.; Schoen, I.; Vogel, V. Improved Side Chain Dynamics in MARTINI Simulations of Protein-Lipid Interfaces. J. Chem. Theory Comput. 2016, 12, 2446−2458. (64) Hills, R. D., Jr; Lu, L.; Voth, G. A. Multiscale Coarse-Graining of the Protein Energy Landscape. PLoS Comput. Biol. 2010, 6, e1000827. (65) Kar, P.; Gopal, S. M.; Cheng, Y.-M.; Predeus, A.; Feig, M. PRIMO: A Transferable Coarse-Grained Force Field for Proteins. J. Chem. Theory Comput. 2013, 9, 3769−3788. (66) Emperador, A.; Sfriso, P.; Villarreal, M. A.; Gelpí, J. L.; Orozco, M. PACSAB: Coarse-Grained Force Field for the Study of ProteinProtein Interactions and Conformational Sampling in Multiprotein Systems. J. Chem. Theory Comput. 2015, 11, 5929−5938. (67) Haxton, T. K. High-Resolution Coarse-Grained Modeling Using Oriented Coarse-Grained Sites. J. Chem. Theory Comput. 2015, 11, 1244−1254. (68) Rauscher, S.; Gapsys, V.; Gajda, M. J.; Zweckstetter, M.; de Groot, B. L.; Grubmüller, H. Structural Ensembles of Intrinsically Disordered Proteins Depend Strongly on Force Field: A Comparison to Experiment. J. Chem. Theory Comput. 2015, 11, 5513−5524. (69) Henriques, J.; Cragnell, C.; Skepö, M. Molecular Dynamics Simulations of Intrinsically Disordered Proteins: Force Field Evaluation and Comparison with Experiment. J. Chem. Theory Comput. 2015, 11, 3420−3431. (70) Fushman, D.; Ghose, R.; Cowburn, D. The Effect of Finite Sampling on the Determination of Orientational Properties: A Theoretical Treatment with Application to Interatomic Vectors in Proteins. J. Am. Chem. Soc. 2000, 122, 10640−10649. (71) Brooks, B. R.; Brooks, C. L.; Mackerell, A. D.; Nilsson, L.; Petrella, R. J.; Roux, B.; Won, Y.; Archontis, G.; Bartels, C.; Boresch, S.; et al. CHARMM: The Biomolecular Simulation Program. J. Comput. Chem. 2009, 30, 1545−1614. (72) Piana, S.; Lindorff-Larsen, K.; Shaw, D. E. How Robust Are Protein Folding Simulations with Respect to Force Field Parameterization? Biophys. J. 2011, 100, L47−L49. (73) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 435−447. (74) Bjelkmar, P.; Larsson, P.; Cuendet, M. A.; Hess, B.; Lindahl, E. Implementation of the CHARMM Force Field in GROMACS: Analysis of Protein Stability Effects from Correction Maps, Virtual Interaction Sites, and Water Models. J. Chem. Theory Comput. 2010, 6, 459−466. (75) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577−8593. (76) Fiorin, G.; Klein, M. L.; Hénin, J. Using Collective Variables to Drive Molecular Dynamics Simulations. Mol. Phys. 2013, 111, 3345− 3362. (77) Tribello, G. A.; Bonomi, M.; Branduardi, D.; Camilloni, C.; Bussi, G. PLUMED 2: New Feathers for an Old Bird. Comput. Phys. Commun. 2014, 185, 604−613. (78) Kearsley, S. K. On the Orthogonal Transformation Used for Structural Comparisons. Acta Crystallogr., Sect. A: Found. Crystallogr. 1989, 45, 208−210. (79) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 33−38. (80) Walker, J. M. The Proteomics Protocols Handbook; Humana Press: Totowa, NJ, 2005; OCLC 272404489. (81) Fischer, H.; Polikarpov, I.; Craievich, A. F. Average Protein Density Is a Molecular-Weight-Dependent Function. Protein Sci. 2004, 13, 2825−2828. (82) Ulmer, T. S.; Ramirez, B. E.; Delaglio, F.; Bax, A. Evaluation of Backbone Proton Positions and Dynamics in a Small Protein by Liquid Crystal NMR Spectroscopy. J. Am. Chem. Soc. 2003, 125, 9179−9191. (83) Derrick, J. P.; Wigley, D. B. The Third IgG-Binding Domain from Streptococcal Protein G. J. Mol. Biol. 1994, 243, 906−918. (84) Lindorff-Larsen, K.; Best, R. B.; DePristo, M. A.; Dobson, C. M.; Vendruscolo, M. Simultaneous Determination of Protein Structure and Dynamics. Nature 2005, 433, 128−132.

(41) Monticelli, L.; Kandasamy, S. K.; Periole, X.; Larson, R. G.; Tieleman, D. P.; Marrink, S.-J. The MARTINI Coarse-Grained Force Field: Extension to Proteins. J. Chem. Theory Comput. 2008, 4, 819− 834. (42) Delong, S.; Usabiaga, F. B.; Donev, A. Brownian Dynamics of Confined Rigid Bodies. J. Chem. Phys. 2015, 143, 144107. (43) Chevrot, G.; Hinsen, K.; Kneller, G. R. Model-Free Simulation Approach to Molecular Diffusion Tensors. J. Chem. Phys. 2013, 139, 154110. (44) Furry, W. H. Isotropic Rotational Brownian Motion. Phys. Rev. 1957, 107, 7−13. (45) Favro, L. D. Theory of the Rotational Brownian Motion of a Free Rigid Body. Phys. Rev. 1960, 119, 53−62. (46) Diebel, J. Representing Attitude: Euler Angles, Unit Quaternions, and Rotation Vectors; 2006; https://www.astro.rug.nl/software/ kapteyn/?downloads/attitude.pdf (accessed Jan 13, 2017). (47) Coutsias, E. A.; Seok, C.; Dill, K. A. Using Quaternions to Calculate RMSD. J. Comput. Chem. 2004, 25, 1849−1857. (48) Dantam, N. Quaternion Computation; 2014; http://www.neil. dantam.name/note/dantam-quaternion.pdf (accessed Jan 13, 2017). (49) Diamond, R. On the Multiple Simultaneous Superposition of Molecular Structures by Rigid Body Transformations. Protein Sci. 1992, 1, 1279−1287. (50) Tjandra, N.; Feller, S. E.; Pastor, R. W.; Bax, A. Rotational Diffusion Anisotropy ofHuman Ubiquitin from 15N NMR Relaxation. J. Am. Chem. Soc. 1995, 117, 12562−12566. (51) Miles, R. E. On Random Rotations in R3. Biometrika 1965, 52, 636−639. (52) Huntress, W. T. Effects of Anisotropic Molecular Rotational Diffusion on Nuclear Magnetic Relaxation in Liquids. J. Chem. Phys. 1968, 48, 3524−3533. (53) Bruschweiler, R.; Liao, X.; Wright, P. E. Long-Range Motional Restrictions in a Multidomain Zinc-Finger Protein from Anisotropic Tumbling. Science 1995, 268, 886−889. (54) Lange, A.; Hoeller, D.; Wienk, H.; Marcillat, O.; Lancelin, J.-M.; Walker, O. NMR Reveals a Different Mode of Binding of the Stam2 VHS Domain to Ubiquitin and Diubiquitin. Biochemistry 2011, 50, 48−62. (55) Echalier, A.; Trivelli, X.; Corbier, C.; Rouhier, N.; Walker, O.; Tsan, P.; Jacquot, J.-P.; Aubry, A.; Krimm, I.; Lancelin, J.-M. Crystal Structure and Solution NMR Dynamics of a D (Type II) Peroxiredoxin Glutaredoxin and Thioredoxin Dependent: A New Insight into the Peroxiredoxin Oligomerism. Biochemistry 2005, 44, 1755−1767. (56) Swaminathan, R.; Hoang, C. P.; Verkman, A. S. Photobleaching Recovery and Anisotropy Decay of Green Fluorescent Protein GFPS65T in Solution and Cells: Cytoplasmic Viscosity Probed by Green Fluorescent Protein Translational and Rotational Diffusion. Biophys. J. 1997, 72, 1900−1907. (57) Pieper, C. M.; Enderlein, J. Fluorescence Correlation Spectroscopy as a Tool for Measuring the Rotational Diffusion of Macromolecules. Chem. Phys. Lett. 2011, 516, 1−11. (58) Fuhrmans, M.; Sanders, B. P.; Marrink, S.-J.; de Vries, A. H. Effects of Bundling on the Properties of the SPC Water Model. Theor. Chem. Acc. 2010, 125, 335−344. (59) Marrink, S. J.; Tieleman, D. P. Perspective on the Martini Model. Chem. Soc. Rev. 2013, 42, 6801−6822. (60) Yang, F.; Moss, L. G.; Phillips, G. N. The Molecular Structure of Green Fluorescent Protein. Nat. Biotechnol. 1996, 14, 1246−1251. (61) Ianeselli, L.; Zhang, F.; Skoda, M. W. A.; Jacobs, R. M. J.; Martin, R. A.; Callow, S.; Prévost, S.; Schreiber, F. Protein-Protein Interactions in Ovalbumin Solutions Studied by Small-Angle Scattering: Effect of Ionic Strength and the Chemical Nature of Cations. J. Phys. Chem. B 2010, 114, 3776−3783. (62) Lu, C.-Y.; Vanden Bout, D. A. Effect of Finite Trajectory Length on the Correlation Function Analysis of Single Molecule Data. J. Chem. Phys. 2006, 125, 124701. K


Article

The Journal of Physical Chemistry B (85) Improta, S.; Politou, A. S.; Pastore, A. Immunoglobulin-like Modules from Titin I-Band: Extensible Components of Muscle Elasticity. Structure 1996, 4, 323−337. (86) Tong, Y.; Hota, P. K.; Hamaneh, M. B.; Buck, M. Insights into Oncogenic Mutations of Plexin-B1 Based on the Solution Structure of the Rho GTPase Binding Domain. Structure 2008, 16, 246−258. (87) Stein, P. E.; Leslie, A. G. W.; Finch, J. T.; Carrell, R. W. Crystal Structure of Uncleaved Ovalbumin at 1·95 Å Resolution. J. Mol. Biol. 1991, 221, 941−959. (88) Sherawat, M.; Tolan, D. R.; Allen, K. N. Structure of a Rabbit Muscle Fructose-1,6-Bisphosphate Aldolase A Dimer Variant. Acta Crystallogr., Sect. D: Biol. Crystallogr. 2008, 64, 543−550. (89) St-Jean, M.; Lafrance-Vanasse, J.; Liotard, B.; Sygusch, J. High Resolution Reaction Intermediates of Rabbit Muscle Fructose-1,6Bisphosphate Aldolase Substrate Cleavage and Induced Fit. J. Biol. Chem. 2005, 280, 27262−27270. (90) Kabsch, W.; Sander, C. Dictionary of Protein Secondary Structure: Pattern Recognition of Hydrogen-Bonded and Geometrical Features. Biopolymers 1983, 22, 2577−2637. (91) Bussi, G.; Donadio, D.; Parrinello, M. Canonical Sampling Through Velocity Rescaling. J. Chem. Phys. 2007, 126, 014101. (92) Parrinello, M.; Rahman, A. Polymorphic Transitions in Single Crystals: A New Molecular Dynamics Method. J. Appl. Phys. 1981, 52, 7182−7190. (93) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A: At., Mol., Opt. Phys. 1985, 31, 1695−1697.

L


Computing the Rotational Diffusion of ... - ACS Publications

Recommend Documents