Get the integral of the fitted function with an error in a certain range

Dear RooFit experts,

I am fitting sidebands of the dimuon mass with polynomial background function in the mass intervals of 12-25 and 25-70 GeV. Now I want to know a number of expected background events in the signal region 25-30 GeV and its error.

I take integral of the fitted function in the mass range of 25-30 GeV and it works.
My question is: how I can get an error of the integral ?

Thank you, Sasha.

Hi Sasha, thanks for joining the ROOT forum!

I presume you have calculated the integral using RooAbsReal::createIntegral(), obtaining another RooAbsReal that represents the integral, right?

In this case, you can simply call getPropagatedError() on the integral object, passing your fit result object, to get the error of the integral.

Cheers,
Jonas

Hi Jonas,

thanks ! It work, but gave me a strange result. I do

// define model
 RooPolynomial bkg("bkg","background",dimuonmass1sideband,RooArgList(a0,a1,a2));
  RooRealVar norm("norm","bkg normalization",1100.,0.0,3000.0);
  RooAddPdf model("model","model",RooArgList(bkg),RooArgList(norm));
  // fit 
  RooFitResult* r = model.fitTo(data,Range("R1,R2"),Save());
 // Integral of bkg function in the signal region
  dimuonmass1sideband.setRange("signalrange",25,30);
  RooAbsReal* intPolynX = model.createIntegral(dimuonmass1sideband,NormSet(dimuonmass1sideband),Range("signalrange"));
  Double_t integral_undersignal = intPolynX->getVal();
  Double_t err_integral_undersignal = intPolynX->getPropagatedError(*r);
  cout <<" Integral under the signal region = " << integral_undersignal <<" +/-" << err_integral_undersignal << endl;

and print gives:

Integral under the signal region = 0.133652 +/-1.88998

So, I do not understand why error is larger than integral. Is it in % ?

Cheers, Sasha.

Hi Jonas,

Thanks ! I do not know yet how ROOT Forum works, but I put my replay to your answer through Forum.
I got a strange thing: integral error is larger than value. Please, check my replay for details.

Cheers, Sasha.

Hi Sasha,

sorry for the late reply!

Thanks for sharing the code, you are doing things correctly. Just one thing that you need to change (my fault for not mentioning it before): if you want the get the expected number of background events, you have to multiply the integral with the yield, and then do the error propagation. You have to do this, because the norm parameter and the pdf shape parameters are generally correlated. If you are not propagating the uncertainty up to the final yield, your result will be wrong.

  std::unique_ptr<RooAbsReal> intPolynX{model.createIntegral(dimuonmass1sideband,NormSet(dimuonmass1sideband),Range("signalrange"))};

  RooProduct yieldPolynX{"yieldPolynX", "yieldPolynX", {norm, *intPolynX}};

  double integralUndersignal = yieldPolynX.getVal();
  double errIntegralUndersignal = yieldPolynX.getPropagatedError(*r);
  std::cout << " Integral under the signal region = " << integralUndersignal
            << " +/-" << errIntegralUndersignal << std::endl;

Concerning your surprisingly large integral value: apparently the uncertainty on your PDF parameters is so large that your pdf can become negative within the parameter uncertainties. Have you inspected your fit result? My guess is that for some of the polynomial coefficients, there will be a huge uncertainty. And this is probably because your are trying to fit a 3rd order polynomial, which is not very stable because the parameter are very correlated (remember that the coefficients you pass to the RooPolynomial constructor already start with the first order, because the zeroth coefficient is fixed by the normalization condition).

To stabilize the fit, I would use the more robust parametrization of the Chebychev polynomials implemented in RooChebychev. You can also consider the Bernstein polynomials implemented in RooBernstein. They are not as orthogonal as the Chebychev polynomials, but sometimes better for fitting because they are positive definite and the PDF can’t become negative.

If you stabilized your fit model, the problem with the large uncertainties should also be gone.

I hope I’m on the right track here and my answer helps!

Cheers,
Jonas

Hi Jonas,

it does not work as expected. I got the same situation as with my original code. Printout says

Integral under the signal region = 0.133652 +/-1.88998

You can try it out with macro I use at

/afs/cern.ch/work/a/anikiten/public/CMSSW_11_0_4/src/13TeV_2017UL_RooFitTest_SR1_FTEST.C
and root file with data at the same dir:

test_2017ul_JECDeepCSV.root

I use interactive root and run my macro with

.L test_2017ul_JECDeepCSV.root
Draw()

Thanks, Sasha.

Hi Sasha,

thanks for sharing the script and the ROOT file.

Indeed, your problem is not the one I thought, and I think I was wrong with my initial suggestion of doing the error propagation.

For a mathematical reason that I can’t exactly pin-point my finger on just yet, the error propagation is not valid in this case. An safer alternative that just came to mind now that I see your script without error propagation would be to reformulate your model such that you don’t need it, and one of the parameters is already the background yield in the signal region.

In the RooAbsPdf::fitTo() documentation, you see explained the SumCoefRange() option to define for which range the yield coefficients are defined. You can pass SumCoefRange("signalrange") to your fitTo call, and the norm parameter will be directly the normalization you are looking for. Then you can drop all the other integration stuff.

I hope it’s fine for you now, otherwise feel free to follow up again!

Cheers,
Jonas

Dear Jonas,

I have just read second part of your mail about possible unstabillities of the fit.
I will check tomorrow and let you know.

But when I have plotted fitted function with the errors the error bar looks fine.

Thanks a lot for your help.

Hi Jonas,

you suggestion works ! I will stay with it.

Thank you very much for all your help; very much appreciated !
Cheers, Sasha.

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.