[SciPy-user] Puzzle with Monte Carlo error estimation

Iain Day iain at day-online.org.uk.invalid
Sun Nov 12 15:03:06 CST 2006


Iain Day wrote:
> Hi,
> 
> I've got a puzzle doing some Monte Carlo stuff to get standard 
> deviations for some fitting parameters I've calculated. The code is 
> below. The problem is that p1 is changing each time optimize.leastsq is 
> called in the loop and I can't see why. Any ideas?
> 
> Thanks,
> 
> Iain

Sorry, there is a typo here:

> for i in range(num_mc):
>      mock_data = buildup(p1, fix0, time_points) + sigma * 
> rand(len(time_points))
>      p2, success = optimize.leastsq(residuals, p1, args=(fix0, ts, expdata))

replace line above with

       p2, success = optimize.leastsq(residuals, p1, args=(fix0, pts, 
mock_data))




>      mocked_parameters[i] = p2
> 
> mc_errors = std(mocked_parameters, axis=0)
> 
> 
> def buildup(varpars, fixpars, time):
>      amplitude, tau, offset = varpars
>      buildup = amplitude * (1.0 - exp(- time / tau)) + offset
>      return buildup
> 
> def residuals(varpars, fixpars, time, expdata):
>     err = expdata - buildup(varpars, fixpars, time)
>     return err



More information about the SciPy-user mailing list