Example of how to fit linear and nonlinear models to data with ROOT

fit_with_root.py

"""
Example of how to fit linear and nonlinear models to data with ROOT [1],
the most widely used data analysis / fitting package for physicists.

It would be great if statsmodels [2] could also fit non-linear models
and if the documentation could be improved to be understandable by physicists.

[1] http://root.cern.ch
[2] http://statsmodels.sourceforge.net
"""

import numpy as np
import ROOT

# Some example data
# http://www.itl.nist.gov/div898/strd/nls/data/daniel_wood.shtml
x = np.array([1.309, 1.471, 1.490, 1.565, 1.611, 1.680])
y = np.array([2.138, 3.421, 3.597, 4.340, 4.882, 5.660])
data = ROOT.TGraph(len(x), x, y)

# Define and fit a linear model (first degree polynomial)
model = ROOT.TF1('linear_model', 'pol1')
data.Fit(model)
model.Print()

# Define and fit a non-linear model
model = ROOT.TF1('nonlinear_model', '[0] * pow(x, [1])')
data.Fit(model)
model.Print()

"""
Here's the output you get when running the script:

****************************************
Minimizer is Linear
Chi2                      =    0.0746053
NDf                       =            4
p0                        =      -10.427   +/-   0.720062    
p1                        =      9.48935   +/-   0.471993    
         linear_model : pol1 Ndim= 1, Npar= 2, Noper= 1
 fExpr[0] = pol1  action = 130 action param = 101 
Optimized expression
 fExpr[0] = pol1          action = 130 action param = 101 
 Par  0                    p0 = -10.427
 Par  1                    p1 = 9.48935

****************************************
Minimizer is Minuit / Migrad
Chi2                      =   0.00431731
NDf                       =            4
Edm                       =   6.8811e-10
NCalls                    =          168
p0                        =     0.768872   +/-   0.01833     
p1                        =      3.86038   +/-   0.0518581   
      nonlinear_model : [0]*pow(x,[1]) Ndim= 1, Npar= 2, Noper= 5
 fExpr[0] = [0]  action = 140 action param = 0 
 fExpr[1] = x  action = 144 action param = 0 
 fExpr[2] = [1]  action = 140 action param = 1 
 fExpr[3] = ^  action = 20 action param = 0 
 fExpr[4] = *  action = 3 action param = 0 
Optimized expression
 fExpr[0] = [0]          action = 146 action param = 0 
 fExpr[1] = x          action = 146 action param = 0 
 fExpr[2] = [1]          action = 146 action param = 1 
 fExpr[3] = ^          action = 20 action param = 0 
 fExpr[4] = *          action = 3 action param = 0 
 Par  0                    p0 = 0.768872
 Par  1                    p1 = 3.86038
"""