rf616_morphing.py File Reference

Detailed Description

Use Morphing in RooFit.

This tutorial shows how to use template morphing inside RooFit. As input we have several Gaussian distributions. The output is one gaussian, with a specific mean value. Since likelihoods are often used within the framework of morphing, we provide a way to estimate the negative log likelihood (nll).

Based on example of Kyle Cranmer https://gist.github.com/cranmer/46fff8d22015e5a26619.

 
 
import ROOT
 
# Number of samples to fill the histograms
n_samples = 1000
 
 
# Kills warning messages
ROOT.RooMsgService.instance().setGlobalKillBelow(ROOT.RooFit.WARNING)
 
 
# morphing as a baseline
def morphing(setting):
    # set up a frame for plotting
    frame1 = x_var.frame()
 
    # define binning for morphing
    bin_mu_x = ROOT.RooBinning(4, 0.0, 4.0)
    grid = ROOT.RooMomentMorphFuncND.Grid(bin_mu_x)
    x_var.setBins(50)
 
    # number of 'sampled' Gaussians, if you change it, adjust the binning properly
    for i in range(5):
        # Create the sampled Gaussian
        workspace.factory(f"Gaussian::g{i}(x, mu{i}[{i}], {sigma})".format(i=i))
 
        # Fill the histograms
        hist = workspace[f"g{i}"].generateBinned([x_var], n_samples * 100)
        # Make sure that every bin is filled and we don't get zero probability
        for i_bin in range(hist.numEntries()):
            hist.add(hist.get(i_bin), 1.0)
 
        # Add the pdf to the workspace, the inOrder of 1 is necessary for calculation of the nll
        # Adjust it to 0 to see binning
        workspace.Import(ROOT.RooHistPdf(f"histpdf{i}", f"histpdf{i}", [x_var], hist, intOrder=1))
 
        # Add the pdf to the grid and to the plot
        grid.addPdf(workspace[f"histpdf{i}"], int(i))
        workspace[f"histpdf{i}"].plotOn(frame1)
 
    # Create the morphing and add it to the workspace
    morph_func = ROOT.RooMomentMorphFuncND("morph_func", "morph_func", [mu_var], [x_var], grid, setting)
 
    # Normalizes the morphed object to be a pdf, set it false to prevent warning messages and gain computational speed up
    morph_func.setPdfMode()
 
    # Creating the morphed pdf
    morph = ROOT.RooWrapperPdf("morph", "morph", morph_func, True)
    workspace.Import(morph)
    workspace["morph"].plotOn(frame1, LineColor="r")
 
    return frame1
 
 
# Define the "observed" data in a workspade
def build_ws(mu_observed, sigma):
    ws = ROOT.RooWorkspace()
    ws.factory(f"Gaussian::gauss(x[-5,15], mu[{mu_observed},0,4], {sigma})".format(mu_observed=mu_observed))
 
    return ws
 
 
# The "observed" data
mu_observed = 2.5
sigma = 1.5
workspace = build_ws(mu_observed, sigma)
x_var = workspace["x"]
mu_var = workspace["mu"]
gauss = workspace.pdf("gauss")
obs_data = gauss.generate(x_var, n_samples)
 
 
# Create the exact negative log likelihood functions for Gaussian model
nll_gauss = gauss.createNLL(obs_data)
 
# Compute the morphed nll
frame1 = morphing(ROOT.RooMomentMorphFuncND.Linear)
 
# TODO: Fix RooAddPdf::fixCoefNormalization(nset) warnings with new CPU backend
nll_morph = workspace["morph"].createNLL(obs_data, EvalBackend="legacy")
 
# Plot the negative logarithmic summed likelihood
frame2 = mu_var.frame(Title="Negative log Likelihood")
nll_gauss.plotOn(frame2, LineColor="b", ShiftToZero=True, Name="gauss")
nll_morph.plotOn(frame2, LineColor="r", ShiftToZero=True, Name="morphed")
 
c = ROOT.TCanvas("rf616_morphing", "rf616_morphing", 800, 400)
c.Divide(2)
c.cd(1)
ROOT.gPad.SetLeftMargin(0.15)
frame1.GetYaxis().SetTitleOffset(1.8)
frame1.Draw()
c.cd(2)
ROOT.gPad.SetLeftMargin(0.15)
frame2.GetYaxis().SetTitleOffset(1.8)
frame2.Draw()
 
 
# Compute the minimum via minuit and display the results
for nll in [nll_gauss, nll_morph]:
    minimizer = ROOT.RooMinimizer(nll)
    minimizer.setPrintLevel(-1)
    minimizer.minimize("Minuit2")
    result = minimizer.save()
    result.Print()

[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf0_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf1_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf2_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf3_Int[x] is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name histpdf4_Int[x] is already in this set
 
  RooFitResult: minimized FCN value: 1862.97, estimated distance to minimum: 2.32702e-05
                covariance matrix quality: Full, accurate covariance matrix
                Status : MINIMIZE=0 
 
    Floating Parameter    FinalValue +/-  Error   
  --------------------  --------------------------
                    mu    2.5399e+00 +/-  4.74e-02
 
 
  RooFitResult: minimized FCN value: 1862.82, estimated distance to minimum: 9.31471e-06
                covariance matrix quality: Full, accurate covariance matrix
                Status : MINIMIZE=0 
 
    Floating Parameter    FinalValue +/-  Error   
  --------------------  --------------------------
                    mu    2.5372e+00 +/-  4.08e-02
 

Date

August 2024

Author

Robin Syring

Definition in file rf616_morphing.py.