Logo ROOT  
Reference Guide
rf612_recoverFromInvalidParameters.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_roofit
3/// \notebook -js
4/// Likelihood and minimization: Recover from regions where the function is not defined.
5///
6/// We demonstrate improved recovery from disallowed parameters. For this, we use a polynomial PDF of the form
7/// \f[
8/// \mathrm{Pol2} = \mathcal{N} \left( c + a_1 \cdot x + a_2 \cdot x^2 + 0.01 \cdot x^3 \right),
9/// \f]
10/// where \f$ \mathcal{N} \f$ is a normalisation factor. Unless the parameters are chosen carefully,
11/// this function can be negative, and hence, it cannot be used as a PDF. In this case, RooFit passes
12/// an error to the minimiser, which might try to recover.
13///
14/// Before ROOT 6.24, RooFit always passed the highest function value that was encountered during the minimisation
15/// to the minimiser. If a parameter is far in a disallowed region, the minimiser has to blindly test various
16/// values of the parameters. It might find the correct values by chance, but often be unable to recover from bad
17/// starting values. Here, we use a model with such bad values.
18///
19/// Starting with ROOT 6.24, the minimiser receives more information. For example, when a PDF is negative,
20/// the magnitude of the "undershoot" is passed to the minimiser. The minimiser can use this to compute a
21/// gradient, which will eventually lead it out of the disallowed region. The steepness of this gradient
22/// can be chosen using RooFit::RecoverFromUndefinedRegions(double). A value of zero is equivalent to RooFit
23/// before ROOT 6.24. Positive values activate the recovery. Values between 1. and 10. were found to be a
24/// good default. If no argument is passed, RooFit uses 10.
25///
26/// \macro_output
27/// \macro_code
28///
29/// \date 10/2020
30/// \author Stephan Hageboeck
31
32#include <RooRealVar.h>
33#include <RooPolynomial.h>
34#include <RooPlot.h>
35#include <RooDataSet.h>
36#include <RooGlobalFunc.h>
37#include <RooFitResult.h>
38#include <RooMsgService.h>
39
40#include <TCanvas.h>
41#include <TLegend.h>
42
44
45 // Create a fit model:
46 // The polynomial is notoriously unstable, because it can quickly go negative.
47 // Since PDFs need to be positive, one often ends up with an unstable fit model.
48 RooRealVar x("x", "x", -15, 15);
49 RooRealVar a1("a1", "a1", -0.5, -10., 20.);
50 RooRealVar a2("a2", "a2", 0.2, -10., 20.);
51 RooRealVar a3("a3", "a3", 0.01);
52 RooPolynomial pdf("pol3", "c + a1 * x + a2 * x*x + 0.01 * x*x*x", x, RooArgSet(a1, a2, a3));
53
54 // Create toy data with all-positive coefficients:
55 std::unique_ptr<RooDataSet> data(pdf.generate(x, 10000));
56
57 // For plotting.
58 // We create pointers to the plotted objects. We want these objects to leak out of the function,
59 // so we can still see them after it returns.
60 TCanvas* c = new TCanvas();
61 RooPlot* frame = x.frame();
62 data->plotOn(frame, RooFit::Name("data"));
63
64 // Plotting a PDF with disallowed parameters doesn't work. We would get a lot of error messages.
65 // Therefore, we disable plotting messages in RooFit's message streams:
68
69
70 // RooFit before ROOT 6.24
71 // --------------------------------
72 // Before 6.24, RooFit wasn't able to recover from invalid parameters. The minimiser just errs around
73 // the starting values of the parameters without finding any improvement.
74
75 // Set up the parameters such that the PDF would come out negative. The PDF is now undefined.
76 a1.setVal(10.);
77 a2.setVal(-1.);
78
79 // Perform a fit:
80 RooFitResult* fitWithoutRecovery = pdf.fitTo(*data, RooFit::Save(),
81 RooFit::RecoverFromUndefinedRegions(0.), // This is how RooFit behaved prior to ROOT 6.24
82 RooFit::PrintEvalErrors(-1), // We are expecting a lot of evaluation errors. -1 switches off printing.
84
85 pdf.plotOn(frame, RooFit::LineColor(kRed), RooFit::Name("noRecovery"));
86
87
88
89 // RooFit since ROOT 6.24
90 // --------------------------------
91 // The minimiser gets information about the "badness" of the violation of the function definition. It uses this
92 // to find its way out of the disallowed parameter regions.
93 std::cout << "\n\n\n-------------- Starting second fit ---------------\n\n" << std::endl;
94
95 // Reset the parameters such that the PDF is again undefined.
96 a1.setVal(10.);
97 a2.setVal(-1.);
98
99 // Fit again, but pass recovery information to the minimiser:
100 RooFitResult* fitWithRecovery = pdf.fitTo(*data, RooFit::Save(),
101 RooFit::RecoverFromUndefinedRegions(1.), // The magnitude of the recovery information can be chosen here.
102 // Higher values mean more aggressive recovery.
103 RooFit::PrintEvalErrors(-1), // We are still expecting a few evaluation errors.
105
106 pdf.plotOn(frame, RooFit::LineColor(kBlue), RooFit::Name("recovery"));
107
108
109
110 // Collect results and plot.
111 // --------------------------------
112 // We print the two fit results, and plot the fitted curves.
113 // The curve of the fit without recovery cannot be plotted, because the PDF is undefined if a2 < 0.
114 fitWithoutRecovery->Print();
115 std::cout << "Without recovery, the fitter encountered " << fitWithoutRecovery->numInvalidNLL()
116 << " invalid function values. The parameters are unchanged." << std::endl;
117
118 fitWithRecovery->Print();
119 std::cout << "With recovery, the fitter encountered " << fitWithRecovery->numInvalidNLL()
120 << " invalid function values, but the parameters are fitted." << std::endl;
121
122 TLegend* legend = new TLegend(0.5, 0.7, 0.9, 0.9);
123 legend->SetBorderSize(0);
124 legend->SetFillStyle(0);
125 legend->AddEntry("data", "Data", "P");
126 legend->AddEntry("noRecovery", "Without recovery (cannot be plotted)", "L");
127 legend->AddEntry("recovery", "With recovery", "L");
128 frame->Draw();
129 legend->Draw();
130 c->Draw();
131}
132
#define c(i)
Definition: RSha256.hxx:101
@ kRed
Definition: Rtypes.h:66
@ kBlue
Definition: Rtypes.h:66
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:56
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
Definition: RooFitResult.h:40
void Print(Option_t *options=nullptr) const override
Print TNamed name and title.
Definition: RooFitResult.h:66
Int_t numInvalidNLL() const
Return number of NLL evaluations with problems.
Definition: RooFitResult.h:91
RooPlot * plotOn(RooPlot *frame, const RooAbsArg &par1, const RooAbsArg &par2, const char *options="ME") const
Add objects to a 2D plot.
Definition: RooFitResult.h:144
static RooMsgService & instance()
Return reference to singleton instance.
StreamConfig & getStream(Int_t id)
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:43
void Draw(Option_t *options=nullptr) override
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:649
RooPolynomial implements a polynomial p.d.f of the form.
Definition: RooPolynomial.h:28
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition: TAttFill.h:39
The Canvas class.
Definition: TCanvas.h:23
This class displays a legend box (TPaveText) containing several legend entries.
Definition: TLegend.h:23
TLegendEntry * AddEntry(const TObject *obj, const char *label="", Option_t *option="lpf")
Add a new entry to this legend.
Definition: TLegend.cxx:317
void Draw(Option_t *option="") override
Draw this legend with its current attributes.
Definition: TLegend.cxx:422
virtual void SetBorderSize(Int_t bordersize=4)
Definition: TPave.h:73
RooCmdArg Save(bool flag=true)
RooCmdArg PrintEvalErrors(Int_t numErrors)
RooCmdArg PrintLevel(Int_t code)
RooCmdArg RecoverFromUndefinedRegions(double strength)
When parameters are chosen such that a PDF is undefined, try to indicate to the minimiser how to leav...
RooCmdArg LineColor(Color_t color)
RooCmdArg Name(const char *name)
Double_t x[n]
Definition: legend1.C:17
void removeTopic(RooFit::MsgTopic oldTopic)