Making an Asimov dataset

Roofit and RooStats
1

Hello,

I wanted to generate an Asimov dataset for my model, corresponding to some choice for the parameter of interest. My current understanding is that this involves setting any nuisance parameters equal to their conditional MLE values (using the observed data), and then reading off the expected values of all the observables. And in the large sample limit, we expect to find that fitting our model to the asimov dataset will recover the parameter of interest value. I was curious if the “ExpectedData()” options of RooFit’s generate method could do this for me. To compare, I made a simple model and tried to make the asimov datasets:

{
   using namespace RooFit;

   RooWorkspace w;

   w.factory("PROD::model(Poisson::pois(n[61,0,100],expr::mean('a*b',a[10,0,100],b[10,0,100])),Gaussian::gaus(b_obs[5,0,100],b,3))");

   w.defineSet("poi","a");
   w.defineSet("np","b");
   w.defineSet("obs","n,b_obs");

   RooDataSet data("data","data",*w.set("obs"));
   data.add(*w.set("obs"));

   //--------------------------------------------------------

   //Method 1. Set np to their conditional MLE values using the observed data, then calculate the expected observables
   w.var("a")->setVal(4);w.var("a")->setConstant(true);
   w.pdf("model")->fitTo(data);
   w.var("a")->setConstant(false);

   RooDataSet asimov("asimov","asimov",*w.set("obs"));
   w.var("n").setVal(w.function("mean")->getVal());
   w.var("b_obs").setVal(w.function("b")->getVal());
   asimov.add(*w.set("obs"));

   //check the asimov dataset will give back a=4
   w.pdf("model")->fitTo(asimov);
   std::cout << "Asimov designed for a=4 gave a-hat = " << w.var("a")->getVal() << std::endl;

   //--------------------------------------------------------

   //Method 2. Use ExpectedData option of "generate"
   w.var("a")->setVal(4);
   RooDataSet* asimov2 = w.pdf("model")->generate(*w.set("obs"),1,false,true,"",true/*This is expected data option*/);
   RooDataSet* asimov3 = w.pdf("model")->generate(*w.set("obs"),NumEvents(1),ExpectedData()); //this generated 20000 entries!??
   //check if this gives back a=4
   w.pdf("model")->fitTo(*asimov2);
   std::cout << "Asimov2 designed for a=4 gave a-hat = " << w.var("a")->getVal() << std::endl;
}

The results I get are 3.98 for my first method, but 3.4 for the second method. The effort to use ExpectedData() didn’t make sense at all.

Is there a trick to using roofit to build asimov datasets like this? The frustrating part at the moment is that I will have to copy over all the expected values (means in the poissons and gaussians that make up my model) to the observables - it would be nice if Roofit is able to do this automatically, which is what I hoped ExpectedData could do?

Thanks

2

Hi,
The function you should use to make an Asimov data set which correctly takes into account the global observables is the one in the RooStats::AsymptoticCalculator, MakeAsimovData

root.cern.ch/root/html/RooStats_ … AsimovData

Here is an example on how works on your model (you have considered blobs as an extra observable and not global observable, i.e. it is part of the data set):

   ModelConfig mc(&w); 
   mc.SetPdf("model");
   mc.SetObservables(*w.set("obs"));
   mc.SetParametersOfInterest(*w.set("poi"));
   mc.SetNuisanceParameters(*w.set("np"));


   RooAbsData * asimov = AsymptoticCalculator::MakeAsimovData(data, mc, RooArgSet(*w.var("a"),*w.var("b")), RooArgSet());

Lorenzo