Hi @AlkaidCheng!
Sorry for the late reply. Before writing more about your specific question, here are some good references about constraint terms and global observables:
- The documentation of RooAbsPdf::createNLL(), especially the explanation for the
Constran()
and GlobalObservables()
command arguments
- An excellent forum post by @stephanh about global observables and constrained parameters
- For experts: the source code of the
createConstraintTerm
function that is also used by createNLL
to double check what really happens
-
RooFit tutorial on parameter constraints
Now, what is the right choice for you? Let’s look at your propositions again, noting what are used as default constraint parameters and global observables if you don’t pass them explicitly:
ModelConfig *mc = ...;
RooAbsData *data = ...;
RooAbsPdf *pdf=mc->GetPdf();
std::unique_ptr<RooArgSet> allParams{pdf->getParameters(*data)};
# choice 1
RooAbsReal* nll=pdf->createNLL(*data,
Constrain(*mc->GetNuisanceParameters()),
GlobalObservables(*mc->GetGlobalObservables()),
...);
# choice 2
RooAbsReal* nll=pdf->createNLL(*data,
// Constrained parameters not specified, so all parameters
// are constrained if there are constraint terms for them.
// It's like:
// Constrain(*allParams),
GlobalObservables(*mc->GetGlobalObservables()),
...);
# choice 3
RooAbsReal* nll=pdf->createNLL(*data,
Constrain(*mc->GetNuisanceParameters()),
// GlobalObservables not specified,
// so the constrained parameters are used, i.e. it's like:
// GlobalObservables(*mc->GetNuisanceParameters()),
...);
Remember that the global observables are the variables that the constraint term is integrated over to determine its normalization. So not setting these global observables correctly is wrong if your constraints are not symmetric under the exchange of global observables and nuisance parameters! However, most constraints in the wild are Gaussian, so often you don’t notice that mistake, which is why you sill see your choice 3 in the wild even though it is conceptually wrong.
Choice 1 and 2 are both correct. In one you are just passing the set of constrained parameters explicitly, in the other not. Passing them explicitly can be useful if you study how removing some of the constraints will affect your fit result.
Just for you to try out, I have created also an example that showcases how you get a different fit result in a model with a Poisson constraint if you forget to set the global observables correctly:
void constraintsOptionsDemo()
{
using namespace RooFit;
// Construct a Gaussian pdf
RooRealVar x("x", "x", -10, 10);
RooRealVar m("m", "m", 5, 0, 10);
RooRealVar s("s", "s", 2, 0.1, 10);
RooGaussian gauss("gauss", "gauss(x,m,s)", x, m, s);
// To reset parameters after fit to not bias results of the next fit
auto resetParameters = [&](){
m.setVal(5);
m.setError(0);
s.setVal(2);
s.setError(0);
};
// Generate small dataset for our tests
std::unique_ptr<RooDataSet> d{gauss.generate(x, 50)};
// Create the global observables, which is constant (not a parameter)
RooRealVar mGlobalObs("mGlobalObs", "mGlobalObs", 4, 0, 10);
mGlobalObs.setConstant();
// Construct Gaussian constraint pdf on parameter m
RooPoisson constraint("constraint", "constraint", mGlobalObs, m);
// Alternatively, you can try a Gaussian constraint, and you will see that
// not defining the correct global observable doesn't matter here because
// the integral of a Gaussian is symmetric under exchange of x and my.
// RooGaussian constraint("constraint", "constraint", mGlobalObs, m, RooConst(2));
// Multiply constraint with pdf
RooProdPdf model("modec", "model with constraint", RooArgSet(gauss, constraint));
// Create NLLs with different options
std::unique_ptr<RooAbsReal> nll1{model.createNLL(*d)};
std::unique_ptr<RooAbsReal> nll2{model.createNLL(*d, Constrain(m))};
std::unique_ptr<RooAbsReal> nll3{model.createNLL(*d, Constrain(m), GlobalObservables(mGlobalObs))};
// Do some fitting now, to show which differences in the NLL are not just an
// absolute offset with no influence on the fit result.
using R = std::unique_ptr<RooFitResult>;
R r1{model.fitTo(*d, Save(), PrintLevel(-1))};
resetParameters();
R r2{model.fitTo(*d, Constrain(m), Save(), PrintLevel(-1))};
resetParameters();
R r3{model.fitTo(*d, Constrain(m), GlobalObservables(mGlobalObs), Save(), PrintLevel(-1))};
resetParameters();
// Print the likelihood values
std::cout << nll1->getVal() << std::endl;
std::cout << nll2->getVal() << std::endl;
std::cout << nll3->getVal() << std::endl;
// Print the fit results
r1->Print();
r2->Print();
r3->Print();
}
Here is the output:
103.932
103.932
103.948
RooFitResult: minimized FCN value: 103.81, estimated distance to minimum: 9.02216e-08
covariance matrix quality: Full, accurate covariance matrix
Status : MINIMIZE=0 HESSE=0
Floating Parameter FinalValue +/- Error
-------------------- --------------------------
m 5.0597e+00 +/- 2.74e-01
s 1.9075e+00 +/- 2.06e-01
RooFitResult: minimized FCN value: 103.81, estimated distance to minimum: 8.64109e-08
covariance matrix quality: Full, accurate covariance matrix
Status : MINIMIZE=0 HESSE=0
Floating Parameter FinalValue +/- Error
-------------------- --------------------------
m 5.0597e+00 +/- 2.74e-01
s 1.9075e+00 +/- 2.06e-01
RooFitResult: minimized FCN value: 103.825, estimated distance to minimum: 8.50857e-08
covariance matrix quality: Full, accurate covariance matrix
Status : MINIMIZE=0 HESSE=0
Floating Parameter FinalValue +/- Error
-------------------- --------------------------
m 5.0612e+00 +/- 2.74e-01
s 1.9076e+00 +/- 2.06e-01
I hope this clarifies a few things!
Cheers,
Jonas