[SciPy-user] incorrect variance in optimize.curvefit and leastsq

[Travis E. Oliphant]

Travis E. Oliphant wrote:
> josef.pktd at gmail.com wrote:
>   
>> I just saw the new optimize.curvefit which provides a wrapper around
>> optimize.leastsq
>>
>> optimize.leastsq provides the raw covariance matrix (cov_x). As I
>> mentioned once on the mailing list, this is not the covariance matrix
>> of the parameter estimates. To get those, the raw covariance matrix
>> has to be multiplied by the standard error of the residual. So, the
>> docstring in optimize.curvefit doesn't correspond to the actual
>> calculation.
>>   
>>     
> Thank you for the clarification.   I had forgotten your earlier valid 
> concerns.   Help fixing the docstring is appreciated.    If you can 
> figure out how to improve the code, that is even better.   I think it is 
> good to at least report the cov, but the docstring should not mislead.
>   
>> The first parameter is not exactly high precision.
>>
>> The second problem is that, in weighted least squares, the calculation
>> of the standard deviation of the parameter estimates has to take the
>> weights into account. (But I don't have the formulas right now.)
>>   
>> I was looking at this to provide a general non-linear least squares
>> class in stats. But for several calculation, the Jacobian would be
>> necessary. optimize.leastsq only provides cov_x, but I was wondering
>> whether the Jacobian can be calculated from the return of the minpack
>> functions in optimize.leastsq, but I didn't have time to figure this
>> out.
>>   
>>     
> I'm not sure, but it might be.    I would love to spend time on this, 
> but don't have it.   If somebody else can pick up, that would be great.

O.K.  So my desire to spend time on it outweighed my wisdom, and I went 
ahead and looked at the reference linked-to and multipled by the 
necessary scale factor.  I fixed the documentation in leastsq as well.   
   A sanity check on my work would be appreciated. 

I divided by the sum of the weights squared for the weighted case.   I'm 
not sure if this is correct, but it's probably close.  When someone can 
verify the formula that would be great.    Adding a check against the 
test case referred-to would be great.  

-Travis

-- 

Travis Oliphant
Enthought, Inc.
(512) 536-1057 (office)
(512) 536-1059 (fax)
http://www.enthought.com
oliphant at enthought.com