c15-5 (779589), страница 3
Текст из файла (страница 3)
Used with the above routine mrqmin (in turn using mrqcof, covsrt, andgaussj), it fits for the modely(x) =KXk=1" 2 #x − EkBk exp −Gk(15.5.16)which is a sum of K Gaussians, each having a variable position, amplitude, andwidth. We store the parameters in the order B1 , E1 , G1 , B2 , E2 , G2 , . . . , BK ,EK , GK .Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-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 copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).void mrqcof(float x[], float y[], float sig[], int ndata, float a[], int ia[],int ma, float **alpha, float beta[], float *chisq,void (*funcs)(float, float [], float *, float [], int))Used by mrqmin to evaluate the linearized fitting matrix alpha, and vector beta as in (15.5.8),and calculate χ2 .{int i,j,k,l,m,mfit=0;float ymod,wt,sig2i,dy,*dyda;688Chapter 15.Modeling of Data#include <math.h>*y=0.0;for (i=1;i<=na-1;i+=3) {arg=(x-a[i+1])/a[i+2];ex=exp(-arg*arg);fac=a[i]*ex*2.0*arg;*y += a[i]*ex;dyda[i]=ex;dyda[i+1]=fac/a[i+2];dyda[i+2]=fac*arg/a[i+2];}}More Advanced Methods for Nonlinear Least SquaresThe Levenberg-Marquardt algorithm can be implemented as a model-trustregion method for minimization (see §9.7 and ref.
[2]) applied to the special caseof a least squares function. A code of this kind due to Moré [3] can be found inMINPACK [4]. Another algorithm for nonlinear least-squares keeps the secondderivative term we dropped in the Levenberg-Marquardt method whenever it wouldbe better to do so. These methods are called “full Newton-type” methods and arereputed to be more robust than Levenberg-Marquardt, but more complex.
Oneimplementation is the code NL2SOL [5].CITED REFERENCES AND FURTHER READING:Bevington, P.R. 1969, Data Reduction and Error Analysis for the Physical Sciences (New York:McGraw-Hill), Chapter 11.Marquardt, D.W. 1963, Journal of the Society for Industrial and Applied Mathematics, vol. 11,pp. 431–441. [1]Jacobs, D.A.H.
(ed.) 1977, The State of the Art in Numerical Analysis (London: AcademicPress), Chapter III.2 (by J.E. Dennis).Dennis, J.E., and Schnabel, R.B. 1983, Numerical Methods for Unconstrained Optimization andNonlinear Equations (Englewood Cliffs, NJ: Prentice-Hall). [2]Moré, J.J. 1977, in Numerical Analysis, Lecture Notes in Mathematics, vol. 630, G.A.
Watson,ed. (Berlin: Springer-Verlag), pp. 105–116. [3]Moré, J.J., Garbow, B.S., and Hillstrom, K.E. 1980, User Guide for MINPACK-1, Argonne NationalLaboratory Report ANL-80-74. [4]Dennis, J.E., Gay, D.M, and Welsch, R.E. 1981, ACM Transactions on Mathematical Software,vol. 7, pp. 348–368; op. cit., pp.
369–383. [5].Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-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 copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited.
To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).void fgauss(float x, float a[], float *y, float dyda[], int na)y(x; a) is the sum of na/3 Gaussians (15.5.16). The amplitude, center, and width of theGaussians are stored in consecutive locations of a: a[i] = Bk , a[i+1] = Ek , a[i+2] = Gk ,k = 1, ..., na/3. The dimensions of the arrays are a[1..na], dyda[1..na].{int i;float fac,ex,arg;15.6 Confidence Limits on Estimated Model Parameters68915.6 Confidence Limits on Estimated ModelParametersMonte Carlo Simulation of Synthetic Data SetsAlthough the measured parameter set a(0) is not the true one, let us considera fictitious world in which it was the true one. Since we hope that our measuredparameters are not too wrong, we hope that that fictitious world is not too differentfrom the actual world with parameters atrue .
In particular, let us hope — no, let usassume — that the shape of the probability distribution a(i) − a(0) in the fictitiousworld is the same, or very nearly the same, as the shape of the probability distributionSample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-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 copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).Several times already in this chapter we have made statements about the standarderrors, or uncertainties, in a set of M estimated parameters a.
We have given someformulas for computing standard deviations or variances of individual parameters(equations 15.2.9, 15.4.15, 15.4.19), as well as some formulas for covariancesbetween pairs of parameters (equation 15.2.10; remark following equation 15.4.15;equation 15.4.20; equation 15.5.15).In this section, we want to be more explicit regarding the precise meaningof these quantitative uncertainties, and to give further information about howquantitative confidence limits on fitted parameters can be estimated.
The subjectcan get somewhat technical, and even somewhat confusing, so we will try to makeprecise statements, even when they must be offered without proof.Figure 15.6.1 shows the conceptual scheme of an experiment that “measures”a set of parameters. There is some underlying true set of parameters atrue that areknown to Mother Nature but hidden from the experimenter.
These true parametersare statistically realized, along with random measurement errors, as a measured dataset, which we will symbolize as D(0). The data set D(0) is known to the experimenter.He or she fits the data to a model by χ2 minimization or some other technique, andobtains measured, i.e., fitted, values for the parameters, which we here denote a(0) .Because measurement errors have a random component, D(0) is not a uniquerealization of the true parameters atrue . Rather, there are infinitely many otherrealizations of the true parameters as “hypothetical data sets” each of which couldhave been the one measured, but happened not to be.
Let us symbolize theseby D(1) , D(2), . . . . Each one, had it been realized, would have given a slightlydifferent set of fitted parameters, a(1), a(2), . . . , respectively. These parameter setsa(i) therefore occur with some probability distribution in the M -dimensional spaceof all possible parameter sets a.
The actual measured set a(0) is one member drawnfrom this distribution.Even more interesting than the probability distribution of a(i) would be thedistribution of the difference a(i) − atrue . This distribution differs from the formerone by a translation that puts Mother Nature’s true value at the origin. If we knew thisdistribution, we would know everything that there is to know about the quantitativeuncertainties in our experimental measurement a(0) .So the name of the game is to find some way of estimating or approximatingthe probability distribution of a(i) − atrue without knowing atrue and without havingavailable to us an infinite universe of hypothetical data sets..
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
