Bootstrap

Bootstrapping is a method to estimate the statistically distribution of a random variable that depends on a bunch of input random variables. The idea is to sample the population of the input random variables to build the population of the output variable. In the following example, I will fit the curve (x,y) with error yerr, to a third order polynomial and extract the position of the minimum of the curve. Helpful post.

import numpy as np
import scipy as pd

# 1. define the model function with minimum as a fit parameter
def third(x,Re,a,b,c):
  return a*(x-Re)**3.+b*(x-Re)**2.+c
# end def

model = third

# 2. fit to extract initial guess
popt,pcov = op.curve_fit(model,x,y)

# 3. bootstrap
def bootstrap(x,y,yerr,model,p0,nsample=100):
  errfunc = lambda p,x,y: y-third(x,*p)
  param_trace = np.array([p0]*nsample)
  for i in range(1,nsample):
    deviation = np.random.randn(len(yerr))*yerr
    ydata = y + deviation
    param_trace[i] = op.leastsq(errfunc,p0,args=(x,ydata))[0]
  # end for i
return param_trace
# end def bootstrap

params = bootstrap(x,y,yerr,fifth,popt)
Re = params.T[0]
print “Re=%1.2f +- %1.2f” % (Re.mean(),Re.std())