Likelihood and minimization: understanding and customizing error handling in likelihood evaluations
Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, April 17, 2024 at 11:19 AM.
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooArgusBG.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
using namespace RooFit;
Observable
RooRealVar m("m", "m", 5.20, 5.30);
Parameters
RooRealVar m0("m0", "m0", 5.291, 5.20, 5.30);
RooRealVar k("k", "k", -30, -50, -10);
RooArgusBG argus("argus", "argus", m, m0, k);
Sample 1000 events in m from argus
std::unique_ptr<RooDataSet> data{argus.generate(m, 1000)};
input_line_51:2:2: warning: 'data' shadows a declaration with the same name in the 'std' namespace; use '::data' to reference this declaration std::unique_ptr<RooDataSet> data{argus.generate(m, 1000)}; ^
RooPlot *frame1 = m.frame(Bins(40), Title("Argus model and data"));
data->plotOn(frame1);
argus.plotOn(frame1);
input_line_52:3:1: error: reference to 'data' is ambiguous data->plotOn(frame1); ^ input_line_51:2:30: note: candidate found by name lookup is 'data' std::unique_ptr<RooDataSet> data{argus.generate(m, 1000)}; ^ /usr/include/c++/9/bits/range_access.h:318:5: note: candidate found by name lookup is 'std::data' data(initializer_list<_Tp> __il) noexcept ^ /usr/include/c++/9/bits/range_access.h:289:5: note: candidate found by name lookup is 'std::data' data(_Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:299:5: note: candidate found by name lookup is 'std::data' data(const _Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:309:5: note: candidate found by name lookup is 'std::data' data(_Tp (&__array)[_Nm]) noexcept ^
The ARGUS background shape has a sharp kinematic cutoff at m=m0 and is prone to evaluation errors if the cutoff parameter m0 is floated: when the pdf cutoff value is lower than that in data events with m>m0 will have zero probability
Perform unbinned ML fit. Print detailed error messages for up to 10 events per likelihood evaluation. The default error handling strategy is to return a very high value of the likelihood to MINUIT if errors occur, which will force MINUIT to retreat from the problematic area
argus.fitTo(*data, PrintEvalErrors(10));
input_line_53:2:15: error: reference to 'data' is ambiguous argus.fitTo(*data, PrintEvalErrors(10)); ^ input_line_51:2:30: note: candidate found by name lookup is 'data' std::unique_ptr<RooDataSet> data{argus.generate(m, 1000)}; ^ /usr/include/c++/9/bits/range_access.h:318:5: note: candidate found by name lookup is 'std::data' data(initializer_list<_Tp> __il) noexcept ^ /usr/include/c++/9/bits/range_access.h:289:5: note: candidate found by name lookup is 'std::data' data(_Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:299:5: note: candidate found by name lookup is 'std::data' data(const _Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:309:5: note: candidate found by name lookup is 'std::data' data(_Tp (&__array)[_Nm]) noexcept ^
Perform another fit. In this configuration only the number of errors per likelihood evaluation is shown, if it is greater than zero. The EvalErrorWall(false) arguments disables the default error handling strategy and will cause the actual (problematic) value of the likelihood to be passed to MINUIT.
NB: Use of this option is NOT recommended as default strategy as broken -log(L) values can often be lower than 'good' ones because offending events are removed. This may effectively create a false minimum in problem areas. This is clearly illustrated in the second plot
m0.setError(0.1);
argus.fitTo(*data, PrintEvalErrors(0), EvalErrorWall(false));
input_line_54:3:14: error: reference to 'data' is ambiguous argus.fitTo(*data, PrintEvalErrors(0), EvalErrorWall(false)); ^ input_line_51:2:30: note: candidate found by name lookup is 'data' std::unique_ptr<RooDataSet> data{argus.generate(m, 1000)}; ^ /usr/include/c++/9/bits/range_access.h:318:5: note: candidate found by name lookup is 'std::data' data(initializer_list<_Tp> __il) noexcept ^ /usr/include/c++/9/bits/range_access.h:289:5: note: candidate found by name lookup is 'std::data' data(_Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:299:5: note: candidate found by name lookup is 'std::data' data(const _Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:309:5: note: candidate found by name lookup is 'std::data' data(_Tp (&__array)[_Nm]) noexcept ^
Construct likelihood function of model and data
std::unique_ptr<RooAbsReal> nll{argus.createNLL(*data)};
input_line_55:2:51: error: reference to 'data' is ambiguous std::unique_ptr<RooAbsReal> nll{argus.createNLL(*data)}; ^ input_line_51:2:30: note: candidate found by name lookup is 'data' std::unique_ptr<RooDataSet> data{argus.generate(m, 1000)}; ^ /usr/include/c++/9/bits/range_access.h:318:5: note: candidate found by name lookup is 'std::data' data(initializer_list<_Tp> __il) noexcept ^ /usr/include/c++/9/bits/range_access.h:289:5: note: candidate found by name lookup is 'std::data' data(_Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:299:5: note: candidate found by name lookup is 'std::data' data(const _Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:309:5: note: candidate found by name lookup is 'std::data' data(_Tp (&__array)[_Nm]) noexcept ^
Plot likelihood in m0 in range that includes problematic values In this configuration no messages are printed for likelihood evaluation errors, but if an likelihood value evaluates with error, the corresponding value on the curve will be set to the value given in EvalErrorValue().
RooPlot *frame2 = m0.frame(Range(5.288, 5.293), Title("-log(L) scan vs m0, problematic regions masked"));
nll->plotOn(frame2, PrintEvalErrors(-1), ShiftToZero(), EvalErrorValue(nll->getVal() + 10), LineColor(kRed));
frame2->SetMaximum(15);
frame2->SetMinimum(0);
TCanvas *c = new TCanvas("rf606_nllerrorhandling", "rf606_nllerrorhandling", 1200, 400);
c->Divide(2);
c->cd(1);
gPad->SetLeftMargin(0.15);
frame1->GetYaxis()->SetTitleOffset(1.4);
frame1->Draw();
c->cd(2);
gPad->SetLeftMargin(0.15);
frame2->GetYaxis()->SetTitleOffset(1.4);
frame2->Draw();
input_line_56:3:103: error: cannot take the address of an rvalue of type 'EColor' nll->plotOn(frame2, PrintEvalErrors(-1), ShiftToZero(), EvalErrorValue(nll->getVal() + 10), LineColor(kRed)); ^~~~ Error while creating dynamic expression for: nll->plotOn(frame2, PrintEvalErrors(-1), ShiftToZero(), EvalErrorValue(nll->getVal() + 10), LineColor(kRed))
Draw all canvases
%jsroot on
gROOT->GetListOfCanvases()->Draw()