[gmx-users] g_energy averages
Erik Lindahl
lindahl at stanford.edu
Wed Jul 17 20:42:53 CEST 2002
David L. Bostick wrote:
> .. So, if I want to find the standard error in an average value, I just
> take the fluctuation and divide it by sqrt(nsteps)?
>
> i.e.
> If I get from the energy file
> Energy Average RMSD Fluct.
> #Surf*SurfTen 147.805 2270.69 2263.84
>
> and the run was over 250000 steps, the standard error is
>
> 2263.84/sqrt(250000)
>
> even if nstenergy was say 500 steps?
>
Yes and no.
It is only correct if each point is an independent observation. If we
have a membrane where the area fluctuates on a scale of 10ns the area
fluctuations will be very small if we only simulate it for 10ps, even
if we save 100000 points during those 10ps. We can't say anything about
the 'physical' (in contrast to purely statistical) standard error
though, since we don't even cover a single period of the motion.
To get an accurate estimate of the standard error you will have to take
the autocorrelation time into account. Load the curve into xmgrace,
calculate the autocorrelation of the fluctuations, and fit e.g. a single
exponential A*exp(-t/tau) to it - 'tau' is your correlation time.
Assuming the total length of your simulation is T, we can then estimate
the number of independent observations as T/tau, so the standard error
is
<rmsd fluctiations>
s=-------------------
T/tau
In the limit where the observations are completely uncorrelated the
correlation time is the stepsize, so in that case the denominator would
be the total number of steps.
The alternative is to split the trajectory into parts that you *know*
are longer than the correlation time, calculate an average from each of
these, and then the standard error as the fluctuations of these averages
divided by the number of parts. This will only be an upper estimate of
the standard error though, and you will still need an estimate of the
correlation time.
Cheers,
Erik
More information about the gromacs.org_gmx-users
mailing list