15.4 General Linear Least Squares 679 Feeding the matrix v and vector w output by the above program into the following short routine,you easily obtain variances and covariances of the fitted parameters a.The square roots of the variances are the standard deviations of the fitted parameters.The routine straightforwardly implements equation(15.4.20) above,with the convention that singular values equal to zero are recognized as having been edited out of the fit. #include "nrutil.h" void svdvar(float **v,int ma,float w[],float **cvm) To evaluate the covariance matrix cvm[1..ma][1..ma]of the fit for ma parameters obtained by svdfit,call this routine with matrices v[1..ma][1..ma],w[1..ma]as returned from 81 svdfit. int k,j,i; float sum,*wti; -00 from NUMERICAL 18881892 wti=vector(1,ma); for(i=1;i<ma;1++)[ wti[i]=0.0; RECIPES I if(w[i])ti[i]=1.0/(w[i]*w[i]); for(1=1;1<=ma;1++)[ Sum contributions to covariance matrix (15.4.20). for(j=1;j<=1;j++)[ for (sum=0.0,k=1;k<=ma;k++)sum +v[i][k]*v[j][k]*wti[k]; America computer, Press. cvm[j][i]=cvm[i][j]=sum; 9 ART free_vector(wti,1,ma); Programs 2 Examples to dir Be aware that some apparently nonlinear problems can be expressed so that 1992 SCIENTIFIC COMPUTING (ISBN they are linear.For example,an exponential model with two parameters a and b, y(r)=aexp(-bx) (15.4.21) can be rewritten as 10521 logy(x)】=c-bz (15.4.22) Numerical which is linear in its parameters c and b.(Of course you must be aware that such 43106 transformations do not exactly take Gaussian errors into Gaussian errors.) Recipes Also watch out for“non-parameters,.”asin y(x)=aexp(-bx+d) (15.4.23) North Here the parameters a and d are,in fact,indistinguishable.This is a good example of where the normal equations will be exactly singular,and where SVD will find a zero singular value.SVD will then make a"least-squares"choice for setting a balance between a and d (or,rather,their equivalents in the linear model derived by taking the logarithms).However-and this is true whenever SVD gives back a zero singular value-you are better advised to figure out analytically where the degeneracy is among your basis functions,and then make appropriate deletions in the basis set. Here are two examples for user-supplied routines funcs.The first one is trivial and fits a general polynomial to a set of data:15.4 General Linear Least Squares 679 Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). Feeding the matrix v and vector w output by the above program into the following short routine, you easily obtain variances and covariances of the fitted parameters a. The square roots of the variances are the standard deviations of the fitted parameters. The routine straightforwardly implements equation (15.4.20) above, with the convention that singular values equal to zero are recognized as having been edited out of the fit. #include "nrutil.h" void svdvar(float **v, int ma, float w[], float **cvm) To evaluate the covariance matrix cvm[1..ma][1..ma] of the fit for ma parameters obtained by svdfit, call this routine with matrices v[1..ma][1..ma], w[1..ma] as returned from svdfit. { int k,j,i; float sum,*wti; wti=vector(1,ma); for (i=1;i<=ma;i++) { wti[i]=0.0; if (w[i]) wti[i]=1.0/(w[i]*w[i]); } for (i=1;i<=ma;i++) { Sum contributions to covariance matrix (15.4.20). for (j=1;j<=i;j++) { for (sum=0.0,k=1;k<=ma;k++) sum += v[i][k]*v[j][k]*wti[k]; cvm[j][i]=cvm[i][j]=sum; } } free_vector(wti,1,ma); } Examples Be aware that some apparently nonlinear problems can be expressed so that they are linear. For example, an exponential model with two parameters a and b, y(x) = a exp(−bx) (15.4.21) can be rewritten as log[y(x)] = c − bx (15.4.22) which is linear in its parameters c and b. (Of course you must be aware that such transformations do not exactly take Gaussian errors into Gaussian errors.) Also watch out for “non-parameters,” as in y(x) = a exp(−bx + d) (15.4.23) Here the parameters a and d are, in fact, indistinguishable. This is a good example of where the normal equations will be exactly singular, and where SVD will find a zero singular value. SVD will then make a “least-squares” choice for setting a balance between a and d (or, rather, their equivalents in the linear model derived by taking the logarithms). However — and this is true whenever SVD gives back a zero singular value — you are better advised to figure out analytically where the degeneracy is among your basis functions, and then make appropriate deletions in the basis set. Here are two examples for user-supplied routines funcs. The first one is trivial and fits a general polynomial to a set of data: