Re: Calculating Error Estimates [message #12302 is a reply to message #12280] |
Tue, 14 July 1998 00:00   |
Wayne Landsman
Messages: 117 Registered: January 1997
|
Senior Member |
|
|
Before giving my answer David Fanning's question about estimating
errors from
CURVEFIT, I first want to post a related complaint about the procedure
SVDFIT. SVDFIT has an output parameter SIGMA which according to the
documentation gives the "1-sigma error estimates of the returned
parameters." What the documentation doesn't say is that unless you
supply
values in the optional WEIGHTS input keyword, then the values returned
by SIGMA
are likely to be complete nonsense.
In order to compute the uncertainties in the derived parameters, SVDFIT
needs to
know the uncertainties in the individual data points. If the user
does not
supply these, then SVDFIT assume that the uncertainties are all equal
(which is
reasonable), and that they are all equal to 1 unit (which is likely to
be very
unreasonable).
What one should do -- and what David F. should do for his example -- is
to
estimate the uncertainties by comparing the data with the best fit.
This is
done correctly by the procedure LINFIT.
SIGMA = SQRT( (TOTAL ((Y - YFIT)^2) / (NPTS - NPARAMS)))
where NPTS and NPARAMS are respectively the number of data points and
number of
parameters used in the fit.
As an example of the problem with SVDFIT, consider a linear fit with an
uncertainty of 0.01 in the Y coordinate
STIS> x = indgen(10)
STIS> y = x + 0.01*randomn(seed,10)
First, I try a linear fit with LINFIT and print the derived parameters
and
their sigma values
IDL> param = LINFIT(x, y, sigma = sig)
IDL> print,param
-0.00540641 1.00106
IDL> print, sig
0.00654651 0.00122627
Now try a linear fit (a polynomial with 2 parameters) of the same data
using SVDFIT
IDL> param = svdfit(x,y,2,sigma=sig)
IDL> print,param
-0.00540641 1.00106
IDL> print,sig
0.587754 0.110096
Although SVDFIT gives the right parameter values, the associated sigma
values
are much too large. This is because SVDFIT assumed individual errors
of 1
unit, whereas LINFIT computed
yfit = param(1)*x + param(0)
sigma = total( (yfit-y)^2)/8. ))
to estimate individual errors of 0.011 units
Wayne Landsman
landsman@mpb.gsfc.nasa.gov
|
|
|