Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions roofit/roofitcore/inc/RooExtendPdf.h
Original file line number Diff line number Diff line change
Expand Up @@ -48,6 +48,8 @@ class RooExtendPdf : public RooAbsPdf {
double expectedEvents(const RooArgSet* nset) const override ;
std::unique_ptr<RooAbsReal> createExpectedEventsFunc(const RooArgSet* nset) const override;

std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;

RooAbsPdf const& pdf() const { return *_pdf; }

RooAbsReal const &getN() const { return *_n; }
Expand Down
43 changes: 43 additions & 0 deletions roofit/roofitcore/src/RooExtendPdf.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -49,6 +49,8 @@ the nominal integration range \f$ \mathrm{normRegion}[x] \f$.

#include "RooFitImplHelpers.h"

#include <RooFit/Detail/NormalizationHelpers.h>

using std::endl;


Expand Down Expand Up @@ -136,6 +138,47 @@ double RooExtendPdf::expectedEvents(const RooArgSet* nset) const
}


/// RooExtendPdf is self-normalized (selfNormalized() returns true), which means
/// it delegates the normalization of the shape entirely to the wrapped pdf.
/// The fractional correction of the extended term in expectedEvents() /
/// createExpectedEventsFunc() is also evaluated from the wrapped pdf. For a
/// likelihood over a sub-range this means the wrapped pdf must be normalized
/// over the same range as this RooExtendPdf (the fit range), otherwise both the
/// shape term (biasing the shape parameters) and the extended term (collapsing
/// the yield to the fit-range count instead of reinterpreting it to the range
/// requested in the constructor) are wrong.
///
/// Since the normalization range of a pdf is not automatically propagated to
/// its servers, we propagate it explicitly to the wrapped pdf while the
/// computation graph is compiled, mirroring what RooAddPdf does for its
/// component pdfs. This is only done in likelihood mode: for a binned chi2, the
/// per-bin predictions must stay normalized over the full range so that the
/// chi2 remains additive across sub-ranges. See GitHub issue #22959.
std::unique_ptr<RooAbsArg>
RooExtendPdf::compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const
{
// RAII to temporarily set the wrapped pdf's normalization range to the one
// of this RooExtendPdf and reset it after compilation.
struct NormRangeGuard {
NormRangeGuard(RooAbsPdf &pdf, const char *rangeName)
: _pdf{pdf}, _oldNormRange{pdf.normRange() ? pdf.normRange() : ""}
{
_pdf.setNormRange(rangeName);
}
~NormRangeGuard() { _pdf.setNormRange(_oldNormRange.empty() ? nullptr : _oldNormRange.c_str()); }
RooAbsPdf &_pdf;
std::string _oldNormRange;
};

if (!ctx.likelihoodMode()) {
return RooAbsPdf::compileForNormSet(normSet, ctx);
}

NormRangeGuard guard{*_pdf, normRange()};
return RooAbsPdf::compileForNormSet(normSet, ctx);
}


std::unique_ptr<RooAbsReal> RooExtendPdf::createExpectedEventsFunc(const RooArgSet *nset) const
{
const RooAbsPdf& pdf = *_pdf;
Expand Down
59 changes: 59 additions & 0 deletions roofit/roofitcore/test/testTestStatistics.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,7 @@
#include <RooWorkspace.h>
#include <RooFormulaVar.h>
#include <RooGenericPdf.h>
#include <RooAddPdf.h>

#include <TH1D.h>
#include <TMath.h>
Expand Down Expand Up @@ -1212,3 +1213,61 @@ TEST(RooChi2Var, BinnedRangeAdditivityAndNormalization)
EXPECT_NEAR(chi2ExtHi, chi2FuncHi, tol);
EXPECT_NEAR(chi2ExtLoHi, chi2FuncLoHi, tol);
}

// Extended likelihood fit of a RooExtendPdf in a sub-range. This is the
// scenario of the rf204b_extendedLikelihood_rangedFit tutorial. When the
// RooExtendPdf is constructed with a range name, the fitted yield must be
// reinterpreted to that range (here the full range) instead of collapsing to
// the number of events observed in the fit range, and the shape parameters must
// not be biased. Both must agree with the equivalent RooAddPdf model, whose
// coefficients are automatically reinterpreted to the full range. Regression
// test for GitHub issue #22959.
TEST(RooNLLVar, RooExtendPdfRangedFit)
{
using namespace RooFit;
RooHelpers::LocalChangeMsgLevel changeMsgLvl(RooFit::WARNING);

RooRandom::randomGenerator()->SetSeed(42);

RooRealVar x("x", "x", 10, 100);
x.setRange("LEFT", 10, 20);
x.setRange("RIGHT", 60, 100);
x.setRange("FULL", 10, 100);

RooRealVar alphaGen("alpha", "alpha", -0.04, -0.1, -0.0);
RooGenericPdf modelGen("model", "exp(alpha*x)", {x, alphaGen});
std::unique_ptr<RooDataSet> data{modelGen.generate(x, 10000)};

const double nFull = data->sumEntries();
const double nFitRange = std::unique_ptr<RooAbsData>{data->reduce(CutRange("LEFT,RIGHT"))}->sumEntries();
// The fit range only contains a fraction of the events, so this test is only
// meaningful if the two counts are clearly different.
ASSERT_LT(nFitRange, 0.8 * nFull);

// Extended fit with RooExtendPdf, yield N interpreted over the FULL range.
RooRealVar alpha1("alpha", "alpha", -0.04, -0.1, -0.0);
RooGenericPdf model1("model", "exp(alpha*x)", {x, alpha1});
RooRealVar nExpected("N", "Extended term", 5000, 0, 20000);
RooExtendPdf extmodel("extmodel", "Extended model", model1, nExpected, "FULL");
std::unique_ptr<RooFitResult> resExt{extmodel.fitTo(*data, Range("LEFT,RIGHT"), PrintLevel(-1), Save())};

// Equivalent extended fit with a single-component RooAddPdf, whose
// coefficient is automatically reinterpreted to the full range in a ranged
// fit.
RooRealVar alpha2("alpha", "alpha", -0.04, -0.1, -0.0);
RooGenericPdf model2("model", "exp(alpha*x)", {x, alpha2});
RooRealVar nBkg("Nbkg", "Number of background events", 5000, 0, 20000);
RooAddPdf modelsum("modelsum", "", RooArgList{model2}, RooArgList{nBkg});
std::unique_ptr<RooFitResult> resAdd{modelsum.fitTo(*data, Range("LEFT,RIGHT"), PrintLevel(-1), Save())};

// The RooExtendPdf yield must match the RooAddPdf background yield...
EXPECT_NEAR(nExpected.getVal(), nBkg.getVal(), 1e-3 * nBkg.getVal());
// ...and both refer to the full range, not the (much smaller) fit-range count.
EXPECT_NEAR(nExpected.getVal(), nFull, 0.05 * nFull);
EXPECT_GT(nExpected.getVal(), 1.5 * nFitRange);

// The shape parameter must be unbiased, i.e. agree with the RooAddPdf fit and
// with the generated value.
EXPECT_NEAR(alpha1.getVal(), alpha2.getVal(), 1e-4);
EXPECT_NEAR(alpha1.getVal(), -0.04, 5. * alpha1.getError());
}
Original file line number Diff line number Diff line change
Expand Up @@ -83,8 +83,8 @@ void rf204b_extendedLikelihood_rangedFit()
// If the fitted domain is a subset of `FULL`, though, the expected event count is divided by
// \f[
// \mathrm{frac} = \frac{
// \int_{\mathrm{Fit range}} \mathrm{model}(x) \; \mathrm{d}x }{
// \int_{\mathrm{Full range}} \mathrm{model}(x) \; \mathrm{d}x }.
// \int_{\mathrm{Full range}} \mathrm{model}(x) \; \mathrm{d}x }{
// \int_{\mathrm{Fit range}} \mathrm{model}(x) \; \mathrm{d}x }.
// \f]
// `N` will therefore return the count extrapolated to the full range instead of the fit range.
//
Expand Down Expand Up @@ -159,10 +159,9 @@ void rf204b_extendedLikelihood_rangedFit()

// This model will automatically insert the correction factor for the reinterpretation of Nsig and Nnbkg in the full ranges.
//
// When this happens, it reports this with lines like the following:
// When this happens, it reports this with a line like the following:
// ```
// [#1] INFO:Fitting -- RooAbsOptTestStatistic::ctor(nll_modelsum_modelsumData_LEFT) fixing interpretation of coefficients of any RooAddPdf to full domain of observables
// [#1] INFO:Fitting -- RooAbsOptTestStatistic::ctor(nll_modelsum_modelsumData_RIGHT) fixing interpretation of coefficients of any RooAddPdf to full domain of observables
// [#1] INFO:Fitting -- RooAbsPdf::fitTo(modelsum) fixing normalization set for coefficient determination to observables in data
// ```

std::unique_ptr<RooFitResult> r3{modelsum.fitTo(*data, Range("LEFT,RIGHT"), PrintLevel(-1), Save())};
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -72,8 +72,8 @@
# If the fitted domain is a subset of `FULL`, though, the expected event count is divided by
# \f[
# \mathrm{frac} = \frac{
# \int_{\mathrm{Fit range}} \mathrm{model}(x) \; \mathrm{d}x }{
# \int_{\mathrm{Full range}} \mathrm{model}(x) \; \mathrm{d}x }.
# \int_{\mathrm{Full range}} \mathrm{model}(x) \; \mathrm{d}x }{
# \int_{\mathrm{Fit range}} \mathrm{model}(x) \; \mathrm{d}x }.
# \f]
# `N` will therefore return the count extrapolated to the full range instead of the fit range.
#
Expand Down Expand Up @@ -152,10 +152,9 @@
# This model will automatically insert the correction factor for the
# reinterpretation of Nsig and Nnbkg in the full ranges.
#
# When this happens, it reports this with lines like the following:
# When this happens, it reports this with a line like the following:
# ```
# [#1] INFO:Fitting -- RooAbsOptTestStatistic::ctor(nll_modelsum_modelsumData_LEFT) fixing interpretation of coefficients of any RooAddPdf to full domain of observables
# [#1] INFO:Fitting -- RooAbsOptTestStatistic::ctor(nll_modelsum_modelsumData_RIGHT) fixing interpretation of coefficients of any RooAddPdf to full domain of observables
# [#1] INFO:Fitting -- RooAbsPdf::fitTo(modelsum) fixing normalization set for coefficient determination to observables in data
# ```

r3 = modelsum.fitTo(data, Range="LEFT,RIGHT", PrintLevel=-1, Save=True)
Expand Down
Loading