Equilibrium Free Energies from Nonequilibrium MeasurementsUsing Maximum-Likelihood MethodsMichael R. ShirtsDepartment of Chemistry, Stanford University, Stanford, California 94305-5080, USAEric Bair and Giles HookerDepartment of Statistics, Stanford University, Stanford, California 94305-4065, USAVijay S. PandeDepartment of Chemistry, Stanford University, Stanford, California 94305-5080, USA(Received 27 June 2003; published 2 October 2003)We present a maximum likelihood argument for the Bennett acceptance ratio method, and derive asimple formula for the variance of free energy estimates generated using this method. This derivation ofthe acceptance ratio method, using a form of logistic regression, a common statistical technique, allowsus to shed additional light on the underlying physical and statistical properties of the method. Forexample, we demonstrate that the acceptance ratio method yields the lowest variance for any estimatorof the free energy which is unbiased in the limit of large numbers of measurements.DOI: 10.1103/PhysRevLett.91.140601 PACS numbers: 05.10.–a, 05.70.CeIntroduction.—Finding the free energy difference be-tween different states of a physical system is of greatgeneral interest in many scientific fields, including drugdesign [1], basic statistics [2], and even nonperturbativequantum chromodynamics [3]. It is of interest to theexperimental community as well as the theoretical andcomputational communities [4]. Recently, there has beenincreased interest in determining the uncertainty and biasin any attempt to extract free energies from a suitable setof data [2,5–11].We can separate the calculation of precise and accuratefree energy differences into two nonoverlapping prob-lems. First, we must generate a number, n, of statisticallyuncorrelated measurements of the system. Second, wemust extract a free energy estimate from these n mea-surements, along with reliable estimates for the statisticalbias and variances of our estimate.We will assume that weare already in possession of a set of n uncorrelated mea-surements of the proper observable for our method, andaddress only the statistical issues related to the extractionof free energy estimates from these measurements.There are a variety of commonly used methods forfinding the free energy of a physical change in the system.Many of these can be expressed finding the equilibriumfree energy from nonequilibrium work distributions.Thermodynamic perturbation theory (TPT) or free en-ergy perturbation (FEP) estimates free energy differencesby exponentially averaging potential energy differencesbetween a reference state sampled at equilibrium and atarget state [12]. However, FEP can be seen as a specialcase of ‘‘fast growth’’ [5,6] nonequilibrium exponentialwork averaging, as the energy difference is the infinitelyfast adiabatic work of transition between the two states[13,14]. ‘‘Slow growth’’ thermodynamic integration hasbeen shown to have high intrinsic biases and is unreliableas originally implemented [15,16]. However, the ‘‘freeenergies’’ obtained from these simulations are actuallymeasurements from a nonequilibrium work distributionand an ensemble of values can therefore also be used toobtain correct free energies [6]. The rest of this Letterwill focus on the generalized nonequilibrium problem.Assume there are two different equilibrium states de-fined on a phase space by energy functions U0~qq andU1~qq.Let F be the free energy between these states,defined as the log of the ratio of the partition functionsassociated with U0~qq and U1~qq. We can associate a workwith the process of changing energy functions from U0toU1or vice versa while the system is maintained in tem-perature equilibrium with the surroundings. By samplinginitial conditions from equilibrium, we obtain a distribu-tion in either direction of such work values. For infinitelyfast switching, these distributions are simply of U  U1 U0 canonically sampled from the initial state.It has long been known that the exponential average ofequilibrium energy differences between two states yieldsthe free energy difference between the states [12]. Morerecently, Jarzynski demonstrated that distribution of non-equilibrium work values can yield an equilibrium freeenergy by taking the exponential average of the set ofnonequilibrium work values [13]. However, the exponen-tial average of a set of data X fxi; ...;xng, defined as1= lnhexpXi (where   1=kT), is a statisticthat is both inherently noisy and biased, even if the spreadof the data is only moderately larger than kT.Theresults of exponential averaging strongly depend on thebehavior at the tails of the distribution, which, by defini-tion, are not as well sampled as the rest of the distribution.Previous studies have explored and demonstrated thePHYSICAL REVIEW LETTERSweek ending3 OCTOBER 2003VOLUME 91, NUMBER 14140601-1 0031-9007=03=91(14)=140601(4)$20.00  2003 The American Physical Society 140601-1poor behavior of exponential averaging for small samplesizes [7–10].In an examination of free energy estimation betweentwo states sampled at equilibrium, Bennett [17] demon-strated that it is possible to use the information con-tained in both the forward and reverse distributions of thepotential energy difference together in a manner whichwas significantly better than estimates using only energydifference data in one direction. This derivation cantrivially be generalized to the nonequilibrium workcase, replacing U with the nonequilibrium work [14].Bennett, in the FEP case, and Crooks, in the general case,showed that the equation:exp FhfWiFhfW expWiR(1)is true for any function fW, where we define the mea-surement from the initial state to the final state as the‘‘forward’’ direction (denoted by the subscript F)andthe measurement from the final state to the initial stateas the ‘‘reverse’’ direction (denoted by the subscript R).Bennett then minimized the statistical variance withrespect to this function fW to find that fW1 nF=nRexpW  F1minimizes the variance in thisfree energy estimate, where nFand nRare the number ofsimulations in the forward and reverse direction, respec-tively. The free energy difference F can easily be foundby iterative methods [17]. Although this method is knownand referenced in the literature, it