Including efficiencies in simultaneous extended likelihood fit with constraints

[ilam]

Hi everyone,

I have two data sets (DS1 and DS2) that I am fitting simultaneously. Each data set has different run time lengths and different selection cuts, and therefore different efficiencies. Say there is a background that is common between both of them, for example, the solar neutrino flux. The solar neutrino flux is well-measured and so I know how many solar neutrino events I expect in each data set. The PDFs of the solar neutrinos are obtainable from Monte Carlo simulations. The events generated from these simulations are subjected to the same selection cuts unique to each data set.

I would like to constrain the expected number of solar neutrinos in my fit as a single parameter common to both data sets. In RooFit, I can do this by including a single Gaussian constrain term in my likelihood. However, due to the fact that each data set has different efficiencies and live times (assume these values are known), I am not sure how to include the efficiencies and live times into RooFit in Extended() mode.

One way I was considering is to fit for an overall event rate rather than the number of events, and constrain the expected solar neutrino rate instead of the expected number of solar neutrino events. That way, the live time and selection efficiencies are included. But I am not sure how to do this using RooFit. Is it possible?

In short, I have the parameter n_solar that is floated in both DS1 and DS2 but each DS1 and DS2 have their own efficiencies and livetimes such that after taken into account, would return the expected solar neutrino rate (which is constrained). Any comments and guidance is appreciated!

[StephanH]

Hi Ilam,

Just take the rate as a common parameter and from that calculate the expected number of events per dataset. I guess something like this should do it:

RooRealVar neutrinoRate ...
RooRealVar eff1
RooRealVar eff2
RooLinearVar events1("events1", "events1", neutrinoRate, eff1, RooConst(0.));
RooLinearVar events2("events2", "events2", neutrinoRate, eff2, RooConst(0.));

RooGaussian constraint(..., targetRate, neutrinoRate, sigma)

[ilam]

Hi StephanH,

Thanks for your response! I created a working example based on your suggestion by modifying tutorial rf_202. I attach the code here test_simul_eventrate.cpp (8.4 KB) . I think it is working, as it is returning a rate that I gave it. But I want to make sure if I understood your suggestion. Does this look right? The plot output from my code is here. The left plot is the first dataset. The right plot is the second dataset. (

).

Also, in the code, I generate some fake data from the model to perform the tests. But I am unsure of how the routine determines how many events to generate in this context. In the usual case, the RooRealVar passed is the number of events (nsig, nbkg), but in the context of this code, they don’t get passed around. About 145 events seem to be generated in the code but I am not sure where RooFit gets this number from.

[StephanH]

Hi @ilam,

I only looked at the coefficients of the RooAddPdf, as that’s what will determine how many events are generated. You implemented

events_bkg = bkgrate * eff1
events_bkg_ds2 = bkgrate * eff1_ds2

I see that they have a common rate, but different efficiencies. If I’m not mistaken, that’s what you asked for. If you want to check the values of such variables yourself, I suggest to print them.

root [1] events_bkg.Print("T")
0x7fff9c160af8 RooLinearVar::events_bkg = 20 [Auto,Clean] 
  0x7fff9c161c10/VS RooRealVar::bkgrate = 100
  0x7fff9c161440/V- RooRealVar::eff1 = 0.2
  0x1c42210/V- RooConstVar::0 = 0

So yes, you see that it’s 0.2 * 100 = 20 in your case.

[ilam]

Hi @StephanH,

Sorry for the late reply. Yes, that was what I wanted. And yes, those numbers makes sense.

Thanks again!