Changes in protein dynamical parameters derived from quasi-elastic neutron scattering spectra induced by coherent scattering: A molecular dynamics simulation study

Protein dynamics at the sub-nanosecond timescale and the Å length-scale has widely been studied using quasi-elastic neutron scattering (QENS). In almost all QENS studies on hydrogenated proteins in D2O buffers, analysis of the spectra after buffer subtraction is carried out under the assumption that the remaining spectra arise from incoherent scattering of proteins while the contribution of coherent scattering is negligible. On the contrary, a study using polarization analysis has shown that the coherent scattering accounts for more than 10% of the total scattering intensity of hydrogenated proteins (Gaspar et al., Biochim. Biophys. Acta 1804:76–82 (2010)). In addition, the effects of coherent scattering on the values of dynamical parameters of proteins obtained by analysis of QENS spectra remain unclear. Here, molecular dynamics (MD) simulation on hen egg white lysozyme was used to investigate this issue. QENS spectra containing only incoherent scattering and those containing both incoherent and coherent scattering were calculated from the MD trajectory. Dynamical parameters were then extracted from the two simulated QENS spectra. Comparison of the resultant dynamical parameters has shown that the error in the values of the dynamical parameters induced by coherent scattering is at most 6%. This error is unlikely to significantly affect the results of QENS studies that investigate the relative changes in protein dynamics caused by different physicochemical conditions such as temperature unless dynamical parameters need to be determined with high precision at the absolute scale.


Introduction
Quasi-elastic neutron scattering (QENS) is an experimental technique that has extensively been used to study dynamical behavior of proteins at the sub-nanosecond timescale and the Å lengthscale [1][2][3][4][5][6].The incoherent neutron scattering cross-section of a hydrogen (H) atom is more than 40 times larger than that of its isotope, deuterium, and any other atoms found in proteins [7], and about half of the total number of atoms constituting protein molecules are H atoms.As a result, experimentally measured QENS spectra are dominated by scattering of H atoms of the protein molecules.Because H atoms are distributed almost uniformly in space, dynamical information studied by QENS reflects the motions of H atoms averaged over all the H atoms in protein molecules.
In typical QENS studies of proteins, solution samples of hydrogenated proteins in D2O buffers are used as specimens.The QENS spectra are represented by the dynamic structure factor S(Q " "⃑ , ω), where Q " "⃑ and ω are the momentum transfer and the energy transfer in the unit of ℏ of neutrons upon scattering, respectively.The spectra arising from protein molecules and their hydration water are obtained by subtracting the spectra of the buffer from those of the protein solution.The resultant spectra Sexp(Q " "⃑ , ω) give almost flat Sexp(Q " "⃑ ), where Sexp(Q " "⃑ , ω) is integrated with regards to ω, and this is expected for incoherent scattering [8,9].Since water molecules show a prominent broad peak arising from their coherent scattering at Q = 1.5-2.0[Å -1 ] [8,9], the lack of such a peak in Sexp( Q " "⃑ ) indicates that the contribution of coherent scattering is negligible in the system of .CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for this this version posted November 14, 2022.; https://doi.org/10.1101/2022.10.27.514007 doi: bioRxiv preprint hydrogenated proteins in D2O, and thus Sexp(Q " "⃑ , ω) contains the dynamical information of protein molecules.Analysis of Sexp(Q " "⃑ , ω) is, therefore, generally carried out with the assumption that it is equal to the incoherent dynamic structure factor Sinc(Q " "⃑ , ω) of protein molecules, for which various theoretical models and phenomenological fitting equations are formulated [1].On the other hand, Gaspar et al. have experimentally shown by polarization analysis, which is a technique that can separate coherent and incoherent scattering [10], that there exists more than 10% of Q-dependent coherent contribution of the hydrogenated protein itself over a wide Q range in Sexp(Q " "⃑ ) [11].
Although incoherent scattering is still the major component of Sexp(Q " "⃑ ), this seminal work has raised a concern about the validity of the dynamical information derived from the QENS study without polarization analysis.Therefore, it is of critical importance to investigate to what extent coherent scattering might affect dynamical parameters of protein molecules extracted from Sexp( Q " "⃑ , ω) compared with those extracted from Sinc(Q " "⃑ , ω).

Molecular dynamics simulation
In this study, one of the trajectories of a hen egg white lysozyme (HEWL) molecule (PDB ID: 1AKI) generated previously was used [12]: The software GROMACS (version 2018.3)[13] was used to perform molecular dynamics (MD) simulation with the AMBER03 force field [14].The heavy water model TIP3P-HW [15] was employed as a solvent molecule.In a cubic box, a HEWL molecule was placed and solvated such that all edges of the simulation box were ≥ 10 Å away from the protein surface, which simulates a protein molecule in solution.Counter ions were added to neutralize the charges of the entire system.van der Waals force and electrostatic interactions were cut off at 8.0 Å [14].
All the exchangeable (labile) hydrogen atoms of HEWL were replaced with deuterium by doubling the mass of the corresponding hydrogen atoms as in previous studies [16][17][18].Whereas this simple procedure cannot explain an experimentally observed molecular behavior occurring at a long time scale (400 ns) [19], significant errors were not seen up to 50 ns from the start of the simulation.As a simulation time, 6 ns was chosen such that it was short enough not to generate large errors in protein dynamics and at the same time long enough to extract dynamical parameters obtained by QENS.
The system was first energy-minimized by the steepest descent method, followed by NVT and NPT equilibrations, which were carried out for 100 ps each with a time step of 2 fs.The resultant structure was subjected to a production run for 6 ns at 300 K and at 1 atm, and the first 5 ns of the trajectory was discarded to extract part of the trajectory which have reached a stable mean square displacement (MSD).Translational and rotational motions of the whole molecule were removed from the MD trajectory using the GROMACS command "gmx trjconv" so that only intramolecular motions are observed.The trajectory was saved every 2 ps.

Calculation of incoherent and coherent dynamic structure factors from the MD trajectory
The incoherent dynamic structure factor Sinc(Q " "⃑ , ω) and the coherent dynamic structure factor Scoh(Q " "⃑ , ω) were calculated from the MD trajectory through the incoherent intermediate scattering function Iinc(Q " "⃑ , t) and the coherent intermediate scattering function Icoh(Q " "⃑ , t), respectively, as follows: where N is the total number of atoms constituting a protein molecule, and b !,!"# and b !,#*+ are the incoherent and coherent scattering lengths of the i-th atom, respectively, which are tabulated in Table 1.r & ""⃑(t) is the position vector of the i-th atom at time t, which is given by a MD trajectory, and 〈 〉 denotes the thermal average.In the case of solution samples, scattering intensity is orientationally averaged since protein molecules are randomly oriented, and the orientational average and the thermal average are independent.The vectors Q " "⃑ and r & ""⃑(t) − r & */ - """""""""⃑(0) are thus treated as their modulus Q = @Q " "⃑ @ and @r & ""⃑(t) − r & */ - """""""""⃑(0)@, respectively.Eqs.1a and 2a then read and Table where I /E0 (Q, t) is an instrumental resolution function, which defines the time-window of observation.In this study, this function was represented by a Gaussian function as, where  IJK is the resolution time, which was set to 10 ps in this study.This corresponds to the energy resolution (FWHM) of 77 μeV.Since I /E0 (Q, t) drops to almost zero at around t = 50 ps, As stated above, Sexp(Q, ω) is the sum of Sinc(Q ω) and Scoh(Q, ω) of protein molecules, i.e., Sexp(Q, ω) = Sinc(Q, ω) + Scoh(Q, ω).Relative contributions of incoherent and coherent scattering to Sexp(Q, ω) were investigated through the scattering intensities Sinc(Q) and Scoh(Q), where Sinc(Q ω) and Scoh(Q, ω) are integrated in the energy direction: S #*+ (Q) corresponds to a scattering curve, which can also be obtained by small-angle neutron scattering.Because Sinc(Q ω) and Scoh(Q, ω) were simulated in the energy range between -0.5 meV and 0.5 meV, the integrations in Eqs.6 and 7 were also conducted in this energy range.

Analysis of the dynamic structure factors
Dynamical parameters such as the residence time and the jump diffusion coefficient, that are typically obtained by QENS measurements, were extracted by fitting the simulated Sexp(Q, ω) and Sinc(Q, ω) with the following phenomenological equation: where C(Q) is a scaling factor between the fit and the simulated spectra, and S N*# (Q, ω) represents .CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for this this version posted November 14, 2022.; https://doi.org/10.1101/2022.10.27.514007 doi: bioRxiv preprint 9 the contribution from local atomic motions in a protein molecule.R(Q, ω) and B(Q) are a resolution function, which is a Fourier transform of Eq. 5, and a background, respectively.L N*# (Q, ω) is a Lorentzian function describing local atomic motions, which is described as where Γ N*# (Q) is the half width at the half maximum (HWHM) of the Lorentzian.The term including δ(ω) is the elastic component and that including provides the detailed features of the corresponding diffusive motions.Γ N*# (Q), which was obtained by the fitting of the simulated QENS spectra using Eq. 8, was found to be well described by the jump-diffusion model [20]: where τ N*# and D N*# are the residence time and the jump-diffusion coefficient of local atomic motions, respectively.
Q-dependence of A(Q), the elastic incoherent structure factor (EISF), provides information on the geometry of local atomic motions.The EISF was fitted by the following equation that describes a diffusion-inside-a-sphere model [21]: . CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for this this version posted November 14, 2022.; https://doi.org/10.1101/2022.10.27.514007 doi: bioRxiv preprint 10 where p < is the immobile fraction of atoms, motions of which are too slow to be resolved at the energy resolution employed, while (1 − p < ) denotes a fraction of atoms resolved at the energy resolution employed.The latter atoms undergo diffusive motions within a sphere of a radius of a. j ) () is the first-order spherical Bessel function of the first kind.

Results and Discussion
First, Sinc(Q) and Scoh(Q) of HEWL calculated from the MD trajectory were compared to investigate the relative fractions of incoherent and coherent scattering in Sexp(Q), which corresponds to the experimentally obtained S(Q) of proteins (Sexp(Q) = Sinc(Q) + Scoh(Q)).The results are shown in Fig. 1 (a).It is seen that Sinc(Q) is independent of Q values, which is consistent with experimental studies [8,9,11] while Scoh(Q) shows a strong Q-dependence with small peaks at around Q = 0.6 and 1.4 Å - 1 .Sexp(Q) thus shows slightly larger intensity below Q = 1.0 Å -1 than above Q = 1.0 Å -1 , which is a tendency also observed in Sexp(Q) obtained from the experimental QENS spectra on another protein troponin [9].Scoh(Q) corresponds to the so-called small-angle neutron scattering (SANS) curve, which is used for structural analysis of proteins, and hence it was compared with a SANS curve of a crystal structure of HEWL calculated using the program CRYSON [22] (Fig. 1 (b)).The resultant SANS curve agrees quite well with the simulated Scoh(Q), suggesting the validity of the calculation procedure used in this study.Fig. 1 (c) shows the relative fractions of incoherent and coherent scattering in the simulated Sexp(Q).It is found that the fraction of incoherent scattering depends on Q values, as reported previously using polarization analysis [11], and it reaches 80-90% as Q .CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for this this version posted November 14, 2022.; https://doi.org/10.1101/2022.10.27.514007 doi: bioRxiv preprint increases.Accordingly, the fraction of coherent scattering is 40% at Q = 0.4 [Å -1 ], above which analysis of QENS spectra is generally carried out, and then decreases to 10-20% at higher Q values.
Next, the shapes of the simulated Sinc(Q, ω) and Scoh(Q, ω) were directly compared, which has not previously been reported, neither experimentally nor theoretically.For this purpose, Scoh(Q, ω) was normalized such that its peak intensity at ω = 0 [meV] is equal to that of Sinc(Q, ω).The results are shown in Fig. 2 for Q = 0.5, 1.0, and 1.5 [Å -1 ].It is seen that both spectra have essentially almost the same shape.On the other hand, when one takes a closer look at these spectra as shown in the insets of Fig. 2, the degree of broadening of the spectra is slightly different between Sinc(Q, ω) and Scoh(Q, ω) and this difference shows Q-dependence.To investigate the effects of these coherent scattering contributions on dynamical parameters derived from QENS measurements, analyses of the simulated Sexp(Q, ω) and Sinc(Q, ω) were carried out as a next step.
The simulated Sexp(Q, ω) and Sinc(Q, ω) were fitted using the phenomenological equation of Eq. 8.A fitting example is shown in Fig. 3, which shows that both spectra are reproduced in good quality by this equation.The resultant Q 2 -dependences of the widths of the Lorentzian (Γloc(Q)) and Qdependences of the EISFs are summarized in Fig. 4. The Γloc(Q) was found to be fitted by the jumpdiffusion model (Eq.10) in both cases as shown in Fig. 4 (a) and (c).Moreover, both EISFs were also reproduced by the diffusion-inside-a-sphere model (Eq.11) as shown in Fig. 4 (b) and (d Sexp(Q, ω) contains Scoh(Q, ω) in addition to Sinc(Q, ω), slight differences in the broadening of the QENS spectra between Sinc(Q, ω) and Scoh(Q, ω), which was shown in Fig. 2, resulted in almost identical values of the Γloc(Q) and the EISF as shown in Fig. 4 (e) and (f), respectively.Although the relative fraction of coherent scattering is 40% at Q = 0.4 [Å -1 ] (Fig. 1 (c)), this does not affect the Γloc(Q) and the EISF values (Fig. 4 (e) and (f)).This is because the shape of the simulated QENS spectra is very similar between Sinc(Q, ω) and Scoh(Q, ω).
Dynamical parameters obtained by the fitting of the simulated Sexp(Q, ω) and Sinc(Q, ω), i.e., the jump diffusion coefficient (Dloc), the residence time (τloc), the radius of a sphere (a), and the immobile fraction of atoms (p0), are compared in Fig. 5.It is found that the Dloc value increases from 1.54 ± 0.03 × 10 -5 cm 2 /s to 1.63 ± 0.03 × 10 -5 cm 2 /s by the contribution of coherent scattering, and τloc decreases from 11.67 ± 0.07 ps to 11.42 ± 0.07 ps.The a value increases from 1.389 ± 0.004 Å to 1.438 ± 0.003 Å while the p0 value increases from 0.722 ± 0.001 to 0.737 ± 0.001 by introducing coherent scattering.These results indicate that while the immobile fraction of atoms slightly increases, the contribution of coherent scattering to QENS spectra apparently leads to slightly enhanced protein dynamics.
Using modern neutron spectrometers such as BL02 DNA at J-PARC [23], dynamical parameters can now be estimated with a precision of 1-2% [9].On the other hand, the differences in dynamical .CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for this this version posted November 14, 2022.; https://doi.org/10.1101/2022.10.27.514007 doi: bioRxiv preprint parameters caused by coherent scattering discussed in the previous paragraph are 5.8%, 2.2%, 3.5%, and 2.1% for Dloc, τloc, a, and p0, respectively.Thus, differences in some dynamical parameters induced by coherent scattering are beyond the precision of experiments and hence it is essential to bear in mind that dynamical parameters extracted from QENS spectra measured without polarization analysis might be varied by several percent compared with those that would be obtained from purely incoherent QENS spectra, which do not contain the contribution of coherent scattering.The present study thus suggests that in the case where dynamical parameters need to be determined at the absolute scale with high precision, QENS measurements with polarization analysis will be required as already proposed before [11].On the other hand, in the case where changes in dynamical behavior of proteins induced by temperature [24,25], pressure [26,27], and other biological conditions [9,28], are investigated, the widths of the Lorentzian function and the EISF values generally show much larger changes than those shown in Fig. 4 (e) and (f) as the sample state is varied.In this case, a several percent of differences in dynamical parameters caused by coherent scattering would have negligible effects on the conclusion of QENS studies on proteins.
To the author's knowledge, this is the first study to investigate the possible effects of coherent scattering on the dynamical parameters obtained by QENS spectra on proteins.In summary, analysis of the simulated QENS spectra, Sexp(Q, ω) and Sinc(Q, ω), has shown that despite the nonnegligible contribution of coherent scattering, which accounts for 10-40% depending on the Q-

.Figure 5 .
Figure 5. Summary of dynamical parameters extracted from the simulated QENS spectra.
. CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in 23 44⃑ 2•2/ !44 4⃑(9);/ " 44 4⃑(<)2 〉 + .
. 1 Incoherent and coherent scattering lengths of atoms used in this study.The values are taken from the literature [1].
Sinc(Q, ω) and Scoh(Q, ω) are the time Fourier transforms of I !"# (Q, t) and I #*+ (Q, t), respectively, CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in ) and I #*+ (Q, t) were calculated up to t = 50 ps using the MD trajectory with the time length .
).Although .CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in . CC-BY 4.0 International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in dynamical parameters extracted from the simulated Sinc(Q, ω) are modulated by at most 6% by coherent scattering.However, these differences in dynamical parameter values are unlikely to affect the results of QENS studies to investigate the relative changes in protein dynamics caused by different physicochemical conditions of the sample unless dynamical parameters need to be determined with high precision at the absolute scale.International license perpetuity.It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in