Logo ROOT   6.08/07
Reference Guide
NumberCountingUtils.cxx
Go to the documentation of this file.
1 // @(#)root/roostats:$Id$
2 // Author: Kyle Cranmer 28/07/2008
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 
13 #ifndef RooStats_NumberCountingUtils
15 #endif
16 
17 #ifndef RooStats_RooStatsUtils
18 #include "RooStats/RooStatsUtils.h"
19 #endif
20 
21 // // Without this macro the THtml doc can not be generated
22 // #if !defined(R__ALPHA) && !defined(R__SOLARIS) && !defined(R__ACC) && !defined(R__FBSD)
23 // NamespaceImp(RooStats)
24 // //NamespaceImp(NumberCountingUtils)
25 // #endif
26 
27 //using namespace RooStats;
28 
29 Double_t RooStats::NumberCountingUtils::BinomialExpP(Double_t signalExp, Double_t backgroundExp, Double_t relativeBkgUncert){
30  // Expected P-value for s=0 in a ratio of Poisson means.
31  // Here the background and its uncertainty are provided directly and
32  // assumed to be from the double Poisson counting setup described in the
33  // BinomialWithTau functions.
34  // Normally one would know tau directly, but here it is determiend from
35  // the background uncertainty. This is not strictly correct, but a useful
36  // approximation.
37 
38 
39  //SIDE BAND EXAMPLE
40  //See Eqn. (19) of Cranmer and pp. 36-37 of Linnemann.
41  //150 total events in signalExp region, 100 in sideband of equal size
42  Double_t mainInf = signalExp+backgroundExp; //Given
43  Double_t tau = 1./backgroundExp/(relativeBkgUncert*relativeBkgUncert);
44  Double_t auxiliaryInf = backgroundExp*tau; //Given
45 
46  Double_t P_Bi = TMath::BetaIncomplete(1./(1.+tau),mainInf,auxiliaryInf+1);
47  return P_Bi;
48 
49 /*
50 Now, if instead the mean background level b in the signal region is
51 specified, along with Gaussian rms sigb, then one can fake a Poisson
52 sideband (see Linnemann, p. 35, converted to Cranmer's notation) by
53 letting tau = b/(sigb*sigb) and y = b*tau. Thus, for example, if one
54 has x=150 and b = 100 +/- 10, one then derives tau and y. Then one
55 has the same two lines of ROOT calling BetaIncomplete and ErfInverse.
56 Since I chose these numbers to revert to the previous example, we get
57 the same answer:
58 */
59 /*
60 //GAUSSIAN FAKED AS POISSON EXAMPLE
61 x = 150. //Given
62 b = 100. //Given
63 sigb = 10. //Given
64 tau = b/(sigb*sigb)
65 y = tau*b
66 Z_Bi = TMath::BetaIncomplete(1./(1.+tau),x,y+1)
67 S = sqrt(2)*TMath::ErfInverse(1 - 2*Z_Bi)
68 
69 */
70 
71 }
72 
73 
75  // Expected P-value for s=0 in a ratio of Poisson means.
76  // Based on two expectations, a main measurement that might have signal
77  // and an auxiliarly measurement for the background that is signal free.
78  // The expected background in the auxiliary measurement is a factor
79  // tau larger than in the main measurement.
80 
81  Double_t mainInf = signalExp+backgroundExp; //Given
82  Double_t auxiliaryInf = backgroundExp*tau; //Given
83 
84  Double_t P_Bi = TMath::BetaIncomplete(1./(1.+tau),mainInf,auxiliaryInf+1);
85 
86  return P_Bi;
87 
88 }
89 
91  // P-value for s=0 in a ratio of Poisson means.
92  // Here the background and its uncertainty are provided directly and
93  // assumed to be from the double Poisson counting setup.
94  // Normally one would know tau directly, but here it is determiend from
95  // the background uncertainty. This is not strictly correct, but a useful
96  // approximation.
97 
98  Double_t tau = 1./backgroundObs/(relativeBkgUncert*relativeBkgUncert);
99  Double_t auxiliaryInf = backgroundObs*tau; //Given
100 
101 
102  //SIDE BAND EXAMPLE
103  //See Eqn. (19) of Cranmer and pp. 36-37 of Linnemann.
104  Double_t P_Bi = TMath::BetaIncomplete(1./(1.+tau),mainObs,auxiliaryInf+1);
105 
106  return P_Bi;
107 
108 }
109 
110 
112  // P-value for s=0 in a ratio of Poisson means.
113  // Based on two observations, a main measurement that might have signal
114  // and an auxiliarly measurement for the background that is signal free.
115  // The expected background in the auxiliary measurement is a factor
116  // tau larger than in the main measurement.
117 
118  //SIDE BAND EXAMPLE
119  //See Eqn. (19) of Cranmer and pp. 36-37 of Linnemann.
120  Double_t P_Bi = TMath::BetaIncomplete(1./(1.+tau),mainObs,auxiliaryObs+1);
121 
122  return P_Bi;
123 
124 }
125 
126 Double_t RooStats::NumberCountingUtils::BinomialExpZ(Double_t signalExp, Double_t backgroundExp, Double_t relativeBkgUncert) {
127  // See BinomialExpP
128  return RooStats::PValueToSignificance( BinomialExpP(signalExp,backgroundExp,relativeBkgUncert) ) ;
129  }
130 
132  // See BinomialWithTauExpP
133  return RooStats::PValueToSignificance( BinomialWithTauExpP(signalExp,backgroundExp,tau) ) ;
134 }
135 
136 
138  // See BinomialObsP
139  return RooStats::PValueToSignificance( BinomialObsP(mainObs,backgroundObs,relativeBkgUncert) ) ;
140 }
141 
143  // See BinomialWithTauObsP
144  return RooStats::PValueToSignificance( BinomialWithTauObsP(mainObs,auxiliaryObs,tau) ) ;
145 }
Double_t PValueToSignificance(Double_t pvalue)
Definition: RooStatsUtils.h:59
Double_t BinomialWithTauObsP(Double_t nObs, Double_t bExp, Double_t tau)
Double_t BinomialObsP(Double_t nObs, Double_t, Double_t fractionalBUncertainty)
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t tau
Definition: TRolke.cxx:630
Double_t BinomialWithTauObsZ(Double_t nObs, Double_t bExp, Double_t tau)
Double_t BinomialExpZ(Double_t sExp, Double_t bExp, Double_t fractionalBUncertainty)
double Double_t
Definition: RtypesCore.h:55
Double_t BinomialExpP(Double_t sExp, Double_t bExp, Double_t fractionalBUncertainty)
Double_t BinomialWithTauExpP(Double_t sExp, Double_t bExp, Double_t tau)
Double_t BetaIncomplete(Double_t x, Double_t a, Double_t b)
Calculates the incomplete Beta-function.
Definition: TMath.cxx:2068
Double_t BinomialObsZ(Double_t nObs, Double_t bExp, Double_t fractionalBUncertainty)
Double_t BinomialWithTauExpZ(Double_t sExp, Double_t bExp, Double_t tau)