[SciPy-User] confidence interval for leastsq fit

[josef.pktd at gmail.com]

On Thu, Apr 28, 2011 at 8:12 AM, eat <e.antero.tammi at gmail.com> wrote:
>
>
> On Thu, Apr 28, 2011 at 3:01 PM, eat <e.antero.tammi at gmail.com> wrote:
>>
>> Hi
>>
>> On Thu, Apr 28, 2011 at 1:09 PM, <josef.pktd at gmail.com> wrote:
>>>
>>> On Wed, Apr 27, 2011 at 11:08 PM, dima osin <dima.osin at gmail.com> wrote:
>>> > How to calculate confidence interval for scipy.optimize.leastsq fit
>>> > using
>>> > the Student's t distribution and NOT the bootstrapping method?
>>> >
>>> > _______________________________________________
>>> > SciPy-User mailing list
>>> > SciPy-User at scipy.org
>>> > http://mail.scipy.org/mailman/listinfo/scipy-user
>>> >
>>> >
>>>
>>> if you mean the confidence interval for the parameter estimates
>>>
>>> formula and calulation for the covariance matrix of the parameter
>>> estimate is in optimize curvefit
>>> bse = np.sqrt(np.diag(cov)
>>> confint = params + np.array([-1, +1]) * stats.t.ppf(alpha/ 2) *
>>> bse[:,None]
>>
>> According to leastsq documentation (v0.10.dev) " cov_x : ndarray ... This
>> matrix must be multiplied by the residual standard deviation to get the
>> covariance of the parameter estimates – see curve_fit." However curve_fit
>> doc states "pcov : 2d array The estimated covariance of popt. The diagonals
>> provide the variance of the parameter estimate."
>> Shouldn't bse (based on leastsq) be claculated like:
>> bse= sqrt(sqrt(sum(residual** 2)/ (M- N))* dig(cov)) # assumed residual
>> zero mean normally distributed
>
> Oops
> meant to be like:
> bse= sqrt(sum(residual** 2)/ (M- N))* dig(cov) # assumed residual zero mean
> normally distributed

optimize.leastsq  returns cov_x which would be the cov in this
calculation and needs the multiplication.

optimize.curve_fit already does the multiplication with the sum
squared residual, so the cov that curve_fit returns is the actual
parameter covariance, which is the one I meant originally

    cf_params, cf_pcov = optimize.curve_fit(func0, x, y)
    cf_bse = np.sqrt(np.diag(cf_pcov))

Josef


>>
>> Regards,
>> eat
>>>
>>> (from scikits.statsmodels LikelihoodModel confint: )
>>>            lower = self.params[cols] - dist.ppf(1-\
>>>                        alpha/2,self.model.df_resid) * bse[cols]
>>>            upper = self.params[cols] + dist.ppf(1-\
>>>                        alpha/2,self.model.df_resid) * bse[cols]
>>>
>>> there is also a scikits.statsmodels.miscmodels.NonlinearLS which  is
>>> similar to scipy.optimize curvefit but a class with all regression
>>> results using the gradient information. Caution: doesn't have complete
>>> test suite yet, only partially tested on some examples.
>>>
>>> If you mean confidence interval for prediction, then I'm not
>>> completely sure the analogy to the linear models works unchanged.
>>>
>>> Josef
>>> _______________________________________________
>>> SciPy-User mailing list
>>> SciPy-User at scipy.org
>>> http://mail.scipy.org/mailman/listinfo/scipy-user
>>
>
>
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>
>