Rf 6 0 1_Intminuit

'LIKELIHOOD AND MINIMIZATION' RooFit tutorial macro #601

Interactive minimization with MINUIT

Author: Wouter Verkerke (C version)
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, January 19, 2022 at 10:26 AM.

In [ ]:
import ROOT

Setup pdf and likelihood

Observable

In [ ]:
x = ROOT.RooRealVar("x", "x", -20, 20)

Model (intentional strong correlations)

In [ ]:
mean = ROOT.RooRealVar("mean", "mean of g1 and g2", 0)
sigma_g1 = ROOT.RooRealVar("sigma_g1", "width of g1", 3)
g1 = ROOT.RooGaussian("g1", "g1", x, mean, sigma_g1)

sigma_g2 = ROOT.RooRealVar("sigma_g2", "width of g2", 4, 3.0, 6.0)
g2 = ROOT.RooGaussian("g2", "g2", x, mean, sigma_g2)

frac = ROOT.RooRealVar("frac", "frac", 0.5, 0.0, 1.0)
model = ROOT.RooAddPdf("model", "model", [g1, g2], [frac])

Generate 1000 events

In [ ]:
data = model.generate({x}, 1000)

Construct unbinned likelihood of model w.r.t. data

In [ ]:
nll = model.createNLL(data)

Interactive minimization, error analysis

Create MINUIT interface object

In [ ]:
m = ROOT.RooMinimizer(nll)

Activate verbose logging of MINUIT parameter space stepping

In [ ]:
m.setVerbose(True)

Call MIGRAD to minimize the likelihood

In [ ]:
m.migrad()

Print values of all parameters, reflect values (and error estimates) that are back propagated from MINUIT

In [ ]:
model.getParameters({x}).Print("s")

Disable verbose logging

In [ ]:
m.setVerbose(False)

Run HESSE to calculate errors from d2L/dp2

In [ ]:
m.hesse()

Print value (and error) of sigma_g2 parameter, reflects value and error back propagated from MINUIT

In [ ]:
sigma_g2.Print()

Run MINOS on sigma_g2 parameter only

In [ ]:
m.minos({sigma_g2})

Print value (and error) of sigma_g2 parameter, reflects value and error back propagated from MINUIT

In [ ]:
sigma_g2.Print()

Saving results, contour plots

Save a snapshot of the fit result. ROOT.This object contains the initial fit parameters, final fit parameters, complete correlation matrix, EDM, minimized FCN , last MINUIT status code and the number of times the ROOT.RooFit function object has indicated evaluation problems (e.g. zero probabilities during likelihood evaluation)

In [ ]:
r = m.save()

Make contour plot of mx vs sx at 1,2, sigma

In [ ]:
frame = m.contour(frac, sigma_g2, 1, 2, 3)
frame.SetTitle("Contour plot")

Print the fit result snapshot

In [ ]:
r.Print("v")

Change parameter values, plotting

At any moment you can manually change the value of a (constant) parameter

In [ ]:
mean.setVal(0.3)

Rerun MIGRAD,HESSE

In [ ]:
m.migrad()
m.hesse()
frac.Print()

Now fix sigma_g2

In [ ]:
sigma_g2.setConstant(True)

Rerun MIGRAD,HESSE

In [ ]:
m.migrad()
m.hesse()
frac.Print()

c = ROOT.TCanvas("rf601_intminuit", "rf601_intminuit", 600, 600)
ROOT.gPad.SetLeftMargin(0.15)
frame.GetYaxis().SetTitleOffset(1.4)
frame.Draw()

c.SaveAs("rf601_intminuit.png")

Draw all canvases

In [ ]:
from ROOT import gROOT 
gROOT.GetListOfCanvases().Draw()