0README atom-fluctuations |-- 0NOTE obvious |-- 0README obvious |-- yout184-dget.Data trajectory data: chi-2 of BPTI Tyr-35 |-- yout184.s R commands for analysis of trajectory |-- yout184-out.s same, but to generate ps and txt files |-- yout184-out.ps plots generated by analsysis `-- yout184-out.txt text generated by analsysis TRAJECTORY AND MODEL FOR ANALYSIS A molecular dynamics simulation for a protein at equilibrium generates a trajectory for a set of variables sufficient to describe the motions of the molecule. These stochastic motions, driven by the ambient thermal energy, are fluctuations about the equilibrium values. The set of variables can be the atom coordinates, from which at each time point one can calculate trajectories for particular bond lengths, bond angles, and dihedral angles. We consider fluctuation of the dihedral angle chi-2 of tyrosine 35 of bovine pancreatic trypsin inhibitor (BPTI). These motions correpond to rotations of the tyrosine ring about an axis approximately coincident with the C-beta--C-gamma bond. Tyr-35 is largely buried in the protein matrix. We can view it as moving in a potential well the shape of which is determined by interactions of the ring atoms with atoms of residues adjacent in tbe protein structure. The thermal motions of the adjacent atoms provide both the random force driving the ring motions and the viscous force opposing its motion (relaxation). This picture suggests that the fluctuations of chi-2 of Tyr-35 may be modeled as those of an harmonically bound Brownian particle. The harmonic motion, described by a torsional force constant, or equivalently, a frequency of oscillation, is determined by the shape of the potential well. The Brownian motion of the ring is determined by thermal motions of the protein matrix and is described by a friction constant, in effect modeling the protein matrix as a solvent in which the ring is embedded. The trajectory of chi-2 of Tyr-35 of BPTI analyzed below is from a simulation at 100K (please do not ask, why not 300K?). Details of the simulation are in file "0NOTE". Part of the analysis is fit of the data with a tethered Brownian oscillator model. The quality of the fit is adequate. It is perhaps remarkable that a simple 2-parameter model can capture the essential features of Brownian motion of a ring buried in the protein matrix. R LANGUAGE COMPUTATIONS AND RESULTS The file "yout184.s" has a set of commands to analyze a trajectory for the dihedral angle chi-2 of Tyr-35 of BPTI: read in the chi-2 trajectory as a vector of angular fluctuations (degrees from the mean value) for each time point, saved in an S/R dump format (dput/dget); construct an equal-length vector of the time of simulation for each time point (picoseconds from start of run); two plots of the trajectory: full trajectory, 16384 data points, 32.768 ps; first 500 points, 1 ps; Fourier transform the signal and normalize; two plots, of modulus = Mod(fft(yout184)/16484) full transform, 8192 frequency bins, F-hat and F-hat*; (2nd half of plot is redundant) first 750 points of the transform, F-hat; construct the time correlation function: C(t) = Re(fft(Conj(fft(yout184)/16484)*fft(yout184)/16484,inverse=TRUE)) fit the tethered Brownian oscillator model to the data for C(t) (first 500 points), using the nonlinear least squares algorithm of nls(); the data are from a simulation for 100 K: mean chi-2 = -3.845215e-05 variance chi-2 = 19.49047 std.dev. chi-2 = 4.414802 Parvseval's theorem - size signal vector invariant to change basis var(yout184) = sum(yout184^2)/16384 = sum(Mod(a)^2)) = C(0) = 19.48929 first 20 values of C(t); note C(0) = variance [1] 19.489285 18.009616 15.247437 13.633616 13.706762 14.740587 16.006862 [8] 16.270025 14.439963 11.713520 10.733692 12.100012 13.886611 14.502295 [15] 13.789323 11.949852 9.803482 9.191256 10.793810 12.592394 summary of fit to C(t) Formula: Ct ~ C * exp(-time/tau) * (cos(sqrt(XX/C - 1/tau^2) * time) + 1/(tau * sqrt(XX/C - 1/tau^2)) * sin(sqrt(XX/C - 1/tau^2) * time)) Parameters: Estimate Std. Error t value Pr(>|t|) tau 0.185765 0.006422 28.92 <2e-16 *** C 12.316130 0.184796 66.65 <2e-16 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 Residual standard error: 1.664 on 498 degrees of freedom Correlation of Parameter Estimates: tau C 0.2501 torsional spring constant, from fit value C(0): K = 52.77597 Kcal/mol freq. of torsional oscillation, from K: w0 = 16.93435 ps^-1 freq. of Brownian model oscillation, from C & tau: w1 = 16.05596 ps^-1 ANALYSIS Plots 1 and 2 of file "yout184-out.ps" and the averages listed first in file "yout184-out.txt" (see above values) describe the trajectory of chi-2 of BPTI, from a simulation at 100K. The tyrosine ring undergoes Brownian motion with an RMS deviation of 4.4 degrees. Inspection of plot 2 indicates the relaxation time for a large fluctuation is between .1 and .2 ps. The Fourier transform of the signal is shown in plots 3 and 4. The strongest contributions are at the frequencies: 130, 75, 30, 16, and 7 ps^-1 (corresponding to 4300, 2500, 1000, 525, and 230 cm^-1). Below 4 ps^-1, the transform shows white noise (all frequencies with the same amplitude), consistent with Brownian motion of the ring, driven by very fast random collisions of the ring with its environment (the time scale of the collisional events is much shorter and well separated from that of the slow Brownian diffusion, so that ring motions can be modeled as driven by a delta-function random force). The time correlation function C(t) is obtained directly from the Fourier transform, as the inverse transform of the convolute: C(t) = F^-1(F(chi-2(t))* * F(chi-2(t))) where F() is the Fourier transform operator, F()* its complex conjugate, and F^-1() the inverse. The tethered Brownian oscillator model was fit to C(t): C(t) = C(0) * exp( -t/tau) * ( cos( w1 * t ) + 1/( tau * w1 ) * sin( w1 * t ) where w1 = sqrt(K/I - 1/tau^2), K=kT/C(0), and I is the moment of inertia of the ring. The two free parameters are tau, a relaxation time, and C(0), the zero-time correlation function. Plot 5 shows C(t) vs time, with the dashed line calculated from the best-fit parameters (tau=0.186, C(0)=12.3). The fit is good up to t ~ 0.6 ps; for longer times the motion is less strongly damped than predicted by the model. The fit value of C(0) is one third less than the value calculated from the trajectory. We may view the model as describing the principal features of the ring motion, while recognizing that the motion has additional components. As noted, it is remarkable that a simple two-parameter model can account for the essential features of so complex a system as a tyrosine ring moving within an asymmetric and dynamic protein matrix. The relaxation time tau is consistent with the decay of large fluctuations seen in the chi-2 trajectory (plot 2). From the estimate for the parameter C(0) of the fit (12.3 degrees^2), we obtain the torsional force constant for oscillation of the tyrosine ring: K = RT/1000*1/C(0) = 53 Kcal/mol which corresponds to a frequency of oscillation (calculated for a moment of inertia I = 7.7 10^15 gm-cm^2), w0 = sqrt(K/I) = 16.9 ps^-1 This oscillator frequency is close to the Fourier component of 16 ps^-1. The frequency, w1 = 16.1 ps^-1, of the periodic term of the tethered oscillator model is approximately equal to w0, i.e., 1/tau^2 << w0. The torsional force constant, K, corresponds to a soft motion, of much lower frequency than for bond stretching or bond bending, reflecting the deformability of the protein matrix and the weakness of non-covalent interactions. John Rupley rupley@u.arizona.edu -or- jar@rupley.com 30 Calle Belleza, Tucson AZ 85716 - (520) 325-4533; fax - (520) 325-4991 Dept. Biochemistry & Molecular Biophysics, Univ. Arizona, Tucson AZ 85721