rf601_intminuit

Likelihood and minimization: interactive minimization with MINUIT

Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, November 30, 2022 at 11:23 AM.

In [ ]:
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooProdPdf.h"
#include "RooAddPdf.h"
#include "RooMinimizer.h"
#include "RooFitResult.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "TH1.h"
using namespace RooFit;

Setup pdf and likelihood

Observable

In [ ]:
RooRealVar x("x", "x", -20, 20);

Model (intentional strong correlations)

In [ ]:
RooRealVar mean("mean", "mean of g1 and g2", 0);
RooRealVar sigma_g1("sigma_g1", "width of g1", 3);
RooGaussian g1("g1", "g1", x, mean, sigma_g1);

RooRealVar sigma_g2("sigma_g2", "width of g2", 4, 3.0, 6.0);
RooGaussian g2("g2", "g2", x, mean, sigma_g2);

RooRealVar frac("frac", "frac", 0.5, 0.0, 1.0);
RooAddPdf model("model", "model", RooArgList(g1, g2), frac);

Generate 1000 events

In [ ]:
RooDataSet *data = model.generate(x, 1000);

Construct unbinned likelihood of model w.r.t. data

In [ ]:
RooAbsReal *nll = model.createNLL(*data);

Interactive minimization, error analysis

Create MINUIT interface object

In [ ]:
RooMinimizer m(*nll);

Activate verbose logging of MINUIT parameter space stepping

In [ ]:
m.setVerbose(true);

Call MIGRAD to minimize the likelihood

In [ ]:
m.migrad();

Print values of all parameters, that reflect values (and error estimates) that are back propagated from MINUIT

In [ ]:
model.getParameters(x)->Print("s");

Disable verbose logging

In [ ]:
m.setVerbose(false);

Run HESSE to calculate errors from d2L/dp2

In [ ]:
m.hesse();

Print value (and error) of sigma_g2 parameter, that reflects value and error back propagated from MINUIT

In [ ]:
sigma_g2.Print();

Run MINOS on sigma_g2 parameter only

In [ ]:
m.minos(sigma_g2);

Print value (and error) of sigma_g2 parameter, that reflects value and error back propagated from MINUIT

In [ ]:
sigma_g2.Print();

Saving results, contour plots

Save a snapshot of the fit result. This object contains the initial fit parameters, the final fit parameters, the complete correlation matrix, the EDM, the minimized FCN , the last MINUIT status code and the number of times the RooFit function object has indicated evaluation problems (e.g. zero probabilities during likelihood evaluation)

In [ ]:
RooFitResult *r = m.save();

Make contour plot of mx vs sx at 1,2,3 sigma

In [ ]:
RooPlot *frame = m.contour(frac, sigma_g2, 1, 2, 3);
frame->SetTitle("Minuit contour plot");

Print the fit result snapshot

In [ ]:
r->Print("v");

Change parameter values, floating

At any moment you can manually change the value of a (constant) parameter

In [ ]:
mean = 0.3;

Rerun MIGRAD,HESSE

In [ ]:
m.migrad();
m.hesse();
frac.Print();

Now fix sigma_g2

In [ ]:
sigma_g2.setConstant(true);

Rerun MIGRAD,HESSE

In [ ]:
m.migrad();
m.hesse();
frac.Print();

new TCanvas("rf601_intminuit", "rf601_intminuit", 600, 600);
gPad->SetLeftMargin(0.15);
frame->GetYaxis()->SetTitleOffset(1.4);
frame->Draw();

Draw all canvases

In [ ]:
%jsroot on
gROOT->GetListOfCanvases()->Draw()