Molecular Dynamics of a Hydrated Collagen Peptide - ACS Publications

Nov 16, 2016 - to as the magic angle (MA) effect and has been attributed to the restricted rotational .... the collagen peptide and the edge of the wa...
0 downloads 0 Views 4MB Size
Article pubs.acs.org/JPCB

Molecular Dynamics of a Hydrated Collagen Peptide: Insights into Rotational Motion and Residence Times of Single-Water Bridges in Collagen Monique C. Tourell†,‡ and Konstantin I. Momot*,†,‡ †

School of Chemistry, Physics and Mechanical Engineering and ‡Institute of Health and Biomedical Innovation, Queensland University of Technology (QUT), GPO Box 2434, Brisbane, Queensland 4001, Australia S Supporting Information *

ABSTRACT: Magnetic resonance transverse spin relaxation time constants (T2) of water protons in ordered collagenous tissues are dependent on the orientation of the tissue relative to the static magnetic field. This dependence is commonly referred to as the magic angle (MA) effect and has been attributed to the restricted rotational motion of icelike water bridges in the hydrated triple-helix collagen molecule. Understanding of the molecular mechanism of the MA effect is important for clinical and research applications of magnetic resonance spectroscopy and imaging to tissues, such as articular cartilage, tendons, and ligaments. In this work, we have used molecular dynamics simulations to investigate the subnanosecond time scale dynamics of single-water bridges in a model collagen peptide. We ascertain the residence times and the patterns of restricted rotational motion of water molecules. The key findings are strongly anisotropic rotation patterns of water molecules at bridge sites and a dynamic, rather than icelike, nature of the single-water bridges within the individual triple-helix collagen molecule.



INTRODUCTION Transverse relaxation anisotropy, or the magic angle (MA) effect, is a commonly observed phenomenon in magnetic resonance imaging (MRI) of ordered collagenous tissues, such as tendon,1,2 articular cartilage,3 and the intervertebral disk.4 The transverse spin relaxation time constants (T2) of the water 1 H nuclei in these tissues are strongly dependent on the orientation of the tissue relative to the applied static magnetic field B0. Although this dependence can complicate the interpretation of T2-weighted MR images,5 it can also provide valuable insight into the tissue microstructure.6−8 T2 anisotropy is usually attributed to the restricted rotation of water molecules associated with oriented collagen fibers, which results in residual dipolar couplings. Rotational modulation of the dipolar interactions between nuclear spins renders the time-dependent part of the spin Hamiltonian stochastic, resulting in a loss of spin-state coherence and the consequent transverse spin relaxation. The dipolar Hamiltonian HDD describes the interaction between the spin pairs; for two like spins I and S (e.g., two water protons), the secular approximation of HDD is given by9,10 HDD =

μ0 ℏγ 2 4πrIS 3

(1 − 3 cos2 θIS)(3IzSz − I·S)

vector and B0, and I and S are the spin operators with z components given by Iz and Sz. In isotropic solution, dipolar interactions between spin pairs are spatially averaged to zero on the time scale of the MRI experiment by rapid molecular tumbling. In collagenous tissues, water bound to collagen molecules undergoes restricted rotational motion, and the interproton vector is averaged to some residual value corresponding to the ensemble-average alignment. The contribution of HDD to spin relaxation is maximized when the director of the interproton vectors in the ensemble is aligned with B0 (⟨θIS⟩ = 0), resulting in reduced T2 values. Conversely, the HDD is minimized, and T2 values are maximized, when the director of the interproton vectors between spin pairs makes the MA ⟨θIS⟩ = θMA ≈ 54.7° with B0. This dependence has been observed experimentally in aligned collagenous tissues.11−13 In particular, maximum T2 values in tendon occur when the collagen fibers make an angle close to θMA with B0.2,12,13 This indicates that the anisotropic component of spin relaxation scales as (1 − 3 cos2 θ)2, where θ is the angle between the collagen fibers and B0, and that the director of the interproton vectors in the ensemble is coaxial with the axis of the collagen fibers.

(1)

Received: August 23, 2016 Revised: November 15, 2016 Published: November 16, 2016

where γ is the gyromagnetic ratio of the spins, rIS is the distance between spins I and S, θIS is the angle between the interproton © XXXX American Chemical Society

A

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

Article

The Journal of Physical Chemistry B

Figure 1. Nine single-water bridge sites, labeled from B1 (bridge site 1) to B9 (bridge site 9), in the imino-poor regions of 1BKV. The first chain is repeated for clarity.

molecules around the bridge axis, while allowing the protons that do not directly participate in the bridge bonds to undergo chemical exchange; (2) the time-averaged intramolecular proton−proton vectors over all bridge sites are coaxial with the axis of the triple-helix collagen molecule; and (3) intramolecular 1H−1H dipolar interaction within the bridge water molecule imparts a stochastic phase shift (relative to the bulk water protons) to the exchangeable proton, which cannot be reversed by spin-echo or Carr−Purcell−Meiboom−Gill pulse sequences and therefore results in transverse spin relaxation. Single-water bridges have been observed in molecular dynamics (MD) simulations of collagen peptides18,19 and in a larger crystal structure.20 However, in these studies, the water bridges were described as dynamic, not icelike. These previous studies were concerned with the dynamics of the collagen molecule and did not investigate the rotational dynamics of the water during the bridge lifetimes. In the present work, we

Fullerton and Rahal have proposed that the primary source of anisotropic water rotation in collagenous tissues, resulting in the MA effect, are single and double-water bridges between positive N−H and negative CO groups of the triple-helix collagen backbone.1 These bridges, as well as longer bridges in the larger hydration shell, can be directly observed in X-ray crystallography of collagen peptides.14 Additionally, differential scanning calorimetry (DSC) studies suggest their behavior to be icelike.15,16 Fullerton has recently provided a detailed review of the role of icelike water bridges as the source of the MA effect in aligned collagenous tissues.17 The definition of an icelike water bridge, used by Fullerton in this context, indicates the hydrogen bonds between the collagen atoms and the water molecule (the molecule making the water bridge) remain intact for at least the duration of the MRI experiment; this implies the hydrogen bond lifetimes to be of the order of seconds or longer. The three essential elements of Fullerton’s model are: (1) icelike water bridges that restrict the rotation of water B

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

Article

The Journal of Physical Chemistry B

Figure 2. Distances between collagen bridge atoms and the closest water molecules at bridge site 1: (a, b) distances between the oxygen in the closest water molecule and the amide hydrogen, dN−H···O; (c, d) distances between both hydrogens in the closest water molecule and the carboxyl oxygen, dCO···H. The histograms in (b) and (d) compare the distance distributions for the closest water molecule to the bridge site (solid histograms) to the distance distribution when a strict bridge is formed (black lines). In the strict bridge distributions, the probability density at distances greater than 2.4 Å (the hydrogen bond cutoff distance) is zero.

the position of the backbone atoms in the collagen molecule. The system was equilibrated at 310 K for 500 ps with the backbone constraint k = 1 kcal mol−1 Å−2. Constraints were then slowly reduced to 0.1 kcal mol−1 Å−2 and the system reequilibrated for a further 500 ps. Harmonic constraints were applied to the positions of the backbone atoms of the 1BKV molecule throughout the equilibration and the simulation to maintain a fiberlike state and prevent any unwinding in the imino-poor region and movement of the collagen molecule in the simulation volume. Both these complications were observed in previous simulations of hydrated 1BKV by Ravikumar and Hwang.18 The final four production runs applied constraints of 0.1, 1, 5, and 10 kcal mol−1 Å−2 to the collagen backbone. Longrange electrostatic interactions were treated using the particle mesh Ewald method. A cutoff distance of 10 Å was used with switching functions applied to nonbonded interactions from 8 Å. The pressure (1 atm) and the temperature (310 K) were controlled using the Nosé−Hoover method30 with Langevin dynamics used to control fluctuations31 and the Langevin thermostat,32 respectively. Covalent bonds between hydrogens and heavy atoms were constrained using SHAKE,33 and hydrogen bonds in water were constrained using SETTLE.34 A simulation time step of 1 fs was used, and the total simulation time was 10 ns for each of the four production runs. Atomic positions were recorded every 0.1 ps (100 time steps). Trajectory Analysis. Analysis of the simulation trajectories was carried out using in-house C++ and Mathematica (Wolfram Inc) code. We investigated five different water hydrogen bonding states at each bridge site. The first two states do not use strict criteria to determine hydrogen bonds: (1) the closest water molecule to the bridge site; and (2) the closest hydrogen (within the closest water molecule) to the carboxyl oxygen bridge atom. For each simulation time frame, the water molecule closest to each of the nine bridge sites was determined as the molecule with the smallest quantity

performed explicit MD simulations of water hydrating a synthetic collagen peptide containing a segment of type III collagen. We characterize the rotational motion and residence times of water molecules at and near the single-water bridge sites to gain insights into the molecular basis of the MA effect. The applications of this work also extend to other MR imaging modalities used to characterize collagenous tissues, such as diffusion imaging, where knowledge of the residence times of water on a collagen molecule is important for understanding the role of hydrogen bonding in the restriction of translational diffusion of water in these tissues.21



METHODS MD Simulations. A synthetic collagen peptide containing a region of human type III collagen (PDB ID: 1BKV22) was used for the simulations. This collagen peptide has been used in previous MD simulations,23 including one predominantly concerned with water bridges18 as well as in a coarse-grained Langevin dynamics study of the translation diffusion of water in model articular cartilage.24 The peptide 1BKV comprises three identical helices containing 30 residues each with 9 single-water bridge sites between them (Figure 1). We refer the bridge sites from 1 to 9 in the same manner as in Ravikumar and Hwang18 Presimulation setup of the collagen peptide was done using visual molecular dynamics.25 Acetic acid molecules, used for crystallization, were removed from the PDB structure and the missing hydrogen atoms were added. The collagen helix was aligned along the z axis of the simulation volume by applying a rotation to the atomic coordinates. The resulting molecule was solvated in a water box of TIP3P water molecules. The box was constructed such that there was a distance of ∼15 Å between the collagen peptide and the edge of the water box in all directions. The MD simulations were performed using NAMD26 with the CHARMM all-atom force field param2227 that was altered to include the parameters for hydroxyproline.28 This parameter set has been successfully used in MD simulations of 1BKV18,23 and other collagen peptides.29 After energy minimization, the system was gradually heated to 310 K at the rate of 10 K/ps, with a constraint constant of k = 1 kcal mol−1 Å−2 applied to

d 2 = dN − H ··· O2 + dC = O ··· H 2

(2)

where dN−H···O is the distance from the water oxygen to the amino hydrogen bridge atom. The quantity dCO···H is the C

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

Article

The Journal of Physical Chemistry B Table 1. Residence Times and Occupancies for the Five Different Bond States across all the Nine Bridge Sitesa eq 2

strict bonds

(1)

closest water

(2)

closest hydrogen

(3)

N−H···O

(4)

CO···H

(5)

bridge

mean residence time (ps) time-weighted median residence time (ps) mean residence time (ps) time-weighted median residence time (ps) occupancy mean residence time (ps) time-weighted median residence time (ps) occupancy mean residence time (ps) time-weighted median residence time (ps) occupancy mean residence time (ps) time-weighted median residence time (ps)

1

2

3

4

5

6

7

8

9

18.79 132.7

9.39 78.5

6.79 101.2

3.69 36.9

4.65 37.4

2.96 21.3

9.30 63.6

23.81 144.7

17.92 140.0

8.54 45.3

5.32 36.1

4.35 46.0

2.34 18.4

2.67 15.0

1.76 9.8

4.10 20.2

5.63 22.0

4.43 17.9

0.80 0.86 1.6

0.75 1.05 2.2

0.69 0.56 1.0

0.82 1.21 2.4

0.83 1.19 2.4

0.82 1.13 2.2

0.88 1.64 3.2

0.92 2.18 4.5

0.88 1.72 3.7

0.89 2.56 6.1

0.91 2.36 5.3

0.87 2.40 6.2

0.76 1.29 2.9

0.77 1.38 3.3

0.72 1.18 3.0

0.82 1.70 4.3

0.81 1.55 3.8

0.77 1.24 3.1

0.72 0.68 1.2

0.69 0.76 1.4

0.60 0.48 0.8

0.61 0.72 1.2

0.63 0.72 1.2

0.57 0.70 1.2

0.72 0.92 1.6

0.74 0.94 1.8

0.67 0.72 1.3

a Occupancy refers to the total fraction of the simulation time spent in the given state; residence times are defined as the expectation values of the continuous time spent in the state. Residence times and occupancies for the single strict bonds between a water molecule and the carboxyl oxygen or amide hydrogen include times when a strict bridge was formed. By definition, the occupancies for the closest water and closest hydrogen states are one for all nine bridges.

Figure 2a. Similarly, the distance between both hydrogens in the closest water molecules and the carboxyl oxygen for bridge site 1 is shown in Figure 2b. The distributions of distances in Figure 2c,d are shown for both the closest water molecule and for the cases when a strict water bridge is formed. For all nine bridge sites, the maxima of the distance distributions for the two different states (closest water or strict bridge) were within 0.1 Å of each other. At all of the sites, the mean distance between the water oxygens and amide hydrogen was larger than the mean distance between the closest hydrogen and the carboxyl oxygen. In the strict bridge case, the mean distance between the water oxygen and the amino hydrogen ranged from 2.02 to 2.11 Å across the nine bridge sites. The mean distance between the water hydrogen participating in the strict bridge and the carboxyl oxygen ranged from 1.90 to 1.95 Å across the nine bridge sites. Detailed data on the distance distributions between water molecules and the collagen bridge atoms can be found in the Supporting Information. Occupancy and Residence Times. Table 1 shows the occupancy and residence times across the nine bridge sites for the five different bond states. The distribution of residence times was asymmetric with a maximum at shorter times and a long tail at larger times. In an effort to characterize the residence time effectively, Table 1 records both the numberweighted mean residence time and the time-weighted median residence time for all states. Nonstrict Hydrogen Bonding Criteria. The residence times of the closest water molecule differed significantly between the bridge sites: generally, the sites near the edge of the imino-poor region (i.e., bridges 1−3 and 7−9) had mean and timeweighted median residence times larger than those at the center of the imino-poor region (i.e., bridges 4−6). Distributions of the time-weighed residence times of the closest water molecule for the nine bridge sites are shown in Figure 3. Although a particular water molecule may remain the closest water molecule to a bridge site for upward of 50 ps, the closest hydrogen atom, within the closest water molecule, to the carboxyl oxygen changed more frequently. The distributions of time-weighted residence times for the closest hydrogen are also

distance from a water hydrogen to the carboxyl oxygen bridge atom and was calculated for both hydrogen atoms in each water molecule; the smallest value was used in eq 2. After identification of the closest water molecule, the closest hydrogen atom was identified as the hydrogen atom within the closest water molecule with the smallest dCO···H. For the other three hydrogen-bonding states, a cutoff distance of 2.4 Å18 was applied to the coordinates of the closest water molecule to identify strict hydrogen bonds: (3) single strict bond between the closest water oxygen and the collagen amino hydrogen, dN−H···O ≤ 2.4; (4) single strict bond between the closest water hydrogen and the collagen carboxyl oxygen, dCO···H ≤ 2.4; and (5) a strict water bridge, dN−H···O ≤ 2.4 and dCO···H ≤ 2.4. Note that states (2)−(5) are not obtained by examination of each water molecule in the simulation but only for the water molecule that is in state (1), that is, the water that exhibits the lowest value of parameter d2 in eq 2. Whereas state (5) will include all instances of a strict water bridge formation, states (3) and (4) may not include all of the single hydrogen bonds between the bridge atoms and water molecules. It is possible that other water molecules form hydrogen bonds with the collagen bridge atoms if the closest water molecule to the bridge site does not form a strict bridge.



RESULTS Comparison of the Harmonic Constraints. No direct correlation between the value of the harmonic constraint constant used in the simulation and any of the physical parameters of interest (residence times, occupancy, etc.) was found. We used the results from the simulations performed with the constraint constant k = 1 kcal mol−1 Å−2 for the analysis presented here. The data for each of the four harmonic constraints used (k = 0.1, 1, 5, 10 kcal mol−1 Å−2) can be found in the Supporting Information. Distance to Bridge Atoms. The distance, throughout the 10 ns simulation, between the oxygens in the closest water molecules and the amide hydrogen for bridge site 1 is shown in D

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

Article

The Journal of Physical Chemistry B

distributions of the residence times, weighted by time, are shown in Figure 5a. The residence time of the strict hydrogen

Figure 3. Residence times, weighted by time, of the closest water molecule (gray) to the bridge sites. The residence times of a particular hydrogen atom, within the closest water molecule, as the closest hydrogen to the carboxyl oxygen are shown in blue. Boxes represent the 25−75th percentiles, whiskers indicate the minimum and maximum values, and the white line is the median value (recorded in Table 1).

shown in Figure 3. The change of the closest hydrogen atom to the carboxyl oxygen occurred as a “flip” between the current closest hydrogen in the water atom and the other; an example is shown in Figure 4. The time scale of this transition was smaller Figure 5. (a) Distributions of the residence times, weighted by time, for the strict bond states: strict hydrogen bonds between the water hydrogen and the carboxyl oxygen (blue, left); strict hydrogen bonds between the water oxygen and the amino hydrogen (red, middle); and the strict bridge state (gray, right). The boxes represent the 25−75th percentiles, the whiskers indicate the minimum and maximum values, and the white line in each box plot represents the median value (recorded in Table 1). (b) Fraction of the total simulation time, occupancy, spent by the closest water molecules making a strict bridge (gray), forming a single hydrogen bond with the amino hydrogen only (red), and forming a single hydrogen bond with the carboxyl oxygen only (blue). The sum of the three occupancies is smaller than one because there are times when the closest water molecule forms no strict hydrogen bonds with either of the collagen bridge atoms.

Figure 4. Distance between the hydrogen atoms in a closest water molecule to bridge site 3 and collagen carboxyl oxygen; this was the longest time spent by a water molecule near a bridge site in the simulation. The distances for the different hydrogen atoms are shown in lighter and darker blue to emphasize the switching of the closest hydrogen to the carboxyl oxygen. The dotted line is at 2.4 Å, the cutoff distance used for strict hydrogen bonding.

bond between a water hydrogen and the carboxyl oxygen was larger than the residence time of the bond between the water oxygen and the amide hydrogen group for each of the bridges 1−7, but not for bridges 8 and 9. Figure 5b compares the strict bridge occupancy to the fraction of the simulation time spent in a strict hydrogen bond with either the carboxyl oxygen or the amino hydrogen atom (but not both). The occupancies for each strict bridge state ranged from 57 to 74% across the nine bridge sites (Table 1). One or more strict hydrogen bond between the closest water molecule and the collagen bridge atoms was observed at least 96% of the simulation time for each of the bridge sites. Water Molecule Energetics. To investigate the mechanism responsible for constraining the water molecules at the bridge sites, the potential energy of the water molecules at the sites was calculated. For each of the nine sites, the ten water molecules that spent the longest time as the closest water molecule to the site were identified, and the energy of these water molecules during their time at the bridge site was calculated. Figure 6a shows the distribution of energies of these ten water molecules across the nine different sites. The bulk water energy distribution was obtained by averaging over 50 water molecules from a 40 Å TIP3P water box simulation, 900

than the 0.1 ps output of the MD simulation, and a detailed analysis of the motion was unavailable to us. The mean and time-weighted median residence times of the closest hydrogen to the carboxyl oxygen did not follow the same trend across the bridge sites as the residence time of the closest water molecule. This is particularly obvious in the case of bridge sites 8 and 9, which have the longest time-weighted median residence times for the closest water to the bridge site (∼140 ps) but have the closest hydrogen residence times similar to those of bridges 4 and 5 (∼20 ps). Strict Hydrogen Bonding Criteria. The mean residence time for the strict water bridges was subpicosecond for each of the nine bridge sites and did not have the same dependence on bridge site, as observed in the closest water molecule residence times. As expected, the lifetimes of the strict single hydrogen bonds were longer than the residence times for the water bridges, but the former were still shorter than 3 ps. The E

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

Article

The Journal of Physical Chemistry B

the bridge vector was taken as the vector that connects the water oxygen to the hydrogen that was closest to the CO group. The mobile proton O−H vector connects the water oxygen to the other water hydrogen (which does not directly participate in the bridge). The distributions of angle θ made by the bridge O−H vector with the mean bridge O−H vector are shown in Figure 7 (left column) for four example bridge sites. These distributions were similar between the bridge sites. The maximum angle (taken as the 98th percentile) made by a bridge O−H vector while in the strict bridge state ranged from θ = 35.3 to 37.9° across the nine bridge sites. Without the strict bridge requirement, this angle ranged from θ = 67.9 to 80.2° across the nine sites. Conversely, the rotation of the mobile proton O−H vector, around the mean bridge axis, varied significantly depending on the bridge site. The distributions of rotational angle ϕ of the mobile proton O−H vector around the mean bridge O−H vector are also shown in Figure 7 (middle column). To compare the rotations between sites, the data are presented such that the vector connecting the z axis (collagen axis) and the mean position of the water oxygen molecule at each site is at ϕ = 0. It is clear from these distributions that the mobile proton O−H vector is not able to fully rotate about the mean bridge axis but spends the majority of time oriented away from the collagen molecule. The correlation function of sin2(ϕ) was calculated for each separate bridge trajectory longer than 1 ps (10 time points) and for trajectories where the same hydrogen remained as the closest hydrogen to the CO group for longer than 1 ps. The individual correlation functions for each trajectory varied dramatically for a given site. The average of the individual trajectory correlation functions Cϕ(τ) was calculated for each bridge site, and examples when no strict bridge criteria are imposed are shown in Figure 7 (right column). Note that the error bars in the right column of Figure 7 are standard errors, not standard deviations. The correlation functions plateau at a value of ∼0.5 or higher, indicating strongly restricted rotational motion. Allowed H−H Vector Rotations. The intramolecular H−H vector of a water molecule near a bridge site was defined as the vector from the hydrogen atom closest to the carboxyl oxygen atom to the other water hydrogen atom. The mean angle of the intramolecular H−H vector with the axis of the collagen fiber, the z axis, is shown in Table 2 for the nine different sites. The average of all H−H vectors across the total simulation time is the time-averaged proton−proton alignment vector ⟨H−H⟩ (Figure 8). The angle between ⟨H−H⟩ and the z axis, or collagen axis, was 15.9° when only H−H vectors from water molecules making a strict water bridge were considered and 13.2° when no strict criteria were imposed. Examples of the intramolecular H−H vector orientations throughout the simulation are shown in Figure 9 for the same bridge sites as those in Figure 7. The mobile proton O−H rotations and resulting intramolecular H−H vector rotations of bridge groups 1, 5−9 were consistent across all four harmonic constraints used (Supporting Information). The rotational motion of the water molecules at bridge sites 2−4 did change across the four harmonic constraint values used in the simulations. In addition to the two bridge bonds, water molecules at bridge sites 2−4 are capable of making a third bond. The water hydrogen that is not directly participating in the water bridge can become hydrogen bonded to the oxygen in the O−H side chain of the THR amino acid residue; the resulting H−H vector rotations are shown in Figure 10 for the

Figure 6. (a) Distributions of the potential energy of water molecules in a bulk water simulation (blue; labeled w) and the energy of water molecules at each of the nine bridge sites (gray). The distributions of the energy across the nine bridges were calculated from the ten water molecules that spent the longest time (as the closest water molecule) to each bridge site. The boxes represent the 25th to 75th percentiles, the whiskers indicate the minimum and maximum values, and the white line in each box plot represents the median value (b) The relationship between residence time and mean energy of a water molecule (during its time as the closest water to one of the nine bridge sites) for the ten water molecules with the longest trajectories at each bridge site. Results for other constraint values are shown in the Supporting Information.

ps in length (other parameters as presented in the Method section). The mean energy value of water molecules at the nine bridge sites ranged from 0.98 to 1.17 of the mean bulk water energy value (−9.63 kcal mol−1). The mean energy of bridge sites in the middle of the imino-poor region (sites 3−6) ranged from 12 to 13% lower than the bulk water value. The mean energy of water molecules at the edges of the imino-poor region (bridge sites 1, 2 and 7−9) deviated from the bulk water value by 1− 8%. Despite this, the bridge sites at the edge of the imino-poor region tend to confine water molecules for longer (Table 1). The relationship between the mean energy of the 10 water molecules during the time at the bridge site and the length of the trajectory at the bridge site is shown in Figure 6b. Water molecules with lower mean energies spent less time at the bridge sites than those water molecules with mean energies closer to that of bulk water. Restricted Rotational Motion. Bridge and Mobile O−H Vectors. The rotational motion of the water molecules during their time at the bridge site was investigated by separating the two O−H vectors into the bridge O−H vector and the mobile proton O−H vector, which were defined as follows. In the case of a strict bridge, the bridge O−H vector was the vector that connects the water oxygen to the hydrogen that was bonded to the CO group. In the case where no strict bridge is formed, F

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

Article

The Journal of Physical Chemistry B

Figure 7. Rotational motion of the O−H vectors in the water molecules associated with bridge sites 1, 3, 7, and 8 (top row to bottom row). The data are rotated such that the vector connecting the z axis (collagen axis) and the mean position of the water oxygen molecule at each site is at ϕ = 0. The left column shows the angle θ made by the bridge O−H vector with the mean bridge O−H vector. The middle column shows the angle ϕ of the mobile O−H vector around the mean bridge O−H vector axis. The solid histograms represent the relevant distribution when no strict hydrogen bond criterion is imposed, and the dark line is the distribution when only strict bridges are considered. Insets in the left column show the distribution of O−H bridge (orange) and O−H mobile (green) vectors in the strict water bridge state on the surface of a sphere. The water oxygen is positioned at the center of the sphere, the black arrows point in the direction of the mean O−H bridge vector, and the black line on the surface of the sphere indicates ϕ = 0. The right column shows the average correlation function Cϕ(τ) of the mobile O−H vector around the mean bridge O−H vector when no strict hydrogen bond criteria are imposed. Only trajectories where a particular hydrogen in the closest water remained the closest hydrogen to the carboxyl oxygen for longer than 1 ps were included in the calculation of Cϕ(τ). Error bars are standard errors, not standard deviations.

Table 2. Angle between the z Axis (Collagen Axis) and the Average Intramolecular H−H Vectors at each of the Nine Bridge Sites throughout the 10 ns Simulationa bridge closest hydrogen

1

2

3

4

5

6

7

8

9

⟨H−H⟩

70.4 70.7

46.0 48.2

46.9 49.8

33.9 34.3

61.6 61.4

63.3 61.8

62.1 61.1

66.1 64.4

66.4 65.7

15.9 13.2

The last column refers to the angle between the z axis and the average H−H vector across all nine bridge sites ⟨H−H⟩ (see Figure 8 for illustration).

a

1BKV collagen peptide. Fullerton proposes that these sites, along with other water bridges, are responsible for the MA effect in aligned collagenous tissues due to a combination of fast exchange between the bridge and bulk water protons and restricted rotational motion.16 Although MD simulations do not allow a direct simulation of the proton exchange between bridge and bulk water molecules, the simulations presented here enable us to investigate two key features of Fullerton’s model: the icelike nature of bridge bonds and the rotational motion of the water molecules near the bridge sites. In this work, we followed the distance-only geometric definition of a hydrogen bond used by other researchers.18 This differs from the normally used definition that is based on both the distance and the angular criteria. The less strict, distance-

strict bond case. As a result of this bond, the motion of the water molecule at bridge sites 2−4 was more restricted compared to the motion of water at the other bridge sites, and the mean intramolecular H−H vector was more parallel with the axis of the collagen molecule (Table 2).



DISCUSSION We have previously proposed a schematic molecular model of the MA effect that explains the orientational dependence of the T2 values of water protons in ordered collagenous tissues.35 In the present study, we investigated the molecular origins of this effect at the atomistic level. We used MD simulations to examine the detailed dynamics of the water molecules in singlewater bridge sites in a hydrated, individual molecule of the G

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

Article

The Journal of Physical Chemistry B

Figure 8. Average H−H vectors when no strict bond criteria are imposed for the nine different bridge sites; coloring indicates bridge number as in Figure 1. Left: average H−H at their positions along the collagen molecule; collagen residues participating in the bridges are colored according to bridge. Right: average H−H vectors positioned at the origin; gray circles indicate the labeled angle (in degrees) from the z axis as a visual guide. The black vector is the mean H−H vector across all bridge sites ⟨H−H⟩.

only, definition was used in this study because the positional constraints placed on the collagen backbone have the ability to limit the accessible hydrogen bond configurations. The use of the distance-only definition alleviated this artificial constraint; it also enabled us to directly compare our results with the hydrogen bond lifetimes obtained when the collagen backbone was not constrained.18 Although the lifetimes corresponding to different definitions of a hydrogen bond may differ, this difference can be expected to be sufficiently small not to affect the “icelike” (or non-icelike) nature of the water bridges in the collagen molecule. Previous simulations by Ravikumar et al.18 investigated lifetimes of water bridges in the 1BKV collagen peptide; they used a cutoff distance of 2.4 Å to define hydrogen bonds and saved molecular coordinates every 0.5 ps. They reported mean lifetimes between 1 and 2.5 ps, with occupancies ranging from 10 to 70% across the nine bridge sites at 300 K. Our results for strict water bridge mean residence times at 310 K are comparable, given the frequency of our output was 0.1 ps, which potentially lends itself to lower mean lifetimes but larger recorded occupancies. However, our results also demonstrate that the residence time of the closest water molecule to a bridge site was much larger than the residence time of the strict water bridge (Table 1) with continuous trajectories longer than 100 ps at some bridge sites. Rather than a quick succession of many different water molecules (every 1−2 ps), this suggests that one particular water molecule remains close to a particular bridge site (for tens or hundreds of picoseconds), constantly breaking

Figure 9. Examples of the distributions of the angles between the allowed H−H vectors and the collagen axis for water molecules in the bridge position (solid blue histograms) and in a strict water bridge state (solid black lines): (a) bridge site 1, (b) bridge site 3, (c) bridge site 7, and (d) bridge site 8. Insets for each panel show the distribution of the H−H vectors in the strict water bridge state on the surface of a sphere; the black arrows indicate the direction of the collagen axis (and by definition the z axis).

and reforming water bridges according to the strict definition of the bridge. The mean energy of water molecules at the bridge sites was typically lower by between 0.12 and 1.65 kcal mol−1 than the mean bulk water energy (−9.36 kcal mol−1), indicating a weak energetic confinement. However, bridge sites where this energy difference was the greatest (in the middle of the iminopoor region) also had the lowest residence times for water molecules. This suggests that confinement of the water molecules at the bridge sites for these longer periods of time (upward of 100 ps) was predominately due to steric H

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

Article

The Journal of Physical Chemistry B

reported in their study, simulation times were 500 ps in length and they report several distinct water molecules at the bridge sites during these times. This suggests that the lifetime of the water molecule at the bridge site was of the order of hundreds of picoseconds. A similar lifetime (∼100 ps) was observed by Fu et al.20 for dangling hydrogen bonds between the water oxygen and N−H groups on a collagen peptide simulated in a larger crystal structure. Fu et al. also observed a water bridge in their simulations and stated that the lifetimes were of the order of a few nanoseconds. However, it is unclear if this is the lifetime of individual water molecules making the water bridge or the lifetime of the availability of the bridge site to water molecules; the bridge site did collapse to a direct hydrogen bond between CO and N−H groups at some time point in their simulation. Although the residence times are longer when the strict hydrogen bond criteria are relaxed, the water bridge hydrogen bonds in the MD simulations cannot be considered icelike. This presents an apparent discrepancy with experimental evidence, which does suggest that the bridge bonds are icelike. Fullerton and Cameron16 used DSC to record the enthalpy of denaturation of rat tendon at different hydrations. At hydrations of less than four water molecules per amino acid triplet (comprising one single-water bridge and one doublewater bridge16), the calculated enthalpy required to remove one water molecule was 1.045 kcal mol−1. According to Fullerton, this is consistent with the enthalpy of fusion of pure water ice (1.436 kcal mol−1), given that the four water molecules in the water bridges form a total of six icelike bonds, not the eight formed in bulk water ice.16 Similarly, Boryskina et al.39 used IR spectroscopy to calculate the average enthalpy of hydration of a water molecule as 3.9 kJ mol−1 (0.93 kcal mol−1) in the combined Langmuir and Henry sorption layers (comprising ∼3.2 water molecules per amino acid triplet) in a collagen molecule. A recent MRI study by Tadimalla and Momot40 also seems to agree with the ice bridge hypothesis. The authors studied the effect of H2O−2H2O replacement to monitor the effect on transverse spin relaxation rates (R2 = 1/T2) in articular cartilage. The gyromagnetic ratio γ of deuterium (2H) is ∼7 times weaker than that of a proton nucleus. If the anisotropic contribution to R2, R2A, is the result of 1H−1H interactions before deuterium replacement and 2H−1H interactions after deuterium replacement, then R2A should decrease with increasing 2H2O concentrations. Instead, the authors found R2A to be independent of the 2H2O concentration. Fullerton notes that if the bridge protons were icelike, they would not exchange with the deuterium population and there would be no change to the dipolar fields experienced by the exchanging proton population.17 However, the 2H2O replacement in Tadimalla and Momot’s study took place over 12−24 h. For the protons participating in the water bridges to remain unavailable for exchange, the hydrogen bond lifetime would have to be of the order of the 2H2O replacement time. Theoretically, the MA effect could still be achieved through stochastic phase shifts of exchangeable protons, according to Fullerton’s model, without icelike bridges or strict hydrogen bonds. For this to be possible, two requirements must be met: (1) the restricted rotation of the intramolecular H−H vector in water molecules associated with the bridge sites such that the time-averaged vector across all bridge sites ⟨H−H⟩ is coaxial with the collagen molecule; and (2) fast exchange between protons in water molecules associated with the bridge sites and

Figure 10. Distribution of the angles between the H−H vector and the collagen axis for water molecules that make a strict bridge at bridge site 3. The solid gray histograms show the angle distribution when the mobile proton is (a) not hydrogen bonded to the THR residue (>2.4 Å from the oxygen of the O−H group in the THR residue) and (b) hydrogen bonded to the THR residue (