[SciPy-user] negative values in diagonal of covariance matrix

Pavlo Shchelokovskyy shchelokovskyy@gmail....
Thu Dec 11 07:09:14 CST 2008


Hi all,

I'm a moderately new user of scipy, trying to make some curve-fitting
with it. I wanted to use cov matrix output of leastsq to estimate
errors of fitted parameters, but stumbled upon strange discrepancy
(for one particular dataset):

On Linux (Fedora 8), using Python 2.5.1, numpy 1.2.0 and scipy 0.6.0

>>> from scipy import *
>>> from scipy import optimize
>>> y = asarray([217, 182, 162, 170, 255])
>>> x = linspace(0, y.size - 1, y.size)
>>> gauss = lambda p, x: p[0] + p[1] * exp(-(x-p[2])**2/(2*p[3]**2))
>>> errfunc = lambda p, x, y: y - gauss(p, x)
>>> pinit = (max(y), -y.ptp(), argmin(y), y.size/4)
>>> fit,cov,info,mesg,success = optimize.leastsq(errfunc, pinit, args=(x,y), full_output = 1)
>>> fit
array([  3.01602865e+05,  -3.01444487e+05,   1.83283239e+00, 8.87273494e+01])
>>> cov
array([[ -2.27903048e+16,   2.27903047e+16,  -4.72378743e+05, -3.35454486e+12],
       [  2.27903047e+16,  -2.27903046e+16,   4.72378733e+05, 3.35454485e+12],
       [ -4.72378950e+05,   4.72378947e+05,   6.38886491e-05, -6.95317302e+01],
       [ -3.35454486e+12,   3.35454485e+12,  -6.95316881e+01, -4.93761332e+08]])

I know that the fit is very bad (I can reject it afterwards in my
procedure), but what bothers me the most are the negative numbers on
the diagonal of the covariance matrix - as far as I understand, there
shouldn't be any... especially when the same script run on Windows XP,
using Python 2.5.2, numpy 1.2.0 and scipy 0.6.0 gives no such values
>>> fit

array([  3.58426106e+05,  -3.58267727e+05,   1.83283417e+00, 9.67300286e+01])

>>> cov

array([[  6.40320190e+14,  -6.40320188e+14,   4.24641077e+04, 8.64500286e+10],

       [ -6.40320188e+14,   6.40320186e+14,  -4.24641083e+04, -8.64500284e+10],

       [  4.24641250e+04,  -4.24641094e+04,   7.64931901e-05, 5.73152924e+00],

       [  8.64500286e+10,  -8.64500284e+10,   5.73152989e+00, 1.16716728e+07]])

Is it my error or misunderstanding somewhere, or is it really a bug in
Linux implementation?

Thanks in advance,


More information about the SciPy-user mailing list