This example compares different approaches to determining parameter uncertainties in weighted unbinned maximum likelihood fits. Performing a weighted unbinned maximum likelihood fits can be useful to account for acceptance effects and to statistically subtract background events using the sPlot formalism. It is however well known that the inverse Hessian matrix does not yield parameter uncertainties with correct coverage in the presence of event weights. Three approaches to the determination of parameter uncertainties are compared in this example:
The example performs the fit of a second order polynomial in the angle cos(theta) [-1,1] to a weighted data set. The polynomial is given by
The two coefficients and and their uncertainties are to be determined in the fit.
The per-event weight is used to correct for an acceptance effect, two different acceptance models can be studied:
The performance of the different approaches to determine parameter uncertainties is compared using the pull distributions from a large number of pseudoexperiments. The pull is defined as , where is the fitted parameter and its uncertainty for pseudoexperiment number i. If the fit is unbiased and the parameter uncertainties are estimated correctly, the pull distribution should be a Gaussian centered around zero with a width of one.
{
TH1 *
haccepted =
new TH1D(
"haccepted",
"Generated events;cos(#theta);#events", 40, -1.0, 1.0);
TH1 *
hweighted =
new TH1D(
"hweighted",
"Generated events;cos(#theta);#events", 40, -1.0, 1.0);
std::array<std::string, 3>
methodLabels{
"Inverse weighted Hessian matrix [SumW2Error(false)]",
"Hessian matrix with squared weights [SumW2Error(true)]",
"Asymptotically correct approach [Asymptotic(true)]"};
};
for (std::size_t i = 0; i < 3; ++i) {
std::string
const &
iLabel = std::to_string(i);
}
constexpr std::size_t
ntoys = 500;
constexpr std::size_t
nstats = 500;
constexpr double c0gen = 0.0;
constexpr double c1gen = 0.0;
std::cout <<
"Running " <<
ntoys * 3 <<
" toy fits ..." << std::endl;
for (std::size_t i = 0; i <
ntoys; i++) {
RooRealVar weight(
"weight",
"weight", 0.0, 1000.0);
ngen *= 2.0 / (23.0 / 15.0);
else
ngen *= 2.0 / (16.0 / 15.0);
bool finished = false;
while (!finished) {
double eff = 1.0;
else
weight = 1.0 / eff;
if (10.0 *
rnd->Rndm() < eff *
pol.getVal())
finished = true;
}
}
};
}
std::cout << "... done." << std::endl;
for (std::size_t i = 0; i <
pullHistos.size(); ++i) {
}
}
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
R__EXTERN TStyle * gStyle
Container class to hold unbinned data.
RooExtendPdf is a wrapper around an existing PDF that adds a parameteric extended likelihood term to ...
static RooMsgService & instance()
Return reference to singleton instance.
RooPolynomial implements a polynomial p.d.f of the form.
Variable that can be changed from the outside.
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
virtual void SetTextSize(Float_t tsize=1)
Set the text size.
1-D histogram with a double per channel (see TH1 documentation)
TH1 is the base class of all histogram classes in ROOT.
This class displays a legend box (TPaveText) containing several legend entries.
Random number generator class based on M.
void SetPadTopMargin(Float_t margin=0.1)
void SetOptStat(Int_t stat=1)
The type of information printed in the histogram statistics box can be selected via the parameter mod...
void SetPadBottomMargin(Float_t margin=0.1)
void SetPaintTextFormat(const char *format="g")
void SetEndErrorSize(Float_t np=2)
Set the size (in pixels) of the small lines drawn at the end of the error bars (TH1 or TGraphErrors).
void SetPadRightMargin(Float_t margin=0.1)
void SetTitleOffset(Float_t offset=1, Option_t *axis="X")
Specify a parameter offset to control the distance between the axis and the axis title.
void SetPadLeftMargin(Float_t margin=0.1)
void SetHistLineColor(Color_t color=1)
void SetTitleSize(Float_t size=0.02, Option_t *axis="X")
void SetHistLineWidth(Width_t width=1)
void SetLabelSize(Float_t size=0.04, Option_t *axis="X")
Set size of axis labels.
void SetOptFit(Int_t fit=1)
The type of information about fit parameters printed in the histogram statistics box can be selected ...
RooCmdArg WeightVar(const char *name="weight", bool reinterpretAsWeight=false)
RooCmdArg AsymptoticError(bool flag)
RooCmdArg SumW2Error(bool flag)
RooCmdArg PrintLevel(Int_t code)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Running 1500 toy fits ...
... done.
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 12.2745
NDf = 12
Edm = 3.84566e-06
NCalls = 53
Constant = 82.0512 +/- 4.52865
Mean = -0.042844 +/- 0.0572716
Sigma = 1.19191 +/- 0.0392558 (limited)
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 8.83881
NDf = 9
Edm = 6.19877e-08
NCalls = 53
Constant = 104.798 +/- 5.86613
Mean = -0.0132039 +/- 0.0432773
Sigma = 0.938314 +/- 0.0320986 (limited)
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 7.28099
NDf = 10
Edm = 6.30801e-07
NCalls = 53
Constant = 103.063 +/- 5.55774
Mean = -0.0233563 +/- 0.0435444
Sigma = 0.954358 +/- 0.0285394 (limited)
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 23.1247
NDf = 14
Edm = 2.27317e-06
NCalls = 53
Constant = 69.5383 +/- 3.95119
Mean = -0.160071 +/- 0.0652133
Sigma = 1.37036 +/- 0.0473756 (limited)
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 7.18474
NDf = 15
Edm = 9.01145e-06
NCalls = 61
Constant = 51.3838 +/- 3.49309
Mean = -0.985792 +/- 0.0760226
Sigma = 1.36366 +/- 0.0605719 (limited)
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 11.9456
NDf = 11
Edm = 4.88765e-08
NCalls = 61
Constant = 94.7614 +/- 5.45747
Mean = -0.0978895 +/- 0.0483575
Sigma = 1.03309 +/- 0.0382346 (limited)
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 20.9563
NDf = 16
Edm = 8.77656e-06
NCalls = 53
Constant = 69.8331 +/- 3.99247
Mean = -0.0767062 +/- 0.0640105
Sigma = 1.36919 +/- 0.0474577 (limited)
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 8.84994
NDf = 10
Edm = 1.85865e-06
NCalls = 53
Constant = 101.665 +/- 5.77203
Mean = -0.0558063 +/- 0.0444665
Sigma = 0.96468 +/- 0.0336688 (limited)
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 9.73174
NDf = 10
Edm = 1.95484e-06
NCalls = 53
Constant = 99.7164 +/- 5.69932
Mean = -0.0653755 +/- 0.0459247
Sigma = 0.982023 +/- 0.0349408 (limited)