Logo ROOT  
Reference Guide
rf601_intminuit.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4##
5## 'LIKELIHOOD AND MINIMIZATION' RooFit tutorial macro #601
6##
7## Interactive minimization with MINUIT
8##
9## \macro_code
10##
11## \date February 2018
12## \authors Clemens Lange, Wouter Verkerke (C version)
13
14
15import ROOT
16
17
18# Setup pdf and likelihood
19# -----------------------------------------------
20
21# Observable
22x = ROOT.RooRealVar("x", "x", -20, 20)
23
24# Model (intentional strong correlations)
25mean = ROOT.RooRealVar("mean", "mean of g1 and g2", 0)
26sigma_g1 = ROOT.RooRealVar("sigma_g1", "width of g1", 3)
27g1 = ROOT.RooGaussian("g1", "g1", x, mean, sigma_g1)
28
29sigma_g2 = ROOT.RooRealVar("sigma_g2", "width of g2", 4, 3.0, 6.0)
30g2 = ROOT.RooGaussian("g2", "g2", x, mean, sigma_g2)
31
32frac = ROOT.RooRealVar("frac", "frac", 0.5, 0.0, 1.0)
33model = ROOT.RooAddPdf("model", "model", [g1, g2], [frac])
34
35# Generate 1000 events
36data = model.generate({x}, 1000)
37
38# Construct unbinned likelihood of model w.r.t. data
39nll = model.createNLL(data)
40
41# Interactive minimization, error analysis
42# -------------------------------------------------------------------------------
43
44# Create MINUIT interface object
45m = ROOT.RooMinimizer(nll)
46
47# Activate verbose logging of MINUIT parameter space stepping
48m.setVerbose(True)
49
50# Call MIGRAD to minimize the likelihood
51m.migrad()
52
53# Print values of all parameters, reflect values (and error estimates)
54# that are back propagated from MINUIT
55model.getParameters({x}).Print("s")
56
57# Disable verbose logging
58m.setVerbose(False)
59
60# Run HESSE to calculate errors from d2L/dp2
61m.hesse()
62
63# Print value (and error) of sigma_g2 parameter, reflects
64# value and error back propagated from MINUIT
65sigma_g2.Print()
66
67# Run MINOS on sigma_g2 parameter only
68m.minos({sigma_g2})
69
70# Print value (and error) of sigma_g2 parameter, reflects
71# value and error back propagated from MINUIT
72sigma_g2.Print()
73
74# Saving results, contour plots
75# ---------------------------------------------------------
76
77# Save a snapshot of the fit result. ROOT.This object contains the initial
78# fit parameters, final fit parameters, complete correlation
79# matrix, EDM, minimized FCN , last MINUIT status code and
80# the number of times the ROOT.RooFit function object has indicated evaluation
81# problems (e.g. zero probabilities during likelihood evaluation)
82r = m.save()
83
84# Make contour plot of mx vs sx at 1,2, sigma
85frame = m.contour(frac, sigma_g2, 1, 2, 3)
86frame.SetTitle("Contour plot")
87
88# Print the fit result snapshot
89r.Print("v")
90
91# Change parameter values, plotting
92# -----------------------------------------------------------------
93
94# At any moment you can manually change the value of a (constant)
95# parameter
96mean.setVal(0.3)
97
98# Rerun MIGRAD,HESSE
99m.migrad()
100m.hesse()
101frac.Print()
102
103# Now fix sigma_g2
104sigma_g2.setConstant(True)
105
106# Rerun MIGRAD,HESSE
107m.migrad()
108m.hesse()
109frac.Print()
110
111c = ROOT.TCanvas("rf601_intminuit", "rf601_intminuit", 600, 600)
112ROOT.gPad.SetLeftMargin(0.15)
113frame.GetYaxis().SetTitleOffset(1.4)
114frame.Draw()
115
116c.SaveAs("rf601_intminuit.png")
void Print(std::ostream &os, const OptionType &opt)