Likelihood and minimization: demonstration of options of the RooFitResult class
Author: Clemens Lange, Wouter Verkerke (C++ version)
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Sunday, February 05, 2023 at 11:16 AM.
from __future__ import print_function
import ROOT
Welcome to JupyROOT 6.29/01
Declare observable x
x = ROOT.RooRealVar("x", "x", 0, 10)
Create two Gaussian PDFs g1(x,mean1,sigma) anf g2(x,mean2,sigma) and their parameters
mean = ROOT.RooRealVar("mean", "mean of gaussians", 5, -10, 10)
sigma1 = ROOT.RooRealVar("sigma1", "width of gaussians", 0.5, 0.1, 10)
sigma2 = ROOT.RooRealVar("sigma2", "width of gaussians", 1, 0.1, 10)
sig1 = ROOT.RooGaussian("sig1", "Signal component 1", x, mean, sigma1)
sig2 = ROOT.RooGaussian("sig2", "Signal component 2", x, mean, sigma2)
Build Chebychev polynomial pdf
a0 = ROOT.RooRealVar("a0", "a0", 0.5, 0.0, 1.0)
a1 = ROOT.RooRealVar("a1", "a1", -0.2)
bkg = ROOT.RooChebychev("bkg", "Background", x, [a0, a1])
Sum the signal components into a composite signal pdf
sig1frac = ROOT.RooRealVar("sig1frac", "fraction of component 1 in signal", 0.8, 0.0, 1.0)
sig = ROOT.RooAddPdf("sig", "Signal", [sig1, sig2], [sig1frac])
Sum the composite signal and background
bkgfrac = ROOT.RooRealVar("bkgfrac", "fraction of background", 0.5, 0.0, 1.0)
model = ROOT.RooAddPdf("model", "g1+g2+a", [bkg, sig], [bkgfrac])
Generate 1000 events
data = model.generate({x}, 1000)
Perform fit and save result
r = model.fitTo(data, Save=True)
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization [#1] INFO:Minimization -- The following expressions will be evaluated in cache-and-track mode: (bkg,sig1,sig2) ********** ** 1 **SET PRINT 1 ********** ********** ** 2 **SET NOGRAD ********** PARAMETER DEFINITIONS: NO. NAME VALUE STEP SIZE LIMITS 1 a0 5.00000e-01 1.00000e-01 0.00000e+00 1.00000e+00 2 bkgfrac 5.00000e-01 1.00000e-01 0.00000e+00 1.00000e+00 3 mean 5.00000e+00 2.00000e+00 -1.00000e+01 1.00000e+01 4 sig1frac 8.00000e-01 1.00000e-01 0.00000e+00 1.00000e+00 5 sigma1 5.00000e-01 2.00000e-01 1.00000e-01 1.00000e+01 6 sigma2 1.00000e+00 4.50000e-01 1.00000e-01 1.00000e+01 ********** ** 3 **SET ERR 0.5 ********** ********** ** 4 **SET PRINT 1 ********** ********** ** 5 **SET STR 1 ********** NOW USING STRATEGY 1: TRY TO BALANCE SPEED AGAINST RELIABILITY ********** ** 6 **MIGRAD 3000 1 ********** FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4. START MIGRAD MINIMIZATION. STRATEGY 1. CONVERGENCE WHEN EDM .LT. 1.00e-03 FCN=1890.61 FROM MIGRAD STATUS=INITIATE 20 CALLS 21 TOTAL EDM= unknown STRATEGY= 1 NO ERROR MATRIX EXT PARAMETER CURRENT GUESS STEP FIRST NO. NAME VALUE ERROR SIZE DERIVATIVE 1 a0 5.00000e-01 1.00000e-01 2.01358e-01 -1.79143e+01 2 bkgfrac 5.00000e-01 1.00000e-01 2.01358e-01 8.01836e+00 3 mean 5.00000e+00 2.00000e+00 2.35352e-01 -3.18732e+02 4 sig1frac 8.00000e-01 1.00000e-01 2.57889e-01 1.67753e+00 5 sigma1 5.00000e-01 2.00000e-01 1.06123e-01 -2.80511e+01 6 sigma2 1.00000e+00 4.50000e-01 1.63378e-01 -2.79262e+00 ERR DEF= 0.5 MIGRAD MINIMIZATION HAS CONVERGED. MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX. COVARIANCE MATRIX CALCULATED SUCCESSFULLY FCN=1885.34 FROM MIGRAD STATUS=CONVERGED 177 CALLS 178 TOTAL EDM=0.000199951 STRATEGY= 1 ERROR MATRIX ACCURATE EXT PARAMETER STEP FIRST NO. NAME VALUE ERROR SIZE DERIVATIVE 1 a0 7.28245e-01 1.10053e-01 4.01494e-03 -5.30986e-04 2 bkgfrac 4.34386e-01 8.18255e-02 1.44990e-03 -2.37203e-02 3 mean 5.03463e+00 3.36192e-02 1.15677e-04 -1.36017e-01 4 sig1frac 7.78347e-01 9.66774e-02 3.39065e-03 -2.60523e-02 5 sigma1 5.23396e-01 4.47433e-02 4.30203e-04 -2.02274e-01 6 sigma2 1.77668e+00 1.13135e+00 3.09140e-03 2.51106e-02 ERR DEF= 0.5 EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 6 ERR DEF=0.5 1.237e-02 -7.217e-03 -1.043e-04 -4.043e-03 1.989e-03 1.052e-01 -7.217e-03 6.757e-03 -1.369e-04 4.830e-03 -1.372e-03 -8.173e-02 -1.043e-04 -1.369e-04 1.130e-03 -2.833e-04 -6.480e-05 3.890e-04 -4.043e-03 4.830e-03 -2.833e-04 9.520e-03 1.344e-03 -2.831e-02 1.989e-03 -1.372e-03 -6.480e-05 1.344e-03 2.002e-03 2.936e-02 1.052e-01 -8.173e-02 3.890e-04 -2.831e-02 2.936e-02 1.322e+00 PARAMETER CORRELATION COEFFICIENTS NO. GLOBAL 1 2 3 4 5 6 1 0.84288 1.000 -0.790 -0.028 -0.373 0.400 0.823 2 0.95650 -0.790 1.000 -0.050 0.602 -0.373 -0.865 3 0.13243 -0.028 -0.050 1.000 -0.086 -0.043 0.010 4 0.87689 -0.373 0.602 -0.086 1.000 0.308 -0.252 5 0.76698 0.400 -0.373 -0.043 0.308 1.000 0.571 6 0.94237 0.823 -0.865 0.010 -0.252 0.571 1.000 ********** ** 7 **SET ERR 0.5 ********** ********** ** 8 **SET PRINT 1 ********** ********** ** 9 **HESSE 3000 ********** COVARIANCE MATRIX CALCULATED SUCCESSFULLY FCN=1885.34 FROM HESSE STATUS=OK 40 CALLS 218 TOTAL EDM=0.000205499 STRATEGY= 1 ERROR MATRIX ACCURATE EXT PARAMETER INTERNAL INTERNAL NO. NAME VALUE ERROR STEP SIZE VALUE 1 a0 7.28245e-01 1.11109e-01 1.60598e-04 4.74047e-01 2 bkgfrac 4.34386e-01 8.36079e-02 2.89981e-04 -1.31608e-01 3 mean 5.03463e+00 3.36219e-02 2.31353e-05 5.27602e-01 4 sig1frac 7.78347e-01 9.69912e-02 6.78131e-04 5.90402e-01 5 sigma1 5.23396e-01 4.51307e-02 8.60406e-05 -1.15419e+00 6 sigma2 1.77668e+00 1.15533e+00 6.18281e-04 -7.22519e-01 ERR DEF= 0.5 EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 6 ERR DEF=0.5 1.261e-02 -7.502e-03 -9.635e-05 -4.154e-03 2.084e-03 1.091e-01 -7.502e-03 7.058e-03 -1.441e-04 4.954e-03 -1.470e-03 -8.595e-02 -9.635e-05 -1.441e-04 1.130e-03 -2.873e-04 -6.383e-05 4.916e-04 -4.154e-03 4.954e-03 -2.873e-04 9.583e-03 1.310e-03 -3.000e-02 2.084e-03 -1.470e-03 -6.383e-05 1.310e-03 2.037e-03 3.075e-02 1.091e-01 -8.595e-02 4.916e-04 -3.000e-02 3.075e-02 1.380e+00 PARAMETER CORRELATION COEFFICIENTS NO. GLOBAL 1 2 3 4 5 6 1 0.84619 1.000 -0.795 -0.026 -0.378 0.411 0.827 2 0.95839 -0.795 1.000 -0.051 0.602 -0.388 -0.871 3 0.13303 -0.026 -0.051 1.000 -0.087 -0.042 0.012 4 0.87775 -0.378 0.602 -0.087 1.000 0.297 -0.261 5 0.77155 0.411 -0.388 -0.042 0.297 1.000 0.580 6 0.94489 0.827 -0.871 0.012 -0.261 0.580 1.000 [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
Summary printing: Basic info plus final values of floating fit parameters
r.Print()
RooFitResult: minimized FCN value: 1885.34, estimated distance to minimum: 0.000205499 covariance matrix quality: Full, accurate covariance matrix Status : MINIMIZE=0 HESSE=0 Floating Parameter FinalValue +/- Error -------------------- -------------------------- a0 7.2825e-01 +/- 1.11e-01 bkgfrac 4.3439e-01 +/- 8.36e-02 mean 5.0346e+00 +/- 3.36e-02 sig1frac 7.7835e-01 +/- 9.70e-02 sigma1 5.2340e-01 +/- 4.51e-02 sigma2 1.7767e+00 +/- 1.16e+00
Verbose printing: Basic info, of constant parameters, and final values of floating parameters, correlations
r.Print("v")
RooFitResult: minimized FCN value: 1885.34, estimated distance to minimum: 0.000205499 covariance matrix quality: Full, accurate covariance matrix Status : MINIMIZE=0 HESSE=0 Constant Parameter Value -------------------- ------------ a1 -2.0000e-01 Floating Parameter InitialValue FinalValue +/- Error GblCorr. -------------------- ------------ -------------------------- -------- a0 5.0000e-01 7.2825e-01 +/- 1.11e-01 <none> bkgfrac 5.0000e-01 4.3439e-01 +/- 8.36e-02 <none> mean 5.0000e+00 5.0346e+00 +/- 3.36e-02 <none> sig1frac 8.0000e-01 7.7835e-01 +/- 9.70e-02 <none> sigma1 5.0000e-01 5.2340e-01 +/- 4.51e-02 <none> sigma2 1.0000e+00 1.7767e+00 +/- 1.16e+00 <none>
Construct 2D color plot of correlation matrix
ROOT.gStyle.SetOptStat(0)
ROOT.gStyle.SetPalette(1)
hcorr = r.correlationHist()
Visualize ellipse corresponding to single correlation matrix element
frame = ROOT.RooPlot(sigma1, sig1frac, 0.45, 0.60, 0.65, 0.90)
frame.SetTitle("Covariance between sigma1 and sig1frac")
r.plotOn(frame, sigma1, sig1frac, "ME12ABHV")
<cppyy.gbl.RooPlot object at 0x10dc2960>
Access basic information
print("EDM = ", r.edm())
print("-log(L) minimum = ", r.minNll())
EDM = 0.00020549915143065197 -log(L) minimum = 1885.343815305069
Access list of final fit parameter values
print("final value of floating parameters")
r.floatParsFinal().Print("s")
final value of floating parameters 1) RooRealVar:: a0 = 0.728245 +/- 0.111109 2) RooRealVar:: bkgfrac = 0.434386 +/- 0.0836079 3) RooRealVar:: mean = 5.03463 +/- 0.0336219 4) RooRealVar:: sig1frac = 0.778347 +/- 0.0969912 5) RooRealVar:: sigma1 = 0.523396 +/- 0.0451307 6) RooRealVar:: sigma2 = 1.77668 +/- 1.15533
Access correlation matrix elements
print("correlation between sig1frac and a0 is ", r.correlation(sig1frac, a0))
print("correlation between bkgfrac and mean is ", r.correlation("bkgfrac", "mean"))
correlation between sig1frac and a0 is -0.3778511951671703 correlation between bkgfrac and mean is -0.05102316430532121
Extract covariance and correlation matrix as ROOT.TMatrixDSym
cor = r.correlationMatrix()
cov = r.covarianceMatrix()
Print correlation, matrix
print("correlation matrix")
cor.Print()
print("covariance matrix")
cov.Print()
correlation matrix covariance matrix 6x6 matrix is as follows | 0 | 1 | 2 | 3 | 4 | ---------------------------------------------------------------------- 0 | 1 -0.7952 -0.02552 -0.3779 0.4111 1 | -0.7952 1 -0.05102 0.6023 -0.3876 2 | -0.02552 -0.05102 1 -0.0873 -0.04206 3 | -0.3779 0.6023 -0.0873 1 0.2966 4 | 0.4111 -0.3876 -0.04206 0.2966 1 5 | 0.8272 -0.8708 0.01245 -0.2609 0.5799 | 5 | ---------------------------------------------------------------------- 0 | 0.8272 1 | -0.8708 2 | 0.01245 3 | -0.2609 4 | 0.5799 5 | 1 6x6 matrix is as follows | 0 | 1 | 2 | 3 | 4 | ---------------------------------------------------------------------- 0 | 0.01261 -0.007502 -9.635e-05 -0.004154 0.002084 1 | -0.007502 0.007058 -0.0001441 0.004954 -0.00147 2 | -9.635e-05 -0.0001441 0.00113 -0.0002873 -6.383e-05 3 | -0.004154 0.004954 -0.0002873 0.009583 0.00131 4 | 0.002084 -0.00147 -6.383e-05 0.00131 0.002037 5 | 0.1091 -0.08595 0.0004916 -0.03 0.03075 | 5 | ---------------------------------------------------------------------- 0 | 0.1091 1 | -0.08595 2 | 0.0004916 3 | -0.03 4 | 0.03075 5 | 1.38
Open ROOT file save save result
f = ROOT.TFile("rf607_fitresult.root", "RECREATE")
r.Write("rf607")
f.Close()
In a clean ROOT session retrieve the persisted fit result as follows: r = gDirectory.Get("rf607")
c = ROOT.TCanvas("rf607_fitresult", "rf607_fitresult", 800, 400)
c.Divide(2)
c.cd(1)
ROOT.gPad.SetLeftMargin(0.15)
hcorr.GetYaxis().SetTitleOffset(1.4)
hcorr.Draw("colz")
c.cd(2)
ROOT.gPad.SetLeftMargin(0.15)
frame.GetYaxis().SetTitleOffset(1.6)
frame.Draw()
c.SaveAs("rf607_fitresult.png")
Info in <TCanvas::Print>: png file rf607_fitresult.png has been created
Draw all canvases
from ROOT import gROOT
gROOT.GetListOfCanvases().Draw()