Logo ROOT  
Reference Guide
rf707_kernelestimation.py
Go to the documentation of this file.
1 ## \file
2 ## \ingroup tutorial_roofit
3 ## \notebook
4 ## Special pdf's: using non-parametric (multi-dimensional) kernel estimation pdfs
5 ##
6 ## \macro_code
7 ##
8 ## \date February 2018
9 ## \authors Clemens Lange, Wouter Verkerke (C++ version)
10 
11 import ROOT
12 
13 
14 # Create low stats 1D dataset
15 # -------------------------------------------------------
16 
17 # Create a toy pdf for sampling
18 x = ROOT.RooRealVar("x", "x", 0, 20)
19 p = ROOT.RooPolynomial("p", "p", x, ROOT.RooArgList(ROOT.RooFit.RooConst(
20  0.01), ROOT.RooFit.RooConst(-0.01), ROOT.RooFit.RooConst(0.0004)))
21 
22 # Sample 500 events from p
23 data1 = p.generate(ROOT.RooArgSet(x), 200)
24 
25 # Create 1D kernel estimation pdf
26 # ---------------------------------------------------------------
27 
28 # Create adaptive kernel estimation pdf. In self configuration the input data
29 # is mirrored over the boundaries to minimize edge effects in distribution
30 # that do not fall to zero towards the edges
31 kest1 = ROOT.RooKeysPdf("kest1", "kest1", x, data1,
32  ROOT.RooKeysPdf.MirrorBoth)
33 
34 # An adaptive kernel estimation pdf on the same data without mirroring option
35 # for comparison
36 kest2 = ROOT.RooKeysPdf("kest2", "kest2", x, data1,
37  ROOT.RooKeysPdf.NoMirror)
38 
39 # Adaptive kernel estimation pdf with increased bandwidth scale factor
40 # (promotes smoothness over detail preservation)
41 kest3 = ROOT.RooKeysPdf("kest1", "kest1", x, data1,
42  ROOT.RooKeysPdf.MirrorBoth, 2)
43 
44 # Plot kernel estimation pdfs with and without mirroring over data
45 frame = x.frame(
46  ROOT.RooFit.Title("Adaptive kernel estimation pdf with and w/o mirroring"),
47  ROOT.RooFit.Bins(20))
48 data1.plotOn(frame)
49 kest1.plotOn(frame)
50 kest2.plotOn(frame, ROOT.RooFit.LineStyle(
51  ROOT.kDashed), ROOT.RooFit.LineColor(ROOT.kRed))
52 
53 # Plot kernel estimation pdfs with regular and increased bandwidth
54 frame2 = x.frame(ROOT.RooFit.Title(
55  "Adaptive kernel estimation pdf with regular, bandwidth"))
56 kest1.plotOn(frame2)
57 kest3.plotOn(frame2, ROOT.RooFit.LineColor(ROOT.kMagenta))
58 
59 # Create low status 2D dataset
60 # -------------------------------------------------------
61 
62 # Construct a 2D toy pdf for sampleing
63 y = ROOT.RooRealVar("y", "y", 0, 20)
64 py = ROOT.RooPolynomial("py", "py", y, ROOT.RooArgList(ROOT.RooFit.RooConst(
65  0.01), ROOT.RooFit.RooConst(0.01), ROOT.RooFit.RooConst(-0.0004)))
66 pxy = ROOT.RooProdPdf("pxy", "pxy", ROOT.RooArgList(p, py))
67 data2 = pxy.generate(ROOT.RooArgSet(x, y), 1000)
68 
69 # Create 2D kernel estimation pdf
70 # ---------------------------------------------------------------
71 
72 # Create 2D adaptive kernel estimation pdf with mirroring
73 kest4 = ROOT.RooNDKeysPdf("kest4", "kest4", ROOT.RooArgList(x, y), data2, "am")
74 
75 # Create 2D adaptive kernel estimation pdf with mirroring and double
76 # bandwidth
77 kest5 = ROOT.RooNDKeysPdf(
78  "kest5", "kest5", ROOT.RooArgList(
79  x, y), data2, "am", 2)
80 
81 # Create a histogram of the data
82 hh_data = ROOT.RooAbsData.createHistogram(
83  data2, "hh_data", x, ROOT.RooFit.Binning(10), ROOT.RooFit.YVar(
84  y, ROOT.RooFit.Binning(10)))
85 
86 # Create histogram of the 2d kernel estimation pdfs
87 hh_pdf = kest4.createHistogram("hh_pdf", x, ROOT.RooFit.Binning(
88  25), ROOT.RooFit.YVar(y, ROOT.RooFit.Binning(25)))
89 hh_pdf2 = kest5.createHistogram("hh_pdf2", x, ROOT.RooFit.Binning(
90  25), ROOT.RooFit.YVar(y, ROOT.RooFit.Binning(25)))
91 hh_pdf.SetLineColor(ROOT.kBlue)
92 hh_pdf2.SetLineColor(ROOT.kMagenta)
93 
94 c = ROOT.TCanvas("rf707_kernelestimation",
95  "rf707_kernelestimation", 800, 800)
96 c.Divide(2, 2)
97 c.cd(1)
98 ROOT.gPad.SetLeftMargin(0.15)
99 frame.GetYaxis().SetTitleOffset(1.4)
100 frame.Draw()
101 c.cd(2)
102 ROOT.gPad.SetLeftMargin(0.15)
103 frame2.GetYaxis().SetTitleOffset(1.8)
104 frame2.Draw()
105 c.cd(3)
106 ROOT.gPad.SetLeftMargin(0.15)
107 hh_data.GetZaxis().SetTitleOffset(1.4)
108 hh_data.Draw("lego")
109 c.cd(4)
110 ROOT.gPad.SetLeftMargin(0.20)
111 hh_pdf.GetZaxis().SetTitleOffset(2.4)
112 hh_pdf.Draw("surf")
113 hh_pdf2.Draw("surfsame")
114 
115 c.SaveAs("rf707_kernelestimation.png")