HybridInstructional.C: example demostrating usage of HybridCalcultor | Roostats tutorials | HybridStandardForm.C: A hypothesis testing example based on number counting with background uncertainty. |
// Example on how to use the HybridCalculatorOriginal class // // Author: Gregory Schott #include "RooRandom.h" #include "RooRealVar.h" #include "RooGaussian.h" #include "RooPolynomial.h" #include "RooArgSet.h" #include "RooAddPdf.h" #include "RooDataSet.h" #include "RooExtendPdf.h" #include "RooConstVar.h" #ifndef __CINT__ // problem including this file with CINT #include "RooGlobalFunc.h" #endif #include "RooStats/HybridCalculatorOriginal.h" #include "RooStats/HybridResult.h" #include "RooStats/HybridPlot.h" void HybridOriginalDemo(int ntoys = 1000) { //***********************************************************************// // This macro show an example on how to use RooStats/HybridCalculatorOriginal // //***********************************************************************// // // With this example, you should get: CL_sb = 0.130 and CL_b = 0.946 // (if data had -2lnQ = -3.0742). You can compare to the expected plot: // http://www-ekp.physik.uni-karlsruhe.de/~schott/roostats/hybridplot_example.png using namespace RooFit; using namespace RooStats; /// set RooFit random seed RooRandom::randomGenerator()->SetSeed(3007); /// build the models for background and signal+background RooRealVar x("x","",-3,3); RooArgList observables(x); // variables to be generated // gaussian signal // RooRealVar sig_mean("sig_mean","",0); // RooRealVar sig_sigma("sig_sigma","",0.8); // RooGaussian sig_pdf("sig_pdf","",x,sig_mean,sig_sigma); RooGaussian sig_pdf("sig_pdf","",x, RooConst(0.0),RooConst(0.8)); RooRealVar sig_yield("sig_yield","",20,0,300); // flat background (extended PDF) // RooRealVar bkg_slope("bkg_slope","",0); // RooPolynomial bkg_pdf("bkg_pdf","",x,bkg_slope); RooPolynomial bkg_pdf("bkg_pdf","", x, RooConst(0)); RooRealVar bkg_yield("bkg_yield","",40,0,300); RooExtendPdf bkg_ext_pdf("bkg_ext_pdf","",bkg_pdf,bkg_yield); // bkg_yield.setConstant(kTRUE); sig_yield.setConstant(kTRUE); // total sig+bkg (extended PDF) RooAddPdf tot_pdf("tot_pdf","",RooArgList(sig_pdf,bkg_pdf),RooArgList(sig_yield,bkg_yield)); /// build the prior PDF on the parameters to be integrated // gaussian contraint on the background yield ( N_B = 40 +/- 10 ie. 25% ) RooGaussian bkg_yield_prior("bkg_yield_prior","",bkg_yield,RooConst(bkg_yield.getVal()),RooConst(10.)); RooArgSet nuisance_parameters(bkg_yield); // variables to be integrated /// generate a data sample RooDataSet* data = tot_pdf.generate(observables,RooFit::Extended()); //***********************************************************************// /// run HybridCalculator on those inputs // use interface from HypoTest calculator by default HybridCalculatorOriginal myHybridCalc(*data, tot_pdf , bkg_ext_pdf , &nuisance_parameters, &bkg_yield_prior); // here I use the default test statistics: 2*lnQ (optional) myHybridCalc.SetTestStatistic(1); //myHybridCalc.SetTestStatistic(3); // profile likelihood ratio myHybridCalc.SetNumberOfToys(ntoys); myHybridCalc.UseNuisance(true); // for speed up generation (do binned data) myHybridCalc.SetGenerateBinned(false); // calculate by running ntoys for the S+B and B hypothesis and retrieve the result HybridResult* myHybridResult = myHybridCalc.GetHypoTest(); if (! myHybridResult) { std::cerr << "\nError returned from Hypothesis test" << std::endl; return; } /// run 1000 toys without gaussian prior on the background yield //HybridResult* myHybridResult = myHybridCalc.Calculate(*data,1000,false); /// save the toy-MC results to file, this way splitting into sub-batch jobs is possible //TFile fileOut("some_hybridresult.root","RECREATE"); //fileOut.cd(); //myHybridResult.Write(); //fileOut.Close(); /// read the results from a file //TFile fileIn("some_hybridresult.root"); //HybridResult* myOtherHybridResult = (HybridResult*) fileIn.Get("myHybridCalc"); /// example on how to merge with toy-MC results obtained in another job //HybridResult* mergedHybridResult = new HybridResult("mergedHybridResult","this object holds merged results"); //mergedHybridResult->Add(myHybridResult); //mergedHybridResult->Add(myOtherHybridResult); /// or //myHybridResult->Add(myOtherHybridResult); /// nice plot of the results HybridPlot* myHybridPlot = myHybridResult->GetPlot("myHybridPlot","Plot of results with HybridCalculatorOriginal",100); myHybridPlot->Draw(); /// recover and display the results double clsb_data = myHybridResult->CLsplusb(); double clb_data = myHybridResult->CLb(); double cls_data = myHybridResult->CLs(); double data_significance = myHybridResult->Significance(); double min2lnQ_data = myHybridResult->GetTestStat_data(); /// compute the mean expected significance from toys double mean_sb_toys_test_stat = myHybridPlot->GetSBmean(); myHybridResult->SetDataTestStatistics(mean_sb_toys_test_stat); double toys_significance = myHybridResult->Significance(); std::cout << "Completed HybridCalculatorOriginal example:\n"; std::cout << " - -2lnQ = " << min2lnQ_data << endl; std::cout << " - CL_sb = " << clsb_data << std::endl; std::cout << " - CL_b = " << clb_data << std::endl; std::cout << " - CL_s = " << cls_data << std::endl; std::cout << " - significance of data = " << data_significance << std::endl; std::cout << " - mean significance of toys = " << toys_significance << std::endl; } HybridOriginalDemo.C:1 HybridOriginalDemo.C:2 HybridOriginalDemo.C:3 HybridOriginalDemo.C:4 HybridOriginalDemo.C:5 HybridOriginalDemo.C:6 HybridOriginalDemo.C:7 HybridOriginalDemo.C:8 HybridOriginalDemo.C:9 HybridOriginalDemo.C:10 HybridOriginalDemo.C:11 HybridOriginalDemo.C:12 HybridOriginalDemo.C:13 HybridOriginalDemo.C:14 HybridOriginalDemo.C:15 HybridOriginalDemo.C:16 HybridOriginalDemo.C:17 HybridOriginalDemo.C:18 HybridOriginalDemo.C:19 HybridOriginalDemo.C:20 HybridOriginalDemo.C:21 HybridOriginalDemo.C:22 HybridOriginalDemo.C:23 HybridOriginalDemo.C:24 HybridOriginalDemo.C:25 HybridOriginalDemo.C:26 HybridOriginalDemo.C:27 HybridOriginalDemo.C:28 HybridOriginalDemo.C:29 HybridOriginalDemo.C:30 HybridOriginalDemo.C:31 HybridOriginalDemo.C:32 HybridOriginalDemo.C:33 HybridOriginalDemo.C:34 HybridOriginalDemo.C:35 HybridOriginalDemo.C:36 HybridOriginalDemo.C:37 HybridOriginalDemo.C:38 HybridOriginalDemo.C:39 HybridOriginalDemo.C:40 HybridOriginalDemo.C:41 HybridOriginalDemo.C:42 HybridOriginalDemo.C:43 HybridOriginalDemo.C:44 HybridOriginalDemo.C:45 HybridOriginalDemo.C:46 HybridOriginalDemo.C:47 HybridOriginalDemo.C:48 HybridOriginalDemo.C:49 HybridOriginalDemo.C:50 HybridOriginalDemo.C:51 HybridOriginalDemo.C:52 HybridOriginalDemo.C:53 HybridOriginalDemo.C:54 HybridOriginalDemo.C:55 HybridOriginalDemo.C:56 HybridOriginalDemo.C:57 HybridOriginalDemo.C:58 HybridOriginalDemo.C:59 HybridOriginalDemo.C:60 HybridOriginalDemo.C:61 HybridOriginalDemo.C:62 HybridOriginalDemo.C:63 HybridOriginalDemo.C:64 HybridOriginalDemo.C:65 HybridOriginalDemo.C:66 HybridOriginalDemo.C:67 HybridOriginalDemo.C:68 HybridOriginalDemo.C:69 HybridOriginalDemo.C:70 HybridOriginalDemo.C:71 HybridOriginalDemo.C:72 HybridOriginalDemo.C:73 HybridOriginalDemo.C:74 HybridOriginalDemo.C:75 HybridOriginalDemo.C:76 HybridOriginalDemo.C:77 HybridOriginalDemo.C:78 HybridOriginalDemo.C:79 HybridOriginalDemo.C:80 HybridOriginalDemo.C:81 HybridOriginalDemo.C:82 HybridOriginalDemo.C:83 HybridOriginalDemo.C:84 HybridOriginalDemo.C:85 HybridOriginalDemo.C:86 HybridOriginalDemo.C:87 HybridOriginalDemo.C:88 HybridOriginalDemo.C:89 HybridOriginalDemo.C:90 HybridOriginalDemo.C:91 HybridOriginalDemo.C:92 HybridOriginalDemo.C:93 HybridOriginalDemo.C:94 HybridOriginalDemo.C:95 HybridOriginalDemo.C:96 HybridOriginalDemo.C:97 HybridOriginalDemo.C:98 HybridOriginalDemo.C:99 HybridOriginalDemo.C:100 HybridOriginalDemo.C:101 HybridOriginalDemo.C:102 HybridOriginalDemo.C:103 HybridOriginalDemo.C:104 HybridOriginalDemo.C:105 HybridOriginalDemo.C:106 HybridOriginalDemo.C:107 HybridOriginalDemo.C:108 HybridOriginalDemo.C:109 HybridOriginalDemo.C:110 HybridOriginalDemo.C:111 HybridOriginalDemo.C:112 HybridOriginalDemo.C:113 HybridOriginalDemo.C:114 HybridOriginalDemo.C:115 HybridOriginalDemo.C:116 HybridOriginalDemo.C:117 HybridOriginalDemo.C:118 HybridOriginalDemo.C:119 HybridOriginalDemo.C:120 HybridOriginalDemo.C:121 HybridOriginalDemo.C:122 HybridOriginalDemo.C:123 HybridOriginalDemo.C:124 HybridOriginalDemo.C:125 HybridOriginalDemo.C:126 HybridOriginalDemo.C:127 HybridOriginalDemo.C:128 HybridOriginalDemo.C:129 HybridOriginalDemo.C:130 HybridOriginalDemo.C:131 HybridOriginalDemo.C:132 HybridOriginalDemo.C:133 HybridOriginalDemo.C:134 HybridOriginalDemo.C:135 HybridOriginalDemo.C:136 HybridOriginalDemo.C:137 HybridOriginalDemo.C:138 HybridOriginalDemo.C:139 HybridOriginalDemo.C:140 HybridOriginalDemo.C:141 HybridOriginalDemo.C:142 HybridOriginalDemo.C:143 HybridOriginalDemo.C:144 |
|