Logo ROOT   6.18/05
Reference Guide
HypoTestInverterResult.h
Go to the documentation of this file.
1// @(#)root/roostats:$Id$
2// Author: Kyle Cranmer, Lorenzo Moneta, Gregory Schott, Wouter Verkerke
3/*************************************************************************
4 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
5 * All rights reserved. *
6 * *
7 * For the licensing terms see $ROOTSYS/LICENSE. *
8 * For the list of contributors see $ROOTSYS/README/CREDITS. *
9 *************************************************************************/
10
11#ifndef ROOSTATS_HypoTestInverterResult
12#define ROOSTATS_HypoTestInverterResult
13
14
15
17
19
20class RooRealVar;
21
22namespace RooStats {
23
24class SamplingDistribution;
25
27
28public:
29
30 /// default constructor
31 explicit HypoTestInverterResult(const char* name = 0);
32
33 /// constructor
34 HypoTestInverterResult( const char* name,
35 const RooRealVar& scannedVariable,
36 double cl ) ;
37
38 HypoTestInverterResult( const HypoTestInverterResult& other, const char* name );
39
40 /// destructor
42
43 /// operator =
45
46 /// remove points that appear to have failed.
47 int ExclusionCleanup();
48
49 /// merge with the content of another HypoTestInverterResult object
50 bool Add( const HypoTestInverterResult& otherResult );
51
52 ///add the result of a single point (an HypoTestRsult)
53 bool Add( Double_t x, const HypoTestResult & result );
54
55 /// function to return the value of the parameter of interest for the i^th entry in the results
56 double GetXValue( int index ) const ;
57
58 /// function to return the value of the confidence level for the i^th entry in the results
59 double GetYValue( int index ) const ;
60
61 /// function to return the estimated error on the value of the confidence level for the i^th entry in the results
62 double GetYError( int index ) const ;
63
64 /// return the observed CLsplusb value for the i-th entry
65 double CLsplusb( int index) const;
66
67 /// return the observed CLb value for the i-th entry
68 double CLb( int index) const;
69
70 /// return the observed CLb value for the i-th entry
71 double CLs( int index) const;
72
73 /// return the observed CLsplusb value for the i-th entry
74 double CLsplusbError( int index) const;
75
76 /// return the observed CLb value for the i-th entry
77 double CLbError( int index) const;
78
79 /// return the observed CLb value for the i-th entry
80 double CLsError( int index) const;
81
82 /// return a pointer to the i^th result object
83 HypoTestResult* GetResult( int index ) const ;
84
85 double GetLastYValue( ) const { return GetYValue( fXValues.size()-1); }
86
87 double GetLastXValue( ) const { return GetXValue( fXValues.size()-1); }
88
89 double GetLastYError( ) const { return GetYError( fXValues.size()-1); }
90
91 HypoTestResult * GetLastResult( ) const { return GetResult( fXValues.size()-1); }
92
93 /// number of entries in the results array
94 int ArraySize() const { return fXValues.size(); };
95
96 int FindIndex(double xvalue) const;
97
98 /// set the size of the test (rate of Type I error) (eg. 0.05 for a 95% Confidence Interval)
99 virtual void SetTestSize( Double_t size ) { fConfidenceLevel = 1.-size; }
100
101 /// set the confidence level for the interval (eg. 0.95 for a 95% Confidence Interval)
102 virtual void SetConfidenceLevel( Double_t cl ) { fConfidenceLevel = cl; }
103
104 /// set CLs threshold for exclusion cleanup function
106
107 /// flag to switch between using CLsb (default) or CLs as confidence level
108 void UseCLs( bool on = true ) { fUseCLs = on; }
109
110 /// query if one sided result
111 bool IsOneSided() const { return !fIsTwoSided; }
112 /// query if two sided result
113 bool IsTwoSided() const { return fIsTwoSided; }
114
115 /// lower and upper bound of the confidence interval (to get upper/lower limits, multiply the size( = 1-confidence level ) by 2
118
119 /// rough estimation of the error on the computed bound of the confidence interval
120 /// Estimate of lower limit error
121 ///function evaluates only a rough error on the lower limit. Be careful when using this estimation
123
124 /// Estimate of lower limit error
125 ///function evaluates only a rough error on the lower limit. Be careful when using this estimation
127
128 /// return expected distribution of p-values (Cls or Clsplusb)
129
131
133
135
136 /// same in terms of alt and null
139 }
141 return GetBackgroundTestStatDist(index);
142 }
143
144 /// get expected lower limit distributions
145 /// implemented using interpolation
146 /// The size for the sampling distribution is given (by default is given by the average number of toy/point)
148
149 /// get expected upper limit distributions
150 /// implemented using interpolation
152
153 /// get Limit value corresponding at the desired nsigma level (0) is median -1 sigma is 1 sigma
154 double GetExpectedLowerLimit(double nsig = 0, const char * opt = "" ) const ;
155
156 /// get Limit value corresponding at the desired nsigma level (0) is median -1 sigma is 1 sigma
157 double GetExpectedUpperLimit(double nsig = 0, const char * opt = "") const ;
158
159
160 double FindInterpolatedLimit(double target, bool lowSearch = false, double xmin=1, double xmax=0);
161
163
164 /// set the interpolation option, linear (kLinear ) or spline (kSpline)
166
168
169private:
170
171
172 double CalculateEstimatedError(double target, bool lower = true, double xmin = 1, double xmax = 0);
173
174 int FindClosestPointIndex(double target, int mode = 0, double xtarget = 0);
175
176 SamplingDistribution* GetLimitDistribution(bool lower ) const;
177
178 double GetExpectedLimit(double nsig, bool lower, const char * opt = "" ) const ;
179
180 double GetGraphX(const TGraph & g, double y0, bool lowSearch, double &xmin, double &xmax) const;
181 double GetGraphX(const TGraph & g, double y0, bool lowSearch = true) const {
182 double xmin=1; double xmax = 0;
183 return GetGraphX(g,y0,lowSearch,xmin,xmax);
184 }
185
186
187protected:
188
190 bool fIsTwoSided; /// two sided scan (look for lower/upper limit)
195 InterpolOption_t fInterpolOption; /// interpolation option (linear or spline)
196
199
201
202 static double fgAsymptoticMaxSigma; /// max sigma value used to scan asymptotic expected p values
203 static int fgAsymptoticNumPoints; /// number of points used to build expected p-values
204
205 std::vector<double> fXValues;
206
207 TList fYObjects; /// list of HypoTestResult for each point
208 TList fExpPValues; /// list of expected sampling distribution for each point
209
210 friend class HypoTestInverter;
213
214 ClassDef(HypoTestInverterResult,5) /// HypoTestInverterResult class
215};
216}
217
218#endif
#define g(i)
Definition: RSha256.hxx:105
double Double_t
Definition: RtypesCore.h:55
#define ClassDef(name, id)
Definition: Rtypes.h:326
char name[80]
Definition: TGX11.cxx:109
float xmin
Definition: THbookFile.cxx:93
float xmax
Definition: THbookFile.cxx:93
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
This class is now deprecated and to be replaced by the HypoTestInverter.
Class to plot an HypoTestInverterResult, result of the HypoTestInverter calculator.
HypoTestInverterResult class holds the array of hypothesis test results and compute a confidence inte...
double GetGraphX(const TGraph &g, double y0, bool lowSearch=true) const
virtual void SetTestSize(Double_t size)
set the size of the test (rate of Type I error) (eg. 0.05 for a 95% Confidence Interval)
void SetInterpolationOption(InterpolOption_t opt)
set the interpolation option, linear (kLinear ) or spline (kSpline)
Double_t LowerLimit()
lower and upper bound of the confidence interval (to get upper/lower limits, multiply the size( = 1-c...
Double_t LowerLimitEstimatedError()
rough estimation of the error on the computed bound of the confidence interval Estimate of lower limi...
double GetYValue(int index) const
function to return the value of the confidence level for the i^th entry in the results
SamplingDistribution * GetLowerLimitDistribution() const
get expected lower limit distributions implemented using interpolation The size for the sampling dist...
SamplingDistribution * GetSignalAndBackgroundTestStatDist(int index) const
get the signal and background test statistic distribution
Double_t UpperLimitEstimatedError()
Estimate of lower limit error function evaluates only a rough error on the lower limit.
HypoTestResult * GetResult(int index) const
return a pointer to the i^th result object
double fLowerLimitError
interpolation option (linear or spline)
HypoTestInverterResult & operator=(const HypoTestInverterResult &other)
operator =
double CalculateEstimatedError(double target, bool lower=true, double xmin=1, double xmax=0)
Return an error estimate on the upper(lower) limit.
HypoTestInverterResult(const char *name=0)
default constructor
int ArraySize() const
number of entries in the results array
double FindInterpolatedLimit(double target, bool lowSearch=false, double xmin=1, double xmax=0)
interpolate to find a limit value Use a linear or a spline interpolation depending on the interpolati...
bool fInterpolateLowerLimit
two sided scan (look for lower/upper limit)
virtual void SetConfidenceLevel(Double_t cl)
set the confidence level for the interval (eg. 0.95 for a 95% Confidence Interval)
std::vector< double > fXValues
number of points used to build expected p-values
double CLsplusbError(int index) const
return the observed CLsplusb value for the i-th entry
double CLsError(int index) const
return the observed CLb value for the i-th entry
double GetGraphX(const TGraph &g, double y0, bool lowSearch, double &xmin, double &xmax) const
return the X value of the given graph for the target value y0 the graph is evaluated using linear int...
double CLbError(int index) const
return the observed CLb value for the i-th entry
void UseCLs(bool on=true)
flag to switch between using CLsb (default) or CLs as confidence level
double GetExpectedUpperLimit(double nsig=0, const char *opt="") const
get Limit value corresponding at the desired nsigma level (0) is median -1 sigma is 1 sigma
InterpolOption_t GetInterpolationOption() const
double GetXValue(int index) const
function to return the value of the parameter of interest for the i^th entry in the results
bool IsOneSided() const
query if one sided result
bool Add(const HypoTestInverterResult &otherResult)
merge with the content of another HypoTestInverterResult object
int FindClosestPointIndex(double target, int mode=0, double xtarget=0)
TList fExpPValues
list of HypoTestResult for each point
double GetExpectedLowerLimit(double nsig=0, const char *opt="") const
get Limit value corresponding at the desired nsigma level (0) is median -1 sigma is 1 sigma
HypoTestResult * GetLastResult() const
static int fgAsymptoticNumPoints
max sigma value used to scan asymptotic expected p values
int ExclusionCleanup()
remove points that appear to have failed.
double CLs(int index) const
return the observed CLb value for the i-th entry
int FindIndex(double xvalue) const
find the index corresponding at the poi value xvalue If no points is found return -1 Note that a tole...
double GetYError(int index) const
function to return the estimated error on the value of the confidence level for the i^th entry in the...
void SetCLsCleanupThreshold(Double_t th)
set CLs threshold for exclusion cleanup function
double CLb(int index) const
return the observed CLb value for the i-th entry
double GetExpectedLimit(double nsig, bool lower, const char *opt="") const
get expected limit (lower/upper) depending on the flag for asymptotic is a special case (the distribu...
SamplingDistribution * GetExpectedPValueDist(int index) const
return expected distribution of p-values (Cls or Clsplusb)
double CLsplusb(int index) const
return the observed CLsplusb value for the i-th entry
SamplingDistribution * GetLimitDistribution(bool lower) const
get the limit distribution (lower/upper depending on the flag) by interpolating the expected p values...
SamplingDistribution * GetUpperLimitDistribution() const
get expected upper limit distributions implemented using interpolation
SamplingDistribution * GetNullTestStatDist(int index) const
same in terms of alt and null
SamplingDistribution * GetAltTestStatDist(int index) const
SamplingDistribution * GetBackgroundTestStatDist(int index) const
get the background test statistic distribution
bool IsTwoSided() const
query if two sided result
HypoTestInverter class for performing an hypothesis test inversion by scanning the hypothesis test re...
HypoTestResult is a base class for results from hypothesis tests.
This class simply holds a sampling distribution of some test statistic.
SimpleInterval is a concrete implementation of the ConfInterval interface.
A Graph is a graphics object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
A doubly linked list.
Definition: TList.h:44
Double_t x[n]
Definition: legend1.C:17
Namespace for the RooStats classes.
Definition: Asimov.h:20