f15-0

Sample page from NUMERICAL RECIPES IN FORTRAN 77: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43064-X) Copyright (C) 1986-1992 by Cambridge University Press.Programs Copyright (C) 1986-1992 by Numerical Recipes Software. Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin g of machine- readable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes booksor CDROMs, v isit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to [email protected] (outside North Amer ica).Chapter 15. Modeling of Data 15.0 Introduction Given a set of observations, one often wants to condense and summarize the databyfittingittoa“model”thatdependsonadjustableparameters. Sometimesthe model is simply a convenientclass of functions, such as polynomials or Gaussians, andthefitsuppliestheappropriatecoefficients. Othertimes,themodel’sparameterscome from some underlying theory that the data are supposed to satisfy; examples are coefficients of rate equations in a complex network of chemical reactions, or orbitalelementsofa binarystar. Modelingcanalsobeusedas akindofconstrained interpolation,whereyouwanttoextendafewdatapointsintoacontinuousfunction, but with some underlying idea of what that function should look like. The basic approach in all cases is usually the same: You choose or design a figure-of-merit function (“merit function,” for short) that measures the agreement between the data and the model with a particular choice of parameters. The meritfunction is conventionally arranged so that small values represent close agreement. The parameters of the model are then adjusted to achieve a minimum in the merit function, yielding best-fit parameters . The adjustment process is thus a problem in minimizationinmanydimensions. ThisoptimizationwasthesubjectofChapter10; however, there exist special, more efficient, methods that are specific to modeling,and we will discuss these in this chapter. Thereareimportantissuesthatgobeyondthemerefindingofbest-fitparameters. Data are generally not exact. They are subject to measurement errors (callednoise in the context of signal-processing). Thus, typical data never exactly fit the model that is being used, even when that model is correct. We need the means to assess whether or not the model is appropriate, that is, we need to test the goodness-of-fit against some useful statistical standard. Weusuallyalsoneedtoknowtheaccuracywithwhichparametersaredetermined by the data set. In other words, we need to know the likely errors of the best-fit parameters. Finally, it is not uncommon in fitting data to discover that the merit function is not unimodal, with a single minimum. In some cases, we may be interested in globalratherthanlocalquestions. Not,“howgoodis this fit?” but rather,“howsure am I that there is not a very much better fit in some corner of parameter space?” As we have seen in Chapter 10, especially §10.9, this kind of problem is generally quite difficult to solve. The important message we want to deliver is that fitting of parameters is not the end-all of parameter estimation. To be genuinely useful, a fitting procedure 650 15.1LeastSquaresas aMaximum LikelihoodEstimator 651Sample page from NUMERICAL RECIPES IN FORTRAN 77: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43064-X) Copyright (C) 1986-1992 by Cambridge University Press.Programs Copyright (C) 1986-1992 by Numerical Recipes Software. Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin g of machine- readable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes booksor CDROMs, v isit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to [email protected] (outside North Amer ica).should provide (i) parameters, (ii) error estimates on the parameters, and (iii) a statistical measure of goodness-of-fit. When the third item suggests that the modelis an unlikely match to the data, then items (i) and (ii) are probably worthless. Unfortunately, many practitioners of parameter estimation never proceed beyond item(i). Theydeemafitacceptableifagraphofdataandmodel“looksgood.” Thisapproachis knownas chi-by-eye . Luckily,its practitionersget what theydeserve. CITED REFERENCES AND FURTHER READING: Bevington, P.R. 1969, Data Reduction and Error Analysis for the Physical Sciences (New York: McGraw-Hill). Brownlee, K.A. 1965, Statistical Theory and Methodology , 2nd ed. (New York: Wiley). Martin, B.R. 1971, Statistics for Physicists (New York: Academic Press). von Mises, R. 1964, Mathematical Theory of Probability and Statistics (New York: Academic Press), Chapter X. Korn, G.A., andKorn, T.M. 1968, Mathematical Handbookfor Scientists andEngineers , 2nded. (New York: McGraw-Hill), Chapters 18–19. 15.1 Least Squares as a Maximum Likelihood Estimator Supposethatwearefitting Ndatapoints (xi,yi)i=1 ,...,N,toamodelthat has Madjustable parameters aj,j =1 ,...,M. The model predicts a functional relationship between the measured independentand dependent variables, y(x)= y(x;a1...a M)( 15.1.1 ) wherethedependenceontheparametersisindicatedexplicitlyontheright-handside. What, exactly, do we want to minimize to get fitted values for the aj’s? The first thing that comes to mind is the familiar least-squares fit, minimizeover a1...a M:N/summationdisplay i=1[yi−y(xi;a1...a M)]2(15.1.2 ) Butwheredoesthiscomefrom? Whatgeneralprinciplesisitbasedon? Theanswer to these questions takes us into the subject of maximum likelihoodestimators . Given a particular data set of xi’s and yi’s, we have the intuitive feeling that some parameter sets a1...a Mare very unlikely — those for which the model function y(x)looksnothinglike thedata—whileothersmaybeverylikely—those thatcloselyresemblethedata. Howcanwequantifythisintuitivefeeling? Howcan we select fitted parameters that are “most likely” to be correct? It is not meaningfulto ask thequestion,“What is the probabilitythat a particularset offitted parameters a 1...a Mis correct?” The reason is that there is no statistical universe of models fromwhich the parametersare drawn. Thereis just one model, the correct one, and a statistical universe of data sets that are drawn from it!