Validation and MC studies: RooMCStudy - Using the randomizer and profile likelihood add-on models
Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Tuesday, March 19, 2024 at 07:17 PM.
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooChebychev.h"
#include "RooAddPdf.h"
#include "RooMCStudy.h"
#include "RooRandomizeParamMCSModule.h"
#include "RooDLLSignificanceMCSModule.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "TH1.h"
#include "TDirectory.h"
using namespace RooFit;
Simulation of signal and background of top quark decaying into 3 jets with background
Observable
RooRealVar mjjj("mjjj", "m(3jet) (GeV)", 100, 85., 350.);
Signal component (Gaussian)
RooRealVar mtop("mtop", "m(top)", 162);
RooRealVar wtop("wtop", "m(top) resolution", 15.2);
RooGaussian sig("sig", "top signal", mjjj, mtop, wtop);
[#0] WARNING:InputArguments -- The parameter 'wtop' with range [-inf, inf] of the RooGaussian 'sig' exceeds the safe range of (0, inf). Advise to limit its range.
Background component (Chebychev)
RooRealVar c0("c0", "Chebychev coefficient 0", -0.846, -1., 1.);
RooRealVar c1("c1", "Chebychev coefficient 1", 0.112, -1., 1.);
RooRealVar c2("c2", "Chebychev coefficient 2", 0.076, -1., 1.);
RooChebychev bkg("bkg", "combinatorial background", mjjj, RooArgList(c0, c1, c2));
Composite model
RooRealVar nsig("nsig", "number of signal events", 53, 0, 1e3);
RooRealVar nbkg("nbkg", "number of background events", 103, 0, 5e3);
RooAddPdf model("model", "model", RooArgList(sig, bkg), RooArgList(nsig, nbkg));
Configure manager to perform binned extended likelihood fits (Binned(),Extended()) on data generated with a Poisson fluctuation on Nobs (Extended())
RooMCStudy *mcs = new RooMCStudy(model, mjjj, Binned(), Silence(), Extended(true),
FitOptions(Extended(true), PrintEvalErrors(-1)));
Add module that randomizes the summed value of nsig+nbkg sampling from a uniform distribution between 0 and 1000
In general one can randomize a single parameter, or a sum of N parameters, using either a uniform or a Gaussian distribution. Multiple randomization can be executed by a single randomizer module
RooRandomizeParamMCSModule randModule;
randModule.sampleSumUniform(RooArgSet(nsig, nbkg), 50, 500);
mcs->addModule(randModule);
Add profile likelihood calculation of significance. Redo each fit while keeping parameter nsig fixed to zero. For each toy, the difference in -log(L) of both fits is stored, as well a simple significance interpretation of the delta(-logL) using Dnll = 0.5 sigma^2
RooDLLSignificanceMCSModule sigModule(nsig, 0);
mcs->addModule(sigModule);
Run 1000 experiments. This configuration will generate a fair number of (harmless) MINUIT warnings due to the instability of the Chebychev polynomial fit at low statistics.
mcs->generateAndFit(500);
[#0] PROGRESS:Generation -- RooMCStudy::run: sample 495 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 490 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 485 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 480 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 475 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 470 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 465 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 460 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 455 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 450 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 445 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 440 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 435 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 430 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 425 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 420 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 415 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 410 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 405 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 400 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 395 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 390 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 385 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 380 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 375 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 370 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 365 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 360 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 355 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 350 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 345 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 340 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 335 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 330 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 325 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 320 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 315 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 310 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 305 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 300 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 295 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 290 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 285 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 280 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 275 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 270 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 265 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 260 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 255 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 250 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 245 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 240 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 235 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 230 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 225 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 220 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 215 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 210 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 205 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 200 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 195 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 190 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 185 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 180 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 175 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 170 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 165 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 160 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 155 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 150 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 145 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 140 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 135 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 130 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 125 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 120 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 115 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 110 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 105 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 100 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 95 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 90 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 85 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 80 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 75 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 70 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 65 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 60 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 55 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 50 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 45 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 40 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 35 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 30 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 25 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 20 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 15 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 10 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 5 [#0] PROGRESS:Generation -- RooMCStudy::run: sample 0 [#0] WARNING:Generation -- Fit parameter 'mtop' does not have an error. A pull distribution cannot be generated. This might be caused by the parameter being constant or because the fits were not run. [#0] WARNING:Generation -- Fit parameter 'wtop' does not have an error. A pull distribution cannot be generated. This might be caused by the parameter being constant or because the fits were not run.
Make some plots
TH1 *dll_vs_ngen = mcs->fitParDataSet().createHistogram("ngen,dll_nullhypo_nsig", AutoBinning(40), AutoBinning(40));
TH1 *z_vs_ngen = mcs->fitParDataSet().createHistogram("ngen,significance_nullhypo_nsig", AutoBinning(40), AutoBinning(40));
TH1 *errnsig_vs_ngen = mcs->fitParDataSet().createHistogram("ngen,nsigerr", AutoBinning(40), AutoBinning(40));
TH1 *errnsig_vs_nsig = mcs->fitParDataSet().createHistogram("nsig,nsigerr", AutoBinning(40), AutoBinning(40));
Draw plots on canvas
TCanvas *c = new TCanvas("rf803_mcstudy_addons2", "rf802_mcstudy_addons2", 800, 800);
c->Divide(2, 2);
c->cd(1);
gPad->SetLeftMargin(0.15);
dll_vs_ngen->GetYaxis()->SetTitleOffset(1.6);
dll_vs_ngen->Draw("box");
c->cd(2);
gPad->SetLeftMargin(0.15);
z_vs_ngen->GetYaxis()->SetTitleOffset(1.6);
z_vs_ngen->Draw("box");
c->cd(3);
gPad->SetLeftMargin(0.15);
errnsig_vs_ngen->GetYaxis()->SetTitleOffset(1.6);
errnsig_vs_ngen->Draw("box");
c->cd(4);
gPad->SetLeftMargin(0.15);
errnsig_vs_nsig->GetYaxis()->SetTitleOffset(1.6);
errnsig_vs_nsig->Draw("box");
Make RooMCStudy object available on command line after macro finishes
gDirectory->Add(mcs);
Draw all canvases
%jsroot on
gROOT->GetListOfCanvases()->Draw()