Statistical Errors


A & T pgs. 191-198

Today I am going to review a basic idea that is central to the course: how to estimate the errors of a simulation.

First let me define a simulation. A simulation has some "state" variables, which we denote by S. In classical mechanics these are the positions and velocities of the particles. In an Ising model they are the spins of the particles. We start the simulation at some initial state S0. We have some process which modifies the state to produce a new state so that we can write Sn+1=T(Sn). We will call the iteration index "time". It is in quotes because it may not refer to real time but a fictitious time or pseudo-time, sometimes called Monte Carlo time.

The iteration function T(S) can either be deterministic (like Newton's equations of motion called Molecular Dynamics) or have a random element to it (called Monte Carlo). In the random case, we say that the new state is sampled from the distribution T(S). What we discuss in this lecture and in much of the course, is the situation where the "dynamics" is ergodic, which we define roughly to mean that after a certain time one loses memory of the initial state, S0, except possibly for certain conserved quantities like energy, momentum etc. This implies there is a correlation or renewal time. For time differences much greater than the correlation time, all properties of the system not explictly conserved are unpredictable as if they were randomly sampled from some distribution.  Ergodicity is often easy to prove for the random transition but usually difficult for the deterministic simulation. We will discuss ergodicity further in the next lecture.

Assuming the system is ergodic one can discuss the state of the system as a probability distribution: Fn(S). (it is the distribution of outcomes that you get if you take a distribution of initial states.) At large time this will approach a unique distribution, the equilibrium state F*(S). In classical statistical mechanics this is the Boltzmann distribution.

The main goal of most simulations is the calculation of a set of properties in equilibrium. The most common example is the internal energy E = <e(S)> where the average means over F*(S). In a simulation we sample e(Sn) for n<T where T is the total length of the run. The curve of e versus T we call the trace. Calculating an estimate of E is easy: just take an average over the run. The only "trick" is to be sure to throw away the part of the trace that too much influenced by the initial conditions. This is called the transient, equilibration portion or warm-up. (Question: if you have a big spike at the beginning of your trace and you neglect to throw it out of the average, how long does it take before the influence of the spike is less than the statistical error? That is, how does the estimated mean depend on the propertiess of the warm up portion?)

The next order of business is to figure how accurate is the estimate of the exact value: namely the estimated error in the estimated mean commonly called the error bar. Simulation results without error bars are only suggestive. Without error bars one has no idea of its significance.  All respectable papers in the literature should  contain estimates of any properties that are important to the conclusion of the paper. All homework exercises must include errors estimates.  Many mistakes are made in estimating errors. You should fully understand both the equations and get an intuitive feel for how to estimate errors. (“eye ball” method) Although we have prepared a code to automatically estimate errors, you should know the formulas and what they mean. On the exam you will have to do it without DataSpork. Often one hears the complaint that one cannot estimate errors because the quantity being measure is a very complicated function of the state. However, it is always possible to get a rough idea of the error bar which is often all that is needed,  by simply rerunning the entire code a few times with different initial conditions and looking at the spread of results. This is called the "eyeball" method. While there are many ways for errors to be underestimated, the eyeball method will not seriously overestimate the errors.

The fundamental theorem on which error estimates is based is the central limit theorem due to Gauss. It says that the average of any statistical processes will be "normally" or have a Gaussian distribution. Further the error bar will equal square root of the  variance of the distribution divided by the number of uncorrelated steps. Later we will discuss in detail the conditions under which this result holds. The correlation time of sequence is the number of steps it takes for a fluctuation to disappear.

There is a terminology which we use to describe the process of estimating the <A> and estimating the error of our estimate. You need to learn what all these terms mean so you can communicate with your fellow scientists about statistical distributions. We have developed a JAVA Analysis Tool,  DataSpork, which does the computations for you discussed here. In the first homework you have to learn how to use it and what the results mean.

The symbols and formulas are given on the pdf formula page.

    1. Mean of A (a). The estimate of the mean value is simply the average over the equilibrated data.
    2. Variance of A ( v ). The mean squared deviation of A(t) about the estimated mean.

See also the detailed page of pdf formulas.


Calendar

©1998-2001 David Ceperley and D.D. Johnson