$ontext Nonlinear Least Squares Regression example Erwin Kalvelagen, nov 2007 Reference: http://www.itl.nist.gov/div898/strd/nls/nls_main.shtml -------------------------------------------------------------------------- Procedure: Nonlinear Least Squares Regression Description: These data are the result of a NIST study regarding dental research in monomolecular adsorption. The response variable is volume, and the predictor variable is pressure. Reference: Misra, D., NIST (1978). Dental Research Monomolecular Adsorption Study. Data: 1 Response Variable (y = volume) 1 Predictor Variable (x = pressure) 14 Observations Lower Level of Difficulty Observed Data Model: Exponential Class 2 Parameters (b1 and b2) y = b1*(1-exp[-b2*x]) + e Starting values Certified Values Start 1 Start 2 Parameter Standard Deviation b1 = 500 250 2.3894212918E+02 2.7070075241E+00 b2 = 0.0001 0.0005 5.5015643181E-04 7.2668688436E-06 Residual Sum of Squares: 1.2455138894E-01 Residual Standard Deviation: 1.0187876330E-01 Degrees of Freedom: 12 Number of Observations: 14 $offtext *----------------------------------------------------------------------------- * data *----------------------------------------------------------------------------- set i /i1*i14/; table data(i,*) y x i1 10.07E0 77.6E0 i2 14.73E0 114.9E0 i3 17.94E0 141.1E0 i4 23.93E0 190.8E0 i5 29.61E0 239.9E0 i6 35.18E0 289.0E0 i7 40.02E0 332.8E0 i8 44.82E0 378.4E0 i9 50.76E0 434.8E0 i10 55.05E0 477.3E0 i11 61.01E0 536.8E0 i12 66.40E0 593.1E0 i13 75.47E0 689.1E0 i14 81.78E0 760.0E0 ; * * extract data * parameter x(i),y(i); x(i) = data(i,'x'); y(i) = data(i,'y'); * * certified values * scalars cb1 'certified value for b1' /2.3894212918E+02/ cb2 'certified value for b2' /5.5015643181E-04/ ce1 'certified std err for b1 ' / 2.7070075241E+00 / ce2 'certified std err for b2 ' / 7.2668688436E-06 / ; *----------------------------------------------------------------------------- * statistical model *----------------------------------------------------------------------------- variables sse 'sum of squared errors' b1 'coefficient to estimate' b2 'coefficient to estimate' ; equations fit(i) 'the non-linear model' obj 'objective' ; obj.. sse =n= 0; fit(i).. y(i) =e= b1*(1-exp[-b2*x(i)]); *----------------------------------------------------------------------------- * first set of initial values *----------------------------------------------------------------------------- b1.l = 500; b2.l = 0.0001; option nlp=nls; model nlfit /obj,fit/; solve nlfit minimizing sse using nlp; display sse.l,b1.l,b2.l; abort$((abs(b1.l-cb1)+abs(b2.l-cb2))>0.0001) "Accuracy problem"; abort$((abs(b1.m-ce1)+abs(b2.m-ce2))>0.0001) "Accuracy problem"; *----------------------------------------------------------------------------- * second set of initial values *----------------------------------------------------------------------------- b1.l = 250; b2.l = 0.0005; solve nlfit minimizing sse using nlp; display sse.l,b1.l,b2.l; abort$((abs(b1.l-cb1)+abs(b2.l-cb2))>0.0001) "Accuracy problem"; abort$((abs(b1.m-ce1)+abs(b2.m-ce2))>0.0001) "Accuracy problem";