ROOT   Reference Guide
rf610_visualerror.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4##
5## Likelihood and minimization: visualization of errors from a covariance matrix
6##
7## \macro_code
8##
9## \date February 2018
10## \authors Clemens Lange, Wouter Verkerke (C++ version)
11
12import ROOT
13
14
15# Setup example fit
16# ---------------------------------------
17
18# Create sum of two Gaussians p.d.f. with factory
19x = ROOT.RooRealVar("x", "x", -10, 10)
20
21m = ROOT.RooRealVar("m", "m", 0, -10, 10)
22s = ROOT.RooRealVar("s", "s", 2, 1, 50)
23sig = ROOT.RooGaussian("sig", "sig", x, m, s)
24
25m2 = ROOT.RooRealVar("m2", "m2", -1, -10, 10)
26s2 = ROOT.RooRealVar("s2", "s2", 6, 1, 50)
27bkg = ROOT.RooGaussian("bkg", "bkg", x, m2, s2)
28
29fsig = ROOT.RooRealVar("fsig", "fsig", 0.33, 0, 1)
31 sig, bkg), ROOT.RooArgList(fsig))
32
33# Create binned dataset
34x.setBins(25)
35d = model.generateBinned(ROOT.RooArgSet(x), 1000)
36
37# Perform fit and save fit result
38r = model.fitTo(d, ROOT.RooFit.Save())
39
40# Visualize fit error
41# -------------------------------------
42
43# Make plot frame
44frame = x.frame(ROOT.RooFit.Bins(40), ROOT.RooFit.Title(
45 "P.d.f with visualized 1-sigma error band"))
46d.plotOn(frame)
47
48# Visualize 1-sigma error encoded in fit result 'r' as orange band using linear error propagation
49# ROOT.This results in an error band that is by construction symmetric
50#
51# The linear error is calculated as
52# error(x) = Z* F_a(x) * Corr(a,a') F_a'(x)
53#
54# where F_a(x) = [ f(x,a+da) - f(x,a-da) ] / 2,
55#
56# with f(x) = the plotted curve
57# 'da' = error taken from the fit result
58# Corr(a,a') = the correlation matrix from the fit result
59# Z = requested significance 'Z sigma band'
60#
61# The linear method is fast (required 2*N evaluations of the curve, N is the number of parameters),
62# but may not be accurate in the presence of strong correlations (~>0.9) and at Z>2 due to linear and
64#
65model.plotOn(frame, ROOT.RooFit.VisualizeError(
66 r, 1), ROOT.RooFit.FillColor(ROOT.kOrange))
67
68# Calculate error using sampling method and visualize as dashed red line.
69#
70# In self method a number of curves is calculated with variations of the parameter values, sampled
71# from a multi-variate Gaussian p.d.f. that is constructed from the fit results covariance matrix.
72# The error(x) is determined by calculating a central interval that capture N% of the variations
73# for each valye of x, N% is controlled by Z (i.e. Z=1 gives N=68%). The number of sampling curves
74# is chosen to be such that at least 100 curves are expected to be outside the N% interval, is minimally
75# 100 (e.g. Z=1.Ncurve=356, Z=2.Ncurve=2156)) Intervals from the sampling method can be asymmetric,
76# and may perform better in the presence of strong correlations, may take
77# (much) longer to calculate
78model.plotOn(
79 frame,
80 ROOT.RooFit.VisualizeError(
81 r,
82 1,
83 ROOT.kFALSE),
84 ROOT.RooFit.DrawOption("L"),
85 ROOT.RooFit.LineWidth(2),
86 ROOT.RooFit.LineColor(
87 ROOT.kRed))
88
89# Perform the same type of error visualization on the background component only.
90# The VisualizeError() option can generally applied to _any_ kind of
91# plot (components, asymmetries, etc..)
92model.plotOn(
93 frame, ROOT.RooFit.VisualizeError(
94 r, 1), ROOT.RooFit.FillColor(
95 ROOT.kOrange), ROOT.RooFit.Components("bkg"))
96model.plotOn(
97 frame,
98 ROOT.RooFit.VisualizeError(
99 r,
100 1,
101 ROOT.kFALSE),
102 ROOT.RooFit.DrawOption("L"),
103 ROOT.RooFit.LineWidth(2),
104 ROOT.RooFit.LineColor(
105 ROOT.kRed),
106 ROOT.RooFit.Components("bkg"),
107 ROOT.RooFit.LineStyle(
108 ROOT.kDashed))
109
110# Overlay central value
111model.plotOn(frame)
112model.plotOn(frame, ROOT.RooFit.Components("bkg"),
113 ROOT.RooFit.LineStyle(ROOT.kDashed))
114d.plotOn(frame)
115frame.SetMinimum(0)
116
117# Visualize partial fit error
118# ------------------------------------------------------
119
120# Make plot frame
121frame2 = x.frame(ROOT.RooFit.Bins(40), ROOT.RooFit.Title(
122 "Visualization of 2-sigma partial error from (m,m2)"))
123
124# Visualize partial error. For partial error visualization the covariance matrix is first reduced as follows
125# ___ -1
126# Vred = V22 = V11 - V12 * V22 * V21
127#
128# Where V11,V12,V21, represent a block decomposition of the covariance matrix into observables that
129# are propagated (labeled by index '1') and that are not propagated (labeled by index '2'), V22bar
130# is the Shur complement of V22, as shown above
131#
132# (Note that Vred is _not_ a simple sub-matrix of V)
133
134# Propagate partial error due to shape parameters (m,m2) using linear and
135# sampling method
136model.plotOn(frame2, ROOT.RooFit.VisualizeError(
137 r, ROOT.RooArgSet(m, m2), 2), ROOT.RooFit.FillColor(ROOT.kCyan))
138model.plotOn(frame2, ROOT.RooFit.Components("bkg"), ROOT.RooFit.VisualizeError(
139 r, ROOT.RooArgSet(m, m2), 2), ROOT.RooFit.FillColor(ROOT.kCyan))
140
141model.plotOn(frame2)
142model.plotOn(frame2, ROOT.RooFit.Components("bkg"),
143 ROOT.RooFit.LineStyle(ROOT.kDashed))
144frame2.SetMinimum(0)
145
146# Make plot frame
147frame3 = x.frame(ROOT.RooFit.Bins(40), ROOT.RooFit.Title(
148 "Visualization of 2-sigma partial error from (s,s2)"))
149
150# Propagate partial error due to yield parameter using linear and sampling
151# method
152model.plotOn(frame3, ROOT.RooFit.VisualizeError(
153 r, ROOT.RooArgSet(s, s2), 2), ROOT.RooFit.FillColor(ROOT.kGreen))
154model.plotOn(frame3, ROOT.RooFit.Components("bkg"), ROOT.RooFit.VisualizeError(
155 r, ROOT.RooArgSet(s, s2), 2), ROOT.RooFit.FillColor(ROOT.kGreen))
156
157model.plotOn(frame3)
158model.plotOn(frame3, ROOT.RooFit.Components("bkg"),
159 ROOT.RooFit.LineStyle(ROOT.kDashed))
160frame3.SetMinimum(0)
161
162# Make plot frame
163frame4 = x.frame(ROOT.RooFit.Bins(40), ROOT.RooFit.Title(
164 "Visualization of 2-sigma partial error from fsig"))
165
166# Propagate partial error due to yield parameter using linear and sampling
167# method
168model.plotOn(frame4, ROOT.RooFit.VisualizeError(
169 r, ROOT.RooArgSet(fsig), 2), ROOT.RooFit.FillColor(ROOT.kMagenta))
170model.plotOn(frame4, ROOT.RooFit.Components("bkg"), ROOT.RooFit.VisualizeError(
171 r, ROOT.RooArgSet(fsig), 2), ROOT.RooFit.FillColor(ROOT.kMagenta))
172
173model.plotOn(frame4)
174model.plotOn(frame4, ROOT.RooFit.Components("bkg"),
175 ROOT.RooFit.LineStyle(ROOT.kDashed))
176frame4.SetMinimum(0)
177
178c = ROOT.TCanvas("rf610_visualerror", "rf610_visualerror", 800, 800)
179c.Divide(2, 2)
180c.cd(1)
182frame.GetYaxis().SetTitleOffset(1.4)
183frame.Draw()
184c.cd(2)
186frame2.GetYaxis().SetTitleOffset(1.6)
187frame2.Draw()
188c.cd(3)
190frame3.GetYaxis().SetTitleOffset(1.6)
191frame3.Draw()
192c.cd(4)
194frame4.GetYaxis().SetTitleOffset(1.6)
195frame4.Draw()
196
197c.SaveAs("rf610_visualerror.png")