Dipolar induced spin-lattice relaxation in the myelin sheath: A molecular dynamics study

Interactions between hydrogen protons of water molecules and macromolecules within the myelin sheath surrounding the axons are a major factor influencing the magnetic resonance (MR) contrast in white matter (WM) regions. In past decades, several studies have investigated the underlying effects and reported a wide range of R1 rates for the myelin associated compartments at different field strengths. However, it was also shown that the experimental quantification of the compartment-specific R1 rates is associated with large uncertainties. The current study therefore investigates the longitudinal relaxation rates within the myelin sheath using a molecular dynamic (MD) simulation. For this purpose, a realistic molecular model of the myelin sheath was employed to determine the dipole-dipole induced R1 relaxation rate of the hydrogen protons at clinically relevant field strengths. The results obtained clearly reflect the spatial heterogeneity of R1 with a increased relaxivity of myelin water due to a reduced molecular mobility near the membrane surface. Moreover, the calculated R1 rates for both myelin water and macromolecules are in excellent agreement with experimental findings from the literature at different field strengths.

from the membrane surface and could thus affect the water dynamics over the whole myelin water space 12 . Since the longitudinal relaxation is direct related to the molecular mobility, it is thought that these restrictions strongly affect the R 1 rate of the myelin water associated molecules 13 .
However, several studies in the past years have shown that an accurate isolation of the MR signal of myelin associated water and a reproducible quantification of its R 1 rate is quite complicated and associated with uncertainties [13][14][15][16][17][18][19][20][21][22][23][24][25] . In addition to the limited accuracy, many of these experiments were performed at different field strengths, which complicates a direct comparison of their findings. On the other hand, in the past decades classical molecular dynamic (MD) simulations have become a reliable and well established tool for modeling the physical properties of molecules and reconstructing experimental data. Moreover, the continuous improvement of molecular force field simulations and the variety of parameterized lipid and protein molecules allows to simulate complex biological structures, i.e. a myelin-alike environment on atomistic scale.
The aim of the current work was therefore the theoretical investigation of the dipole-dipole induced R 1 relaxation rate in the myelin sheath using an MD simulation. For this purpose, we first carried out an MD simulation of a lipid-water system representing a small fraction of the membrane structure inside the myelin sheath (the lack of protein molecules in this model is discussed in the Limitations section). The hydrogen trajectories obtained from this atomistic simulation were employed in the next step to determine the dipolar induced relaxation rate of the 1 H protons inside the myelin sheath at field strengths of 1.5 T, 3 T and 7 T. Based on these results, the average R 1 rates of the lipid membrane and the water compartment were calculated and compared with experimental findings at different field strength. In addition, the cross-relaxation rate between the water and the macromolecules was determined to assess its contribution to the total magnetization transfer.

Methods
Dipolar relaxation in a two-compartment system. The formalism for the calculation of nuclear spin relaxation employed in this study is based on the well-established approach introduced by Bloembergen, Purcell and Pound 26 . In this approach the authors derived the relaxation rates in liquid substances for spin-½ particles by treating the motion of the molecules classically whereas the dipolar interaction Hamiltonian was described quantum mechanically. It is assumed that secondary relaxation mechanisms such as j-coupling, spin-rotation or chemical shift anisotropy of the molecules are negligible. Under these conditions, the relevant interaction between two equivalent spins I and S is given by the dipole-dipole Hamiltonian which can be written in terms of second-rank tensors as where T 2,q are spin operators and Fq are second degree spherical harmonics of order q. The polar coordinates r, θ and ϕ describe the position of spin S relative to I, whereas the main magnetic field, B 0 , defines the z-direction of the coordinate system. Spin operators and spherical harmonics satisfy the relationships = = − , respectively, where the asterisk indicates the complex conjugate, and the latter are given by: A full definition of the spin operators can be found in 27 . For the case of molecules that are subject to thermal motion, the relative position and thus the Hamiltonian becomes time dependent. Fluctuations with spectral components around the Lamor frequency, ω γ = B L 0 , serve as energy source for the longitudinal relaxation process. The relevant frequencies for the spin transitions can be obtained from the real part of the corresponding spectral densities which are given by ( ) where G q is the time autocorrelation function of the q-th spherical harmonic and τ is the lag time. The angular brackets on the right-hand side denote the average over time. For isotropic tumbling particles such as water molecules in the liquid phase, the autocorrelation function decays exponentially on timescales of a few picoseconds and the corresponding spectral density can be described by a simple Lorentzian curve. In macromolecules, however, hydrogen nuclei are bonded in a flexible lattice resulting in an overall motion that is governed by the internal vibrations of the nuclei and the global rotation of the whole molecule. Since both processes occur simultaneous and independently on different time scales, the corresponding correlation functions can be factored and written in a normalized form as where M τ is the global correlation time of the macromolecule, i τ the internal vibrational correlation time, S 2 the order parameter describing the restriction of the motion of individual protons and can take values from 0 to 1 28 . The respective spectral density can then be written as linear combination as where  q is the magnitude of the correlation function given by F F . Following the semiclassical approach from Bloembergen et al., the transition probability per unit time between two two-particle spin states q q 1 2 and ′ ′ q q 1 2 , W, is given to first order by where  is the reduced Planck constant, q = 1,2 is the order of the spherical harmonics and γ I(S) is the gyromagnetic ratio of spin I(S). Using these definitions and the relationship between the transition probabilities and the longitudinal relaxation rate, which can be derived from the kinetic equations of the state population, the R 1 rate for two hydrogen nuclei (I S 1/2 = = ) can be expressed in terms of the spectral densities as 27 Here, μ 0 is the vacuum permeability and γ H is the gyromagnetic ratio of the hydrogen nuclei. Eq. (7) holds true for two nuclei in magnetically equivalent environments and can easily be extended for a many spin system. In order to apply Eq. (7) to the myelin sheath, it must be taken into account that the dynamics of hydrogen nuclei within rigid macromolecules is different from those within the mobile water molecules. The motion-restricted nuclei are subject to static magnetic fields from neighbouring spins and ions which make them insensitive to exciting pulses near the Lamor frequency 29,30 . Such an excitation pulses, for example, would depose energy in the water molecules but does not affect the net spin population within the macromolecules. Due to their different absorption properties, it follows that the relaxation processes of hydrogen nuclei assigned to the water-pool (WP) and hydrogen nuclei assigned to the macromolecule-pool (MP) have to be determined separately. For this purpose, Eq. (7) has to be extended to account for the different environments: Here, j R 1,A is the longitudinal relaxation rate of the j-th proton and N A is the total number of hydrogen protons in pool A (WP or MP), i. e. the hydrogen nuclei are exclusively assigned to either WP or MP. In addition to this compartment-intrinsic induced relaxation, it must be taken into account that the hydrogen nuclei still interact across pools, thereby enhancing their transition probability mutually. Moreover, the unbalanced spin populations due to spectral broadening causes a MT between the two pools. These inter-compartmental interactions, which exclusively involve 1 H protons from the respective other pool, are referred as auto-and cross-relaxation 31 and are given by where N B is the total number of 1 H protons assigned to the pool B. Figure 1 shows how the different relaxation processes can be classified in context of the pool affiliation. In order to determine the average auto-and cross-relaxation rates, the Eqs (7-10) have to be averaged over all hydrogen nuclei within the respective pools. For cross-relaxation, this procedure leads to the expression , which can be interpreted as the transfer of magnetization per time unit from pool B to A. Since the latter expression corresponds to a double sum over N A and N B weighted by the pool size, the average cross-relaxation rates satisfy which holds also true for the average auto-relaxation rates. To compare the calculated relaxation rates with experimental findings, it has to be considered that the influence of cross-relaxation cannot be distinguished from the influence of chemical exchange in conventional MR experiments. Although both mechanisms gover the transfer of magnetization between different pools, they do not change the sum magnetization of the combined system and can thus be decoupled from the auto-and intrinsic relaxation rates. The effective relaxation rates, R ef 1 , for the each of the pool, without the contribution of MT effects, can then be defined as www.nature.com/scientificreports www.nature.com/scientificreports/ It should be mentioned that Eqs (7-10) are expressed in a general form without any assumptions regarding the molecular structure or motion. In the literature it is often assumed that the molecular motion is unrestricted and isotropic. For such a stochastic process, the spherical harmonics yield an exponentially decaying spectral density and the relaxation rates can be expressed in terms of the reorientational correlation time of the water molecules. Since the water molecules are motionally restricted and thus the molecular tumbling is not necessarily isotropic, the general form was used here.

System preparation and molecular dynamics simulation.
To simulate the dynamics of hydrogen atoms inside the myelin sheath, a MD simulation was performed using the GROMACS package 5.1.2 and the CHARMM36 all-atom force field 32,33 . The membrane structure of the myelin sheath was represented by two lipid-leaflets and was modeled according to a realistic myelin alike molecular composition using the CHARMM-GUI web interface 34  This composition corresponds to a phosphor-and galactolipid content of 46% and 27%, respectively, in accordance with histological observations 10 . The space between the leaflets was solvated with 4699 water molecules, which corresponds to a water mass-ratio of 40%. In addition, 8 sodium molecules were added to the water compartment to neutralize the net charge of the sulfatides. The water molecules between the leaflets were represented by the TIP4p-FB model 35 . This water model yields significantly better results for water specific characteristics compared to the standard mTIP3p model and is, as far as we can judge, compatible with the CHARMM36 force field (see supporting information).
Before the production runs, during which the trajectories were saved for the further analysis, the system was thermodynamically equilibrated for a MD simulation time of 50 ns to body temperature (310.15 K) and normal pressure according the following settings: The equations of motion were solved using the leap-frog algorithm with a simulation time step of 2 fs. The cut-off lengths of the Lennard-Jones potential and short-range Coulomb interaction were set to 12 nm using the Verlet-scheme. Long-range electrostatic interactions were calculated using the particle-mesh Ewald summation with a grid spacing of 0.5 Å. During all simulation steps, the temperature and pressure were kept constant and controlled by the Nosé-Hoover and Parrinello-Rahmen coupling with a time constant of τ = 1 T ps and 2 p τ = ps, respectively. The length (but not the angles) of the hydrogen bonds in the membranes was constrained using the LINCS algorithm.
Based on the equilibrated system, the production runs were conducted using the same settings as described above. The run was performed for T MD = 20 ns and a saving frequency of ΔT = 0.5 ps. During the production run, the lipid layers were in the L β -phase (gel phase) and the tail distance fluctuates around 8.2 nm, which is in accordance with the typical periodicity of the myelin sheath 36,37 . After the run, the trajectories were recomposed again to avoid broken molecule structures due to the periodic boundary condition using the subroutines provided by the GROMACS package. Finally, the trajectories of the hydrogen nuclei of the WP and MP were separately extracted. The final configuration after the trajectory conversion and the density profiles of the corresponding membrane components are shown in Fig. 2a,b.
It should be mentioned that the production run was also performed as an NVE (constant particle number, volume and energy) instead of an NPT (constant particle number, pressure, and temperature) ensemble in order to avoid additional molecular fluctuations due to the temperature and pressure couplings. However, both ensemble configurations yielded, within the standard error, similar relaxation rates for the water molecules. In addition, it should also be mentioned that employing the original water model (mTIP3p) of the CHARMM36 force field instead of the TIP4p-FB model would lead to quite different relaxation rates for the water pool (see supporting information). www.nature.com/scientificreports www.nature.com/scientificreports/ Computation of the relaxation rates. Based on the extracted trajectory sets, the compartment-intrinsic relaxation rates of the WP protons were calculated as follows: First, a single hydrogen nucleus was randomly selected and the relative positions to the remaining hydrogen nuclei (N WP = 9397) within the WP were computed for each time step and converted into polar coordinates. Then, from the polar coordinates, the spherical harmonics given by Eq. (2b,c) and the corresponding correlation functions G 1 (τ) and G 2 (τ) were calculated in the interval τ = 0, …, 0.3T MD using fast Fourier transformation. The truncation at τ = 0.3T MD was performed to take into account that the tail of the correlation function is dominated by noise due to the finite trajectory length. Next, the spectral densities γ J B ( ) , of the selected hydrogen nucleus was determined using Eq. (8).
Since the water molecules diffuse several times across the water compartment within T MD = 20 ns, the calculation above yields representative results for the average relaxation rate of the WP but is not well suited for studying the specific influence of the reduced molecular mobility near the membrane surface. For this purpose, the procedure describe above was repeated for ′ = .
T 1 5 MD ns, whereby the starting points for the trajectory analysis were randomly selected from the original data set. The time interval of this short-time analysis was chosen so that the average WP-R 1 rate is comparable with the results where the whole trajectory set was analysed, but the averaged positions of the individual hydrogen nuclei still reflect the location within the WP. In addition to the R 1 rates, the diffusion coefficients of the water molecules were calculated by fitting the Stokes-Einstein relation to the mean square displacement of the oxygen trajectories. This short-term analysis was performed for 2000 randomly selected trajectory parts and yields the diffusion profile shown in Fig. 2c.   Figure 2. Snapshot of the lipid-water system after the production run and the trajectory conversation (a). Density profile of different lipid species inside the membrane as well as water molecules (b). Self-Diffusion profile of water molecules across the water pool (c). The dotted lines indicate the Gibbs-dividing surface and the error bars the 25 th and 75 th percentile.
The determination of the relaxation rates of the MP associated hydrogen nuclei (N MP = 13810) was performed similarly, except for the determination of the spectral densities and the short-time analysis. First, the relative positions between a randomly selected hydrogen nucleus and the remaining nuclei were computed and converted into polar coordinates. To determine the spatial distribution of the R 1 rates within the lipid membrane, the average positions of the selected hydrogen nuclei were stored. Then, the spherical harmonics and the correlation functions were computed using Eqs (2)  were computed using Eqs (6)(7)(8). Due to the long computation time of the fitting procedure, the relaxation rates of the MP-protons were calculated for 1381 random selected hydrogen nuclei, which corresponds to 10% of the total compartment size. These hydrogen nuclei were also used to calculate the auto-and cross-relaxation in the MP. For this purpose, the relative positions between one proton and all hydrogen nuclei from the WP were calculated to determine the correlation functions. Since the correlation function describing the inter-compartmental interactions completely converged within the simulation time, the spectral densities H 2 0 could be directly computed using explicit integration. Next, from the spectral densities, the auto-and cross-relaxation rates were calculated according Eqs (9) and (10). Then, the computed rates were averaged in order to determine auto-and cross-relaxation rates of the WP using Eq. (11). Finally, the effective relaxation rates,  Table 1 contains the averaged relaxation rates calculated for WP and MP at 1.5 T, 3 T and 7 T. As is evident from the results, all rates increase with decreasing field strength. For R WP 1,

Results
, the average relaxation rates decrease from 1.175 Hz at 1.5 T by a factor of 1.6 to 0.712 Hz at 7 T. A similar decrease by a factor of 2.0 and 1.7 is also observed for R MP 1, and R MP auto 1, , respectively. In contrast, for the cross-relaxation rates, a decrease from 0.061 Hz and 0.042 Hz by factor of 3.8 to 0.016 and 0.011 Hz is observed for the WP and MP, respectively. The self-diffusion coefficient of water obtained from the short-term analysis is shown in Fig. 2c. As can be seen in this Figure, diffusion of the water molecules near the membrane surface is significant reduced and reaches a maximum of approximately 2 10 5 ⋅ − cm 2 /s in the centre of the WP. Figure 3 shows the different relaxation rates of the WP and MP as a function of the hydrogen positions and the different field strengths. As it is evident from this Figure, the relaxation rates of the macromolecules increases from R 1, MP ≈ 1Hz at the end of the hydrocarbon tails to approximate 5 Hz for 1.5 T, and 4 Hz for 7 T, in the head groups. For the relaxation rates in the WP shown in Fig. 3b, a parabolic relaxation profile is observed with a minimum of R 1, MP ≈ 0.5 Hz at the middle of the WP and a signification increase of R MP , respectively. Here, the auto-and cross-relaxation rates increase sharply from a few mHz at the tail positions by a factor ≈100 at the head groups, which can be directly attributed to the − r 6 dependency of the dipolar interaction between the WP and MP. The fitting procedure of the effective relaxation rates according Eq.

Discussion
The aim of this work was to determine the relaxation rates inside the myelin sheath and its surrounding water compartment using an MD simulation. As is evident from Table 1 and Fig. 3, the calculated relaxation rates clearly reflect the expected R 1 dispersion in WM structures 38,39 . Moreover, the results obtained demonstrate that the relaxation rate of the myelin water is strongly affected by motional restrictions near the membrane surface. As can be derived from the Stokes-Einstein relation and the relaxation theory, the diffusivity of spherical particles such as water molecules is inverse proportional to the rotational correlation time and thus to the R 1 rate. This inverse relationship becomes also evident when one compares the spatial distribution of R 1 rates in water (Fig. 3a) and the spatial profile of the self-diffusion coefficient (Fig. 2c). The influence of the reduced molecular mobility on the R 1 rate of water is evident over the entire WP and increases the average R 1 rate in the WP significantly in comparison to water in the liquid phase (1.1 Hz vs 0.3 Hz). In the molecular pool a distinct increase of the relaxation rates along the tail structure was observed at all three field strengths. This behaviour can be attributed to the decreased nuclei mobility in the hydrocarbon chains for lipids in the L β -phase 40 and illustrates a spatial inhomogeneity of R 1 within the myelin sheath. www.nature.com/scientificreports www.nature.com/scientificreports/ Relaxation rate. The results of the current work fit well with experimental findings at different field strength 13,[16][17][18][19][20][21][22][23][24][25] . In 17-20 , phospholipid membranes were extracted from lecithin probes and their R 1 rates were examined at different temperatures and field strengths between 0.2 T and 8.45 T. For body temperature, the authors reported R 1 rates between 12 Hz at 0.2 T and 2 Hz at 8.45 T. The function fitted to the effective relaxation rate of the MP as determined theoretically in the current study yields R 1 rates of 12.4 Hz and 2.3 Hz, respectively, in agreement with experimental findings. In contrast to the measurement on membrane model systems, an experimental quantification of R 1 within the myelin WP and thus a comparison with the findings of the current study is more complicated. For field strength of 1.5 T, for example, the reported R 1 rates range between 1.8 HZ and 2.7 HZ (see Table 1 in ref. 41 ). These findings are slightly larger than the calculated R 1 rates of 1.34 Hz at the same field strength in this study. However, it has to be considered that some of these studies employed a two-component (free and bound water) analysis and reported the relaxation rate of the fast decaying component whose relative contribution is between 8-20% as the R 1 of myelin water. This method, however, is not sensitive to the MR invisible macromolecules and as such provides a R 1 rate which contains a signal contribution from the MP due to magnetization transfer. An accurate quantification of R 1 , therefore, includes a proper distinction of the respective molecular compartments , which requires a complex multi-component analysis of the tissue samples. Even if this technique is well understood, the underlying model has several degrees of freedoms or covarying parameters and is therefore associated with systematic uncertainties even in sophisticated experiments. For example, in ref. 21 the authors combined MT experiments with a multi-echo gradient echo measurement to determine the residence time of water molecules in the myelin sheath at 7 T. Even though this study provided important insights in the WM microstructure, it fails to separate the relaxation rates of two pools and reported an intermediated R 1 of 1.84 Hz for the myelin sheath. Given the consistent results of the calculated and experimental R 1 rates for the MP, the reported R 1 rate can be compared with the average relaxation rate obtained in this work. The latter is given by   www.nature.com/scientificreports www.nature.com/scientificreports/ associated compartment at the different field strength are shown in Fig. 4. The R 1 rates which could be clearly assigned to the WP, MP or to the weighted average R ( ) , are illustrated in the corresponding colours. In contrast, the R 1 rates, which were assigned to the WP but might still contain a signal contribution from the MP were plotted as black circles. As is evident from this Figure, the experimental findings of the MP agree well with the calculated R 1 rates. In addition, the experimental findings for the myelin water are constrained by the effective relaxation rates of the WP and MP. This means that an appropriate weighting of the respective R 1 rates as well as the inclusion of magnetization transfer between the pools should allow to reconstruct the experimental data.

Cross-relaxation and magnetization transfer.
Apart from the compartment-specific relaxation rates, the dipolar-induced cross-relaxation rates between WP and MP were calculated. As already mentioned in the introduction, it is still unclear how strong this mechanism contributes to the total magnetization transfer. In saturation-transfer experiments, for example, the magnetization transfer between lecithin/cholesterol and isotopically substituted protons was investigated 42 . The authors concluded that the signal differences between 3 H protons and 1 H protons in the surrounding water compartment could be explained by a dominance of magnetically induced exchange over chemical exchange. In contrast, other experiments in tissue and model systems have shown significant pH effects, suggesting chemical exchange 43,44 . If one compares the calculated cross-relaxation rates with results obtained from the measurement of magnetization transfer, this conclusion is also supported by the current study. The total exchange rate, k WP , from the WP to the MP, for example, can be written as WP WP cross c e 1,

= +
where R ce is the contribution due to chemical exchange. Studies on bovine WM have shown that k WP is in the order of 12 Hz which is approximately two orders of magnitude higher than the few hundred mHz calculated in this study 16 . This result indicates that the dipole-dipole induced cross-relaxation is negligible compared to chemical exchange and thus plays a minor role for the total magnetization transfer. Unfortunately, and in contrast to the physically well defined cross-relaxation, chemical exchange is characterised by complex chemical (and kinetic) processes on nanosecond-millisecond time scales and requires a comprehensive analysis of the bound enthalpy and the electron configurations of the reactants. Since these parameters are not extractable from a classical MD simulation, a more detailed investigation of the chemical exchange contribution is outside the scope of the current work.
Relating the calculated R 1 rates within the myelin sheath to the typical exchange rate of approx. 12 Hz as determined in MT experiments reveals that it might be difficult to extract individual compartment-intrinsic R 1 rates. The magnetisation transfer is fast compared to the relaxation rates in SM and MW. Therefore, the two pools are in rapid exchange so that hydrogen nuclei within the myelin appear as a single magnetisation pool on the longitudinal relaxation time scale.
Limitations. The current study is subject to some limitations that may have an impact on the interpretation of the reported findings. First, it should be noted that the R 1 rates were determined for a perpendicular orientation of the membrane normal relative to the main magnetic field (see Fig. 1). As NMR studies have shown, however, the relaxation rate of macromolecules is affected by the membrane orientation relative to B 0 45,46 . This orientation dependency should also hold true for the R 1 relaxation rates of the MP as reported in Table 1. Second, it should be mentioned that the membrane model employed here did not consider the protein content inside the myelin  16,22,24 were performed at room temperature between 20° C to 24° C. The remaining data were either in vivo experiments or were performed at 39° C (ref. 19 ) and 45° C (ref. 18  www.nature.com/scientificreports www.nature.com/scientificreports/ sheath. Previous studies have shown that the protein fraction makes up between 19% to 27% of the dry mass 47,48 of the myelin sheath and consists mainly of myelin basic (∼30%) and proteolipid proteins (∼50%) 49 . Several studies have shown that some of these proteins are located only in the cytoplasmic space or are surrounded by a unique lipid composition 50,51 . The inclusion of proteins employing a realistic localization within the cellular environment and molecular distribution would therefore significantly increase the simulation run-time. However, the absence of proteins might not signifcantly alter the calculated R 1 rates as can be roughly estimated from the results presented in 18 . Here, the authors investigated the effects of proteins on the R 1 rate of choline lipids. It was demonstrated that the relaxation rate increases by about 38% for a lipid to protein ratio of 114, which approximately corresponds to the ratio within the myelin sheath. It was suggested that the higher R 1 rate can be explained by the restriction of the lipid motion due to the heavier protein molecules. However, in the current study calcium pump proteins (∼110 kD) were added to the choline layers, which are approx. five times heavier than the myelin basic or proteolipid proteins (∼18 kD to 30 kD). It can therefore be assumed that the lighter proteins within the myelin sheath might have less influence on the lipid dynamics and thus on the R 1 of the membrane sheath. This assumption is also supported by a recently published ex vivo study, which concluded that the R 1 relaxation in WM regions is mainly governed by the lipid and not the protein content 52 . Therefore, it can be presumed that a more realistic molecular system including proteins would yield a slightly higher R 1 rate for the solid myelin but should not significantly affect the final interpretations of this work. Finally, the current study only investigated the effect of the 1 H induced dipole-dipole relaxation. Although this should be the most effective relaxation mechanism in biological tissue, it can be assumed that other spin = 1/2 particles such as the phosphorus nuclei and secondary effects such as J-coupling, susceptibility anisotropy or the iron content additionally affect the R 1 rate.

conclusion
In this work we investigated the dipolar-induced R 1 relaxation inside the myelin sheath at several clinical relevant field strengths of 1.5 T, 3 T and 7 T. For this purpose, we employed an atomistic model of the myelin sheath using a realistic lipid composition and a suitable water model. The calculated results are in excellent agreement with experimental findings over a wide range of clinically relevant field strengths. Given the consistent results, the presented model might serve as a basis to investigate further MR properties of the myelin sheath on a microscopic scale. Furthermore, the results obtained suggest that the dipolar induced cross-relaxation plays a minor role for MT and illustrates the increased R 1 rate of water molecules due to the reduced molecular mobility. Especially the accurate determination of the compartment specific relaxation rates shown in Table 1 might be helpful in future studies to assess experimental findings or to reconstruct the MR signal evolution in myelin rich structures.

Data Availability
The simulation files to generate the MD data are available from the corresponding author on request.