[SciPy-user] incorrect variance in optimize.curvefit and leastsq
Fri Feb 13 01:06:13 CST 2009
On Fri, Feb 13, 2009 at 00:08, <firstname.lastname@example.org> wrote:
> On Fri, Feb 13, 2009 at 12:24 AM, Travis E. Oliphant
> <email@example.com> wrote:
>> firstname.lastname@example.org wrote:
>>> On Thu, Feb 12, 2009 at 11:09 PM, Travis E. Oliphant
>>> <email@example.com> wrote:
>>>> firstname.lastname@example.org wrote:
>>>>> I just saw the new optimize.curvefit which provides a wrapper around
>>> the standard deviation of the error can be calculated and the
>>> corrected (this is written for the use from outside of curvefit):
>>> yhat = func(x,popt, popt) # get predicted observations
>>> SSE = np.sum((y-yhat)**2)
>>> sig2 = SSE/(len(y)-len(popt))
>>> ecov = sig2*pcov # this is the variance-covariance matrix of
>>> the parameter estimates
>>> inside curvefit, this work (before the return):
>>> err = func(popt, *args)
>>> SSE = np.sum((err)**2)
>>> sig2 = SSE / (len(ydata) - len(popt))
>>> pcov = sig2 * pcov
> I think for the weighted least squares problem the weights should go
> into the SSE calculation. I didn't find directly the reference, but I
> am somewhat confident that this is correct, from the analogy to the
> model ynew = y*weight where weight_i = 1/sigma_i in the linear case.
> But it's too late today to try to check this.
> SSE = np.sum((err*weight)**2)
Yes. Basically, you want a Chi-squared statistic for the residuals.
scipy.odr does this scaling, for example.
"I have come to believe that the whole world is an enigma, a harmless
enigma that is made terrible by our own mad attempt to interpret it as
though it had an underlying truth."
-- Umberto Eco
More information about the SciPy-user