[SciPy-User] optimize.leastsq error estimation

Joe Philip Ninan indiajoe@gmail....
Fri Mar 22 03:58:06 CDT 2013


Hi Matt,
Thanks for clarifying and also for providing an example.
I think it will be nice if we add the term "reduced chi square" also along
with "reduced variance" in the scipy documentation.
And also an example to show how to estimate the error in the documentation
page.
Thanking you,
Joe


On 22 March 2013 08:08, Matt Newville <newville@cars.uchicago.edu> wrote:

> Hi Joe,
>
> On Thu, Mar 21, 2013 at 3:10 PM, Joe Philip Ninan <indiajoe@gmail.com>
> wrote:
> > Hi,
> > I want to confirm whether i understood the documentation of
> optimize.leastsq
> > correctly.
> > (
> >
> http://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.leastsq.html
> > )
> > For estimating the error in the fitted parameter.
> > The documentation says :
> > "This matrix (cov_x) must be multiplied by the residual variance to get
> the
> > covariance of the parameter estimates"
> > i.e. Multiply the output cov_x with the "reduced chi square". Right? and
> > take sqrt of corresponding diagonal element from covariance matrix.
> > Where "residual chi square" is the sum of {[(fitted function -
> > data)/sigma]**2} /N
> > and where N is degrees of freedom [= len(data)-number of parameters]
> > And sigma is the error in each data point.
> >
> > I mainly want to confirm whether the divide by sigma is required or not.
> > (I couldn't find the definition of "residual variance" in documentation,
> so
> > i assumed it is reduced chi square).
> >
> > Thanking you,
> > Joe
>
> I think you're correct, but there may be some confusion in
> terminology.  The standard errors are the square-root of the diagonal
> elements of (covar * reduced_chi_square), where reduced_chi_square is
>
>    (((data-model)/uncertainty)**2).sum() / (len(ydata)- len(variables))
>
> In (untested!) code,  it might look like:
>
>     bestvals, covar, info, errmsg, ier  = leastsq(objectivefunc, initvals,
>
>   args=args, full_output=1, **kws)
>
>     residual = objectivefunc(bestvals, *args)
>     reduced_chi_square = (residual**2).sum() / (len(data) - len(initvals))
>
>     covar = covar *  reduced_chi_square
>
>     std_error = array([sqrt(covar[i, i]) for i in range(len(bestvals))])
>
> Here, objectivefunc is expected to return the scaled misfit array,
> that is  "(data - model)/uncertainty".
>
> Cheers,
>
> --Matt Newville <newville at cars.uchicago.edu> 630-252-0431
> _______________________________________________
> SciPy-User mailing list
> SciPy-User@scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>



-- 
/---------------------------------------------------------------
"GNU/Linux: because a PC is a terrible thing to waste" -  GNU Generation
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/scipy-user/attachments/20130322/776dd533/attachment.html 


More information about the SciPy-User mailing list