rf705_linearmorph

'SPECIAL PDFS' RooFit tutorial macro #705

Linear interpolation between p.d.f shapes using the 'Alex Read' algorithm

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 Wednesday, November 30, 2022 at 11:24 AM.

In [1]:
import ROOT
Welcome to JupyROOT 6.27/01

Create end point pdf shapes

Observable

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

Lower end point shape: a Gaussian

In [3]:
g1mean = ROOT.RooRealVar("g1mean", "g1mean", -10)
g1 = ROOT.RooGaussian("g1", "g1", x, g1mean, ROOT.RooFit.RooConst(2))

Upper end point shape: a Polynomial

In [4]:
g2 = ROOT.RooPolynomial("g2", "g2", x, [-0.03, -0.001])

Create interpolating pdf

Create interpolation variable

In [5]:
alpha = ROOT.RooRealVar("alpha", "alpha", 0, 1.0)

Specify sampling density on observable and interpolation variable

In [6]:
x.setBins(1000, "cache")
alpha.setBins(50, "cache")

Construct interpolating pdf in (x,a) represent g1(x) at a=a_min and g2(x) at a=a_max

In [7]:
lmorph = ROOT.RooIntegralMorph("lmorph", "lmorph", g1, g2, x, alpha)

Plot interpolating pdf aat various alphas a l p h a

Show end points as blue curves

In [8]:
frame1 = x.frame()
g1.plotOn(frame1)
g2.plotOn(frame1)
Out[8]:
<cppyy.gbl.RooPlot object at 0x891dd70>

Show interpolated shapes in red

In [9]:
alpha.setVal(0.125)
lmorph.plotOn(frame1, LineColor="r")
alpha.setVal(0.25)
lmorph.plotOn(frame1, LineColor="r")
alpha.setVal(0.375)
lmorph.plotOn(frame1, LineColor="r")
alpha.setVal(0.50)
lmorph.plotOn(frame1, LineColor="r")
alpha.setVal(0.625)
lmorph.plotOn(frame1, LineColor="r")
alpha.setVal(0.75)
lmorph.plotOn(frame1, LineColor="r")
alpha.setVal(0.875)
lmorph.plotOn(frame1, LineColor="r")
alpha.setVal(0.95)
lmorph.plotOn(frame1, LineColor="r")
Out[9]:
<cppyy.gbl.RooPlot object at 0x891dd70>
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8a68aa0 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8a68aa0 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8a68aa0 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8a68aa0 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8aa65c0 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8aa65c0 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8aa65c0 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8aa65c0 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0

Show 2D distribution of pdf(x,alpha)

Create 2D histogram

In [10]:
hh = lmorph.createHistogram("hh", x, Binning=40, YVar=dict(var=alpha, Binning=40))
hh.SetLineColor(ROOT.kBlue)
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8b4c930 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x_alpha for nset (x,alpha) with code 0 from preexisting content.

Fit pdf to dataset with alpha=0.8

Generate a toy dataset alpha = 0.8

In [11]:
alpha.setVal(0.8)
data = lmorph.generate({x}, 1000)
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8b4c930 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8ccf7f0 with pdf g1_MORPH_g2_CACHE_Obs[x]_NORM_x for nset (x) with code 0 from preexisting content.

Fit pdf to toy data

In [12]:
lmorph.setCacheAlpha(True)
lmorph.fitTo(data, Verbose=True)
Out[12]:
<cppyy.gbl.RooFitResult object at 0x(nil)>
[#0] PROGRESS:Eval -- RooIntegralMorph::fillCacheObject(lmorph) filling multi-dimensional cache..................................................

[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x8ff84e0 with pdf g1_MORPH_g2_CACHE_Obs[alpha,x]_NORM_x for nset (x) with code 0
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
[#1] INFO:Minimization --  The following expressions have been identified as constant and will be precalculated and cached: (g1,g2)
[#0] WARNING:Minimization -- RooAbsMinimizerFcn::synchronize: WARNING: no initial error estimate available for alpha: using 0.1
 **********
 **    1 **SET PRINT           1
 **********
 **********
 **    2 **SET NOGRAD
 **********
 PARAMETER DEFINITIONS:
    NO.   NAME         VALUE      STEP SIZE      LIMITS
     1 alpha        8.00000e-01  1.00000e-01    0.00000e+00  1.00000e+00
 **********
 **    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         500           1
 **********
 FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4.

prevFCN = 2639.749221   START MIGRAD MINIMIZATION.  STRATEGY  1.  CONVERGENCE WHEN EDM .LT. 1.00e-03
alpha=0.8026, 
prevFCN = 2639.409749  alpha=0.7974, 
prevFCN = 2640.17258  alpha=0.8003, 
prevFCN = 2639.703144  alpha=0.7997, 
prevFCN = 2639.796557   FCN=2639.75 FROM MIGRAD    STATUS=INITIATE        4 CALLS           5 TOTAL
                     EDM= unknown      STRATEGY= 1      NO ERROR MATRIX       
  EXT PARAMETER               CURRENT GUESS       STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  alpha        8.00000e-01   1.00000e-01   2.57889e-01  -5.89245e+01
                               ERR DEF= 0.5
alpha=0.8116, 
prevFCN = 2638.919064  alpha=0.8119, 
prevFCN = 2638.930985  alpha=0.8113, 
prevFCN = 2638.908768  alpha=0.8094, 
prevFCN = 2638.901287  alpha=0.8102, 
prevFCN = 2638.885143  alpha=0.8103, 
prevFCN = 2638.886672  alpha=0.8101, 
prevFCN = 2638.883909  alpha=0.8095, 
prevFCN = 2638.898442  alpha=0.8098, 
prevFCN = 2638.889204  alpha=0.8099, 
prevFCN = 2638.884962  alpha=0.8104, 
prevFCN = 2638.887231  alpha=0.8098, 
prevFCN = 2638.888592  alpha=0.8102, 
prevFCN = 2638.885076  alpha=0.81, 
prevFCN = 2638.883889  alpha=0.8098, 
prevFCN = 2638.888745  alpha=0.81, 
prevFCN = 2638.883589  alpha=0.81, 
prevFCN = 2638.883345  alpha=0.8095, 
prevFCN = 2638.898453  alpha=0.8098, 
prevFCN = 2638.89011  alpha=0.8099, 
prevFCN = 2638.886243  alpha=0.81, 
prevFCN = 2638.884385  alpha=0.81, 
prevFCN = 2638.883475  alpha=0.8101, 
prevFCN = 2638.884457  alpha=0.8099, 
prevFCN = 2638.885715  alpha=0.8101, 
prevFCN = 2638.883988  alpha=0.81, 
prevFCN = 2638.884446  alpha=0.81, 
prevFCN = 2638.883552  alpha=0.81, 
prevFCN = 2638.883372  alpha=0.81, 
prevFCN = 2638.883349  alpha=0.81, 
prevFCN = 2638.883345   MIGRAD FAILS TO FIND IMPROVEMENT
 MIGRAD MINIMIZATION HAS CONVERGED.
 MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
alpha=0.81, 
prevFCN = 2638.883345  alpha=0.8101, 
prevFCN = 2638.883988  alpha=0.81, 
prevFCN = 2638.884446  alpha=0.81, 
prevFCN = 2638.883467  alpha=0.81, 
prevFCN = 2638.883226  alpha=0.81, 
prevFCN = 2638.883369  alpha=0.81, 
prevFCN = 2638.883321   COVARIANCE MATRIX CALCULATED SUCCESSFULLY
 FCN=2638.88 FROM MIGRAD    STATUS=CONVERGED      41 CALLS          42 TOTAL
                     EDM=0.000209046    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                   STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  alpha        8.10021e-01   1.63457e-03   1.74006e-04   3.46989e+00
                               ERR DEF= 0.5
 EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  1    ERR DEF=0.5
  2.672e-06 
alpha=0.81,  **********
 **    7 **SET ERR         0.5
 **********
 **********
 **    8 **SET PRINT           1
 **********
 **********
 **    9 **HESSE         500
 **********

prevFCN = 2638.883345  alpha=0.81, 
prevFCN = 2638.883369  alpha=0.81, 
prevFCN = 2638.883321  alpha=0.8103, 
prevFCN = 2638.886457  alpha=0.8097, 
prevFCN = 2638.89054  alpha=0.8131, 
prevFCN = 2638.992315  alpha=0.8069, 
prevFCN = 2639.023618  alpha=0.8171, 
prevFCN = 2639.367946  alpha=0.8029, 
prevFCN = 2639.374671  alpha=0.8102, 
prevFCN = 2638.885206  alpha=0.8098, 
prevFCN = 2638.887615  alpha=0.8101, 
prevFCN = 2638.883674  alpha=0.81, 
prevFCN = 2638.883555   COVARIANCE MATRIX CALCULATED SUCCESSFULLY
 FCN=2638.88 FROM HESSE     STATUS=OK             13 CALLS          55 TOTAL
                     EDM=0.00230644    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                INTERNAL      INTERNAL  
  NO.   NAME      VALUE            ERROR       STEP SIZE       VALUE   
   1  alpha        8.10021e-01   7.16762e-03   1.80513e-02   6.68797e-01
                               ERR DEF= 0.5
 EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  1    ERR DEF=0.5
  5.138e-05 
alpha=0.81, [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization

Plot fitted pdf and data overlaid

In [13]:
frame2 = x.frame(Bins=100)
data.plotOn(frame2)
lmorph.plotOn(frame2)
Out[13]:
<cppyy.gbl.RooPlot object at 0x8fa2540>
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x91af9f0 with pdf g1_MORPH_g2_CACHE_Obs[alpha,x]_NORM_x for nset (x) with code 0 from preexisting content.

Scan -log(L) vs alpha

Show scan -log(L) of dataset w.r.t alpha

In [14]:
frame3 = alpha.frame(Bins=100, Range=(0.1, 0.9))

Make 2D pdf of histogram

In [15]:
nll = lmorph.createNLL(data)
nll.plotOn(frame3, ShiftToZero=True)

lmorph.setCacheAlpha(False)

c = ROOT.TCanvas("rf705_linearmorph", "rf705_linearmorph", 800, 800)
c.Divide(2, 2)
c.cd(1)
ROOT.gPad.SetLeftMargin(0.15)
frame1.GetYaxis().SetTitleOffset(1.6)
frame1.Draw()
c.cd(2)
ROOT.gPad.SetLeftMargin(0.20)
hh.GetZaxis().SetTitleOffset(2.5)
hh.Draw("surf")
c.cd(3)
ROOT.gPad.SetLeftMargin(0.15)
frame3.GetYaxis().SetTitleOffset(1.4)
frame3.Draw()
c.cd(4)
ROOT.gPad.SetLeftMargin(0.15)
frame2.GetYaxis().SetTitleOffset(1.4)
frame2.Draw()

c.SaveAs("rf705_linearmorph.png")
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x9198a40 with pdf g1_MORPH_g2_CACHE_Obs[alpha,x]_NORM_x for nset (x) with code 0 from preexisting content.
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x95bf880 with pdf g1_MORPH_g2_CACHE_Obs[alpha,x]_NORM_x for nset (x) with code 0 from preexisting content.
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(lmorph) creating new cache 0x95bf880 with pdf g1_MORPH_g2_CACHE_Obs[alpha,x]_NORM_x for nset (x) with code 0 from preexisting content.
Info in <TCanvas::Print>: png file rf705_linearmorph.png has been created

Draw all canvases

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