Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
HybridInstructional.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_roostats
3/// \notebook -js
4/// Example demonstrating usage of HybridCalcultor
5///
6/// A hypothesis testing example based on number counting
7/// with background uncertainty.
8///
9/// NOTE: This example must be run with the ACLIC (the + option ) due to the
10/// new class that is defined.
11///
12/// This example:
13/// - demonstrates the usage of the HybridCalcultor (Part 4-6)
14/// - demonstrates the numerical integration of RooFit (Part 2)
15/// - validates the RooStats against an example with a known analytic answer
16/// - demonstrates usage of different test statistics
17/// - explains subtle choices in the prior used for hybrid methods
18/// - demonstrates usage of different priors for the nuisance parameters
19///
20/// The basic setup here is that a main measurement has observed x events with an
21/// expectation of s+b. One can choose an ad hoc prior for the uncertainty on b,
22/// or try to base it on an auxiliary measurement. In this case, the auxiliary
23/// measurement (aka control measurement, sideband) is another counting experiment
24/// with measurement y and expectation tau*b. With an 'original prior' on b,
25/// called \f$\eta(b)\f$ then one can obtain a posterior from the auxiliary measurement
26/// \f$\pi(b) = \eta(b) * Pois(y|tau*b).\f$ This is a principled choice for a prior
27/// on b in the main measurement of x, which can then be treated in a hybrid
28/// Bayesian/Frequentist way. Additionally, one can try to treat the two
29/// measurements simultaneously, which is detailed in Part 6 of the tutorial.
30///
31/// This tutorial is related to the FourBin.C tutorial in the modeling, but
32/// focuses on hypothesis testing instead of interval estimation.
33///
34/// More background on this 'prototype problem' can be found in the
35/// following papers:
36///
37/// - Evaluation of three methods for calculating statistical significance
38/// when incorporating a systematic uncertainty into a test of the
39/// background-only hypothesis for a Poisson process
40/// Authors: Robert D. Cousins, James T. Linnemann, Jordan Tucker
41/// http://arxiv.org/abs/physics/0702156
42/// NIM A 595 (2008) 480--501
43///
44/// - Statistical Challenges for Searches for New Physics at the LHC
45/// Author: Kyle Cranmer
46/// http://arxiv.org/abs/physics/0511028
47///
48/// - Measures of Significance in HEP and Astrophysics
49/// Authors J. T. Linnemann
50/// http://arxiv.org/abs/physics/0312059
51///
52/// \macro_image
53/// \macro_output
54/// \macro_code
55///
56/// \authors Kyle Cranmer, Wouter Verkerke, and Sven Kreiss
57
58#include "RooGlobalFunc.h"
59#include "RooRealVar.h"
60#include "RooProdPdf.h"
61#include "RooWorkspace.h"
62#include "RooDataSet.h"
63#include "TCanvas.h"
64#include "TStopwatch.h"
65#include "TH1.h"
66#include "RooPlot.h"
67#include "RooMsgService.h"
68
70
74
79
80using namespace RooFit;
81using namespace RooStats;
82
83// ----------------------------------
84// A New Test Statistic Class for this example.
85// It simply returns the sum of the values in a particular
86// column of a dataset.
87// You can ignore this class and focus on the macro below
88class BinCountTestStat : public TestStatistic {
89public:
90 BinCountTestStat(void) : fColumnName("tmp") {}
91 BinCountTestStat(string columnName) : fColumnName(columnName) {}
92
93 virtual Double_t Evaluate(RooAbsData &data, RooArgSet & /*nullPOI*/)
94 {
95 // This is the main method in the interface
96 Double_t value = 0.0;
97 for (int i = 0; i < data.numEntries(); i++) {
98 value += data.get(i)->getRealValue(fColumnName.c_str());
99 }
100 return value;
101 }
102 virtual const TString GetVarName() const { return fColumnName; }
103
104private:
105 string fColumnName;
106
107protected:
109};
110
111
112// ----------------------------------
113// The Actual Tutorial Macro
114
115void HybridInstructional(int ntoys = 6000)
116{
117 double nToysRatio = 20; // ratio Ntoys Null/ntoys ALT
118
119 // This tutorial has 6 parts
120 // Table of Contents
121 // Setup
122 // 1. Make the model for the 'prototype problem'
123 // Special cases
124 // 2. Use RooFit's direct integration to get p-value & significance
125 // 3. Use RooStats analytic solution for this problem
126 // RooStats HybridCalculator -- can be generalized
127 // 4. RooStats ToyMC version of 2. & 3.
128 // 5. RooStats ToyMC with an equivalent test statistic
129 // 6. RooStats ToyMC with simultaneous control & main measurement
130
131 TStopwatch t;
132 t.Start();
133 TCanvas *c = new TCanvas;
134 c->Divide(2, 2);
135
136 // ----------------------------------------------------
137 // P A R T 1 : D I R E C T I N T E G R A T I O N
138 // ====================================================
139 // Make model for prototype on/off problem
140 // $Pois(x | s+b) * Pois(y | tau b )$
141 // for Z_Gamma, use uniform prior on b.
142 RooWorkspace *w = new RooWorkspace("w");
143 w->factory("Poisson::px(x[150,0,500],sum::splusb(s[0,0,100],b[100,0.1,300]))");
144 w->factory("Poisson::py(y[100,0.1,500],prod::taub(tau[1.],b))");
145 w->factory("PROD::model(px,py)");
146 w->factory("Uniform::prior_b(b)");
147
148 // We will control the output level in a few places to avoid
149 // verbose progress messages. We start by keeping track
150 // of the current threshold on messages.
152
153 // ----------------------------------------------------
154 // P A R T 2 : D I R E C T I N T E G R A T I O N
155 // ====================================================
156 // This is not the 'RooStats' way, but in this case the distribution
157 // of the test statistic is simply x and can be calculated directly
158 // from the PDF using RooFit's built-in integration.
159 // Note, this does not generalize to situations in which the test statistic
160 // depends on many events (rows in a dataset).
161
162 // construct the Bayesian-averaged model (eg. a projection pdf)
163 // $p'(x|s) = \int db p(x|s+b) * [ p(y|b) * prior(b) ]$
164 w->factory("PROJ::averagedModel(PROD::foo(px|b,py,prior_b),b)");
165
166 RooMsgService::instance().setGlobalKillBelow(RooFit::ERROR); // lower message level
167 // plot it, red is averaged model, green is b known exactly, blue is s+b av model
168 RooPlot *frame = w->var("x")->frame(Range(50, 230));
169 w->pdf("averagedModel")->plotOn(frame, LineColor(kRed));
170 w->pdf("px")->plotOn(frame, LineColor(kGreen));
171 w->var("s")->setVal(50.);
172 w->pdf("averagedModel")->plotOn(frame, LineColor(kBlue));
173 c->cd(1);
174 frame->Draw();
175 w->var("s")->setVal(0.);
176
177 // compare analytic calculation of Z_Bi
178 // with the numerical RooFit implementation of Z_Gamma
179 // for an example with x = 150, y = 100
180
181 // numeric RooFit Z_Gamma
182 w->var("y")->setVal(100);
183 w->var("x")->setVal(150);
184 std::unique_ptr<RooAbsReal> cdf{w->pdf("averagedModel")->createCdf(*w->var("x"))};
185 cdf->getVal(); // get ugly print messages out of the way
186 cout << "-----------------------------------------" << endl;
187 cout << "Part 2" << endl;
188 cout << "Hybrid p-value from direct integration = " << 1 - cdf->getVal() << endl;
189 cout << "Z_Gamma Significance = " << PValueToSignificance(1 - cdf->getVal()) << endl;
190 RooMsgService::instance().setGlobalKillBelow(msglevel); // set it back
191
192 // ---------------------------------------------
193 // P A R T 3 : A N A L Y T I C R E S U L T
194 // =============================================
195 // In this special case, the integrals are known analytically
196 // and they are implemented in RooStats::NumberCountingUtils
197
198 // analytic Z_Bi
199 double p_Bi = NumberCountingUtils::BinomialWithTauObsP(150, 100, 1);
200 double Z_Bi = NumberCountingUtils::BinomialWithTauObsZ(150, 100, 1);
201 cout << "-----------------------------------------" << endl;
202 cout << "Part 3" << endl;
203 std::cout << "Z_Bi p-value (analytic): " << p_Bi << std::endl;
204 std::cout << "Z_Bi significance (analytic): " << Z_Bi << std::endl;
205 t.Stop();
206 t.Print();
207 t.Reset();
208 t.Start();
209
210 // -------------------------------------------------------------
211 // P A R T 4 : U S I N G H Y B R I D C A L C U L A T O R
212 // =============================================================
213 // Now we demonstrate the RooStats HybridCalculator.
214 //
215 // Like all RooStats calculators it needs the data and a ModelConfig
216 // for the relevant hypotheses. Since we are doing hypothesis testing
217 // we need a ModelConfig for the null (background only) and the alternate
218 // (signal+background) hypotheses. We also need to specify the PDF,
219 // the parameters of interest, and the observables. Furthermore, since
220 // the parameter of interest is floating, we need to specify which values
221 // of the parameter corresponds to the null and alternate (eg. s=0 and s=50)
222 //
223 // define some sets of variables obs={x} and poi={s}
224 // note here, x is the only observable in the main measurement
225 // and y is treated as a separate measurement, which is used
226 // to produce the prior that will be used in this calculation
227 // to randomize the nuisance parameters.
228 w->defineSet("obs", "x");
229 w->defineSet("poi", "s");
230
231 // create a toy dataset with the x=150
232 RooDataSet *data = new RooDataSet("d", "d", *w->set("obs"));
233 data->add(*w->set("obs"));
234
235 // Part 3a : Setup ModelConfigs
236 // -------------------------------------------------------
237 // create the null (background-only) ModelConfig with s=0
238 ModelConfig b_model("B_model", w);
239 b_model.SetPdf(*w->pdf("px"));
240 b_model.SetObservables(*w->set("obs"));
241 b_model.SetParametersOfInterest(*w->set("poi"));
242 w->var("s")->setVal(0.0); // important!
243 b_model.SetSnapshot(*w->set("poi"));
244
245 // create the alternate (signal+background) ModelConfig with s=50
246 ModelConfig sb_model("S+B_model", w);
247 sb_model.SetPdf(*w->pdf("px"));
248 sb_model.SetObservables(*w->set("obs"));
249 sb_model.SetParametersOfInterest(*w->set("poi"));
250 w->var("s")->setVal(50.0); // important!
251 sb_model.SetSnapshot(*w->set("poi"));
252
253 // Part 3b : Choose Test Statistic
254 // ----------------------------------
255 // To make an equivalent calculation we need to use x as the test
256 // statistic. This is not a built-in test statistic in RooStats
257 // so we define it above. The new class inherits from the
258 // RooStats::TestStatistic interface, and simply returns the value
259 // of x in the dataset.
260
262
263 // Part 3c : Define Prior used to randomize nuisance parameters
264 // -------------------------------------------------------------
265 //
266 // The prior used for the hybrid calculator is the posterior
267 // from the auxiliary measurement y. The model for the aux.
268 // measurement is Pois(y|tau*b), thus the likelihood function
269 // is proportional to (has the form of) a Gamma distribution.
270 // if the 'original prior' $\eta(b)$ is uniform, then from
271 // Bayes's theorem we have the posterior:
272 // $\pi(b) = Pois(y|tau*b) * \eta(b)$
273 // If $\eta(b)$ is flat, then we arrive at a Gamma distribution.
274 // Since RooFit will normalize the PDF we can actually supply
275 // $py=Pois(y,tau*b)$ that will be equivalent to multiplying by a uniform.
276 //
277 // Alternatively, we could explicitly use a gamma distribution:
278 // `w->factory("Gamma::gamma(b,sum::temp(y,1),1,0)");`
279 //
280 // or we can use some other ad hoc prior that do not naturally
281 // follow from the known form of the auxiliary measurement.
282 // The common choice is the equivalent Gaussian:
283 w->factory("Gaussian::gauss_prior(b,y, expr::sqrty('sqrt(y)',y))");
284 // this corresponds to the "Z_N" calculation.
285 //
286 // or one could use the analogous log-normal prior
287 w->factory("Lognormal::lognorm_prior(b,y, expr::kappa('1+1./sqrt(y)',y))");
288 //
289 // Ideally, the HybridCalculator would be able to inspect the full
290 // model $Pois(x | s+b) * Pois(y | tau b )$ and be given the original
291 // prior $\eta(b)$ to form $\pi(b) = Pois(y|tau*b) * \eta(b)$.
292 // This is not yet implemented because in the general case
293 // it is not easy to identify the terms in the PDF that correspond
294 // to the auxiliary measurement. So for now, it must be set
295 // explicitly with:
296 // - ForcePriorNuisanceNull()
297 // - ForcePriorNuisanceAlt()
298 // the name "ForcePriorNuisance" was chosen because we anticipate
299 // this to be auto-detected, but will leave the option open
300 // to force to a different prior for the nuisance parameters.
301
302 // Part 3d : Construct and configure the HybridCalculator
303 // -------------------------------------------------------
304
306 ToyMCSampler *toymcs1 = (ToyMCSampler *)hc1.GetTestStatSampler();
307 toymcs1->SetNEventsPerToy(1); // because the model is in number counting form
308 toymcs1->SetTestStatistic(&binCount); // set the test statistic
309 hc1.SetToys(ntoys, ntoys / nToysRatio);
310 hc1.ForcePriorNuisanceAlt(*w->pdf("py"));
311 hc1.ForcePriorNuisanceNull(*w->pdf("py"));
312 // if you wanted to use the ad hoc Gaussian prior instead
313 // ~~~
314 // hc1.ForcePriorNuisanceAlt(*w->pdf("gauss_prior"));
315 // hc1.ForcePriorNuisanceNull(*w->pdf("gauss_prior"));
316 // ~~~
317 // if you wanted to use the ad hoc log-normal prior instead
318 // ~~~
319 // hc1.ForcePriorNuisanceAlt(*w->pdf("lognorm_prior"));
320 // hc1.ForcePriorNuisanceNull(*w->pdf("lognorm_prior"));
321 // ~~~
322
323 // these lines save current msg level and then kill any messages below ERROR
324 RooMsgService::instance().setGlobalKillBelow(RooFit::ERROR);
325 // Get the result
326 HypoTestResult *r1 = hc1.GetHypoTest();
327 RooMsgService::instance().setGlobalKillBelow(msglevel); // set it back
328 cout << "-----------------------------------------" << endl;
329 cout << "Part 4" << endl;
330 r1->Print();
331 t.Stop();
332 t.Print();
333 t.Reset();
334 t.Start();
335
336 c->cd(2);
337 HypoTestPlot *p1 = new HypoTestPlot(*r1, 30); // 30 bins, TS is discrete
338 p1->Draw();
339
340 // -------------------------------------------------------------------------
341 // # P A R T 5 : U S I N G H Y B R I D C A L C U L A T O R
342 // # W I T H A N A L T E R N A T I V E T E S T S T A T I S T I C
343 //
344 // A likelihood ratio test statistics should be 1-to-1 with the count x
345 // when the value of b is fixed in the likelihood. This is implemented
346 // by the SimpleLikelihoodRatioTestStat
347
349 slrts.SetNullParameters(*b_model.GetSnapshot());
350 slrts.SetAltParameters(*sb_model.GetSnapshot());
351
352 // HYBRID CALCULATOR
354 ToyMCSampler *toymcs2 = (ToyMCSampler *)hc2.GetTestStatSampler();
355 toymcs2->SetNEventsPerToy(1);
356 toymcs2->SetTestStatistic(&slrts);
357 hc2.SetToys(ntoys, ntoys / nToysRatio);
358 hc2.ForcePriorNuisanceAlt(*w->pdf("py"));
359 hc2.ForcePriorNuisanceNull(*w->pdf("py"));
360 // if you wanted to use the ad hoc Gaussian prior instead
361 // ~~~
362 // hc2.ForcePriorNuisanceAlt(*w->pdf("gauss_prior"));
363 // hc2.ForcePriorNuisanceNull(*w->pdf("gauss_prior"));
364 // ~~~
365 // if you wanted to use the ad hoc log-normal prior instead
366 // ~~~
367 // hc2.ForcePriorNuisanceAlt(*w->pdf("lognorm_prior"));
368 // hc2.ForcePriorNuisanceNull(*w->pdf("lognorm_prior"));
369 // ~~~
370
371 // these lines save current msg level and then kill any messages below ERROR
372 RooMsgService::instance().setGlobalKillBelow(RooFit::ERROR);
373 // Get the result
374 HypoTestResult *r2 = hc2.GetHypoTest();
375 cout << "-----------------------------------------" << endl;
376 cout << "Part 5" << endl;
377 r2->Print();
378 t.Stop();
379 t.Print();
380 t.Reset();
381 t.Start();
382 RooMsgService::instance().setGlobalKillBelow(msglevel);
383
384 c->cd(3);
385 HypoTestPlot *p2 = new HypoTestPlot(*r2, 30); // 30 bins
386 p2->Draw();
387
388 // -----------------------------------------------------------------------------
389 // # P A R T 6 : U S I N G H Y B R I D C A L C U L A T O R W I T H A N A L T E R N A T I V E T E S T
390 // # S T A T I S T I C A N D S I M U L T A N E O U S M O D E L
391 //
392 // If one wants to use a test statistic in which the nuisance parameters
393 // are profiled (in one way or another), then the PDF must constrain b.
394 // Otherwise any observation x can always be explained with s=0 and b=x/tau.
395 //
396 // In this case, one is really thinking about the problem in a
397 // different way. They are considering x,y simultaneously.
398 // and the PDF should be $Pois(x | s+b) * Pois(y | tau b )$
399 // and the set 'obs' should be {x,y}.
400
401 w->defineSet("obsXY", "x,y");
402
403 // create a toy dataset with the x=150, y=100
404 w->var("x")->setVal(150.);
405 w->var("y")->setVal(100.);
406 RooDataSet *dataXY = new RooDataSet("dXY", "dXY", *w->set("obsXY"));
407 dataXY->add(*w->set("obsXY"));
408
409 // now we need new model configs, with PDF="model"
410 ModelConfig b_modelXY("B_modelXY", w);
411 b_modelXY.SetPdf(*w->pdf("model")); // IMPORTANT
412 b_modelXY.SetObservables(*w->set("obsXY"));
413 b_modelXY.SetParametersOfInterest(*w->set("poi"));
414 w->var("s")->setVal(0.0); // IMPORTANT
415 b_modelXY.SetSnapshot(*w->set("poi"));
416
417 // create the alternate (signal+background) ModelConfig with s=50
418 ModelConfig sb_modelXY("S+B_modelXY", w);
419 sb_modelXY.SetPdf(*w->pdf("model")); // IMPORTANT
420 sb_modelXY.SetObservables(*w->set("obsXY"));
421 sb_modelXY.SetParametersOfInterest(*w->set("poi"));
422 w->var("s")->setVal(50.0); // IMPORTANT
423 sb_modelXY.SetSnapshot(*w->set("poi"));
424
425 // Test statistics like the profile likelihood ratio
426 // (or the ratio of profiled likelihoods (Tevatron) or the MLE for s)
427 // will now work, since the nuisance parameter b is constrained by y.
428 // ratio of alt and null likelihoods with background yield profiled.
429 //
430 // NOTE: These are slower because they have to run fits for each toy
431
432 // Tevatron-style Ratio of profiled likelihoods
433 // $Q_Tev = -log L(s=0,\hat\hat{b})/L(s=50,\hat\hat{b})$
434 RatioOfProfiledLikelihoodsTestStat ropl(*b_modelXY.GetPdf(), *sb_modelXY.GetPdf(), sb_modelXY.GetSnapshot());
435 ropl.SetSubtractMLE(false);
436
437 // profile likelihood where alternate is best fit value of signal yield
438 // $\lambda(0) = -log L(s=0,\hat\hat{b})/L(\hat{s},\hat{b})$
440
441 // just use the maximum likelihood estimate of signal yield
442 // $MLE = \hat{s}$
443 MaxLikelihoodEstimateTestStat mlets(*sb_modelXY.GetPdf(), *w->var("s"));
444
445 // However, it is less clear how to justify the prior used in randomizing
446 // the nuisance parameters (since that is a property of the ensemble,
447 // and y is a property of each toy pseudo experiment. In that case,
448 // one probably wants to consider a different y0 which will be held
449 // constant and the prior $\pi(b) = Pois(y0 | tau b) * \eta(b)$.
450 w->factory("y0[100]");
451 w->factory("Gamma::gamma_y0(b,sum::temp0(y0,1),1,0)");
452 w->factory("Gaussian::gauss_prior_y0(b,y0, expr::sqrty0('sqrt(y0)',y0))");
453
454 // HYBRID CALCULATOR
456 ToyMCSampler *toymcs3 = (ToyMCSampler *)hc3.GetTestStatSampler();
457 toymcs3->SetNEventsPerToy(1);
458 toymcs3->SetTestStatistic(&slrts);
459 hc3.SetToys(ntoys, ntoys / nToysRatio);
460 hc3.ForcePriorNuisanceAlt(*w->pdf("gamma_y0"));
461 hc3.ForcePriorNuisanceNull(*w->pdf("gamma_y0"));
462 // if you wanted to use the ad hoc Gaussian prior instead
463 // ~~~{.cpp}
464 // hc3.ForcePriorNuisanceAlt(*w->pdf("gauss_prior_y0"));
465 // hc3.ForcePriorNuisanceNull(*w->pdf("gauss_prior_y0"));
466 // ~~~
467
468 // choose fit-based test statistic
469 toymcs3->SetTestStatistic(&profll);
470 // toymcs3->SetTestStatistic(&ropl);
471 // toymcs3->SetTestStatistic(&mlets);
472
473 // these lines save current msg level and then kill any messages below ERROR
474 RooMsgService::instance().setGlobalKillBelow(RooFit::ERROR);
475 // Get the result
476 HypoTestResult *r3 = hc3.GetHypoTest();
477 cout << "-----------------------------------------" << endl;
478 cout << "Part 6" << endl;
479 r3->Print();
480 t.Stop();
481 t.Print();
482 t.Reset();
483 t.Start();
484 RooMsgService::instance().setGlobalKillBelow(msglevel);
485
486 c->cd(4);
487 c->GetPad(4)->SetLogy();
488 HypoTestPlot *p3 = new HypoTestPlot(*r3, 50); // 50 bins
489 p3->Draw();
490
491 c->SaveAs("zbi.pdf");
492
493 // -----------------------------------------
494 // OUTPUT (2.66 GHz Intel Core i7)
495 // =========================================
496
497 // -----------------------------------------
498 // Part 2
499 // Hybrid p-value from direct integration = 0.00094165
500 // Z_Gamma Significance = 3.10804
501 // -----------------------------------------
502 //
503 // Part 3
504 // Z_Bi p-value (analytic): 0.00094165
505 // Z_Bi significance (analytic): 3.10804
506 // Real time 0:00:00, CP time 0.610
507
508 // -----------------------------------------
509 // Part 4
510 // Results HybridCalculator_result:
511 // - Null p-value = 0.00115 +/- 0.000228984
512 // - Significance = 3.04848 sigma
513 // - Number of S+B toys: 1000
514 // - Number of B toys: 20000
515 // - Test statistic evaluated on data: 150
516 // - CL_b: 0.99885 +/- 0.000239654
517 // - CL_s+b: 0.476 +/- 0.0157932
518 // - CL_s: 0.476548 +/- 0.0158118
519 // Real time 0:00:07, CP time 7.620
520
521 // -----------------------------------------
522 // Part 5
523 // Results HybridCalculator_result:
524 // - Null p-value = 0.0009 +/- 0.000206057
525 // - Significance = 3.12139 sigma
526 // - Number of S+B toys: 1000
527 // - Number of B toys: 20000
528 // - Test statistic evaluated on data: 10.8198
529 // - CL_b: 0.9991 +/- 0.000212037
530 // - CL_s+b: 0.465 +/- 0.0157726
531 // - CL_s: 0.465419 +/- 0.0157871
532 // Real time 0:00:34, CP time 34.360
533
534 // -----------------------------------------
535 // Part 6
536 // Results HybridCalculator_result:
537 // - Null p-value = 0.000666667 +/- 0.000149021
538 // - Significance = 3.20871 sigma
539 // - Number of S+B toys: 1000
540 // - Number of B toys: 30000
541 // - Test statistic evaluated on data: 5.03388
542 // - CL_b: 0.999333 +/- 0.000149021
543 // - CL_s+b: 0.511 +/- 0.0158076
544 // - CL_s: 0.511341 +/- 0.0158183
545 // Real time 0:05:06, CP time 306.330
546
547 // ----------------------------------
548 // Comparison
549 // ----------------------------------
550 // LEPStatToolsForLHC
551 // https://plone4.fnal.gov:4430/P0/phystat/packages/0703002
552 // Uses Gaussian prior
553 // CL_b = 6.218476e-04, Significance = 3.228665 sigma
554 //
555 // ----------------------------------
556 // Comparison
557 // ----------------------------------
558 // Asymptotic
559 // From the value of the profile likelihood ratio (5.0338)
560 // The significance can be estimated using Wilks's theorem
561 // significance = sqrt(2*profileLR) = 3.1729 sigma
562}
#define c(i)
Definition RSha256.hxx:101
double Double_t
Double 8 bytes.
Definition RtypesCore.h:73
#define ClassDef(name, id)
Definition Rtypes.h:344
@ kRed
Definition Rtypes.h:67
@ kGreen
Definition Rtypes.h:67
@ kBlue
Definition Rtypes.h:67
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
Abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
Container class to hold unbinned data.
Definition RooDataSet.h:32
static RooMsgService & instance()
Return reference to singleton instance.
Plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:43
void Draw(Option_t *options=nullptr) override
Draw this plot and all of the elements it contains.
Definition RooPlot.cxx:596
Same purpose as HybridCalculatorOriginal, but different implementation.
This class provides the plots for the result of a study performed with any of the HypoTestCalculatorG...
HypoTestResult is a base class for results from hypothesis tests.
MaxLikelihoodEstimateTestStat: TestStatistic that returns maximum likelihood estimate of a specified ...
< A class that holds configuration information for a model using a workspace as a store
Definition ModelConfig.h:34
ProfileLikelihoodTestStat is an implementation of the TestStatistic interface that calculates the pro...
TestStatistic that returns the ratio of profiled likelihoods.
TestStatistic class that returns -log(L[null] / L[alt]) where L is the likelihood.
TestStatistic is an interface class to provide a facility for construction test statistics distributi...
ToyMCSampler is an implementation of the TestStatSampler interface.
Persistable container for RooFit projects.
The Canvas class.
Definition TCanvas.h:23
void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0) override
Automatic pad generation by division.
Definition TPad.cxx:1258
Stopwatch class.
Definition TStopwatch.h:28
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
void Stop()
Stop the stopwatch.
void Reset()
Definition TStopwatch.h:52
void Print(Option_t *option="") const override
Print the real and cpu time passed between the start and stop events.
Basic string class.
Definition TString.h:138
RooCmdArg LineColor(TColorNumber color)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition CodegenImpl.h:67
MsgLevel
Verbosity level for RooMsgService::StreamConfig in RooMsgService.
Namespace for the RooStats classes.
Definition CodegenImpl.h:61
double PValueToSignificance(double pvalue)
returns one-sided significance corresponding to a p-value
Ta Range(0, 0, 1, 1)