Logo ROOT   6.08/07
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 
16 #ifndef ROOSTATS_SimpleInterval
18 #endif
19 
21 
22 class RooRealVar;
23 
24 namespace RooStats {
25 
26 class SamplingDistribution;
27 
28 /**
29  HypoTestInverterResult class: holds the array of hypothesis test results and compute a confidence interval.
30  Based on the RatioFinder code available in the RooStatsCms package developed by Gregory Schott and Danilo Piparo
31  Ported and adapted to RooStats by Gregory Schott
32  Some contributions to this class have been written by Matthias Wolf (error estimation)
33 
34  \ingroup Roostats
35 **/
36 
37 
39 
40 public:
41 
42  /// default constructor
43  explicit HypoTestInverterResult(const char* name = 0);
44 
45  /// constructor
46  HypoTestInverterResult( const char* name,
47  const RooRealVar& scannedVariable,
48  double cl ) ;
49 
50  HypoTestInverterResult( const HypoTestInverterResult& other, const char* name );
51 
52  /// destructor
53  virtual ~HypoTestInverterResult();
54 
55  /// operator =
57 
58  /// remove points that appear to have failed.
59  int ExclusionCleanup();
60 
61  /// merge with the content of another HypoTestInverterResult object
62  bool Add( const HypoTestInverterResult& otherResult );
63 
64  ///add the result of a single point (an HypoTestRsult)
65  bool Add( Double_t x, const HypoTestResult & result );
66 
67  /// function to return the value of the parameter of interest for the i^th entry in the results
68  double GetXValue( int index ) const ;
69 
70  /// function to return the value of the confidence level for the i^th entry in the results
71  double GetYValue( int index ) const ;
72 
73  /// function to return the estimated error on the value of the confidence level for the i^th entry in the results
74  double GetYError( int index ) const ;
75 
76  /// return the observed CLsplusb value for the i-th entry
77  double CLsplusb( int index) const;
78 
79  /// return the observed CLb value for the i-th entry
80  double CLb( int index) const;
81 
82  /// return the observed CLb value for the i-th entry
83  double CLs( int index) const;
84 
85  /// return the observed CLsplusb value for the i-th entry
86  double CLsplusbError( int index) const;
87 
88  /// return the observed CLb value for the i-th entry
89  double CLbError( int index) const;
90 
91  /// return the observed CLb value for the i-th entry
92  double CLsError( int index) const;
93 
94  /// return a pointer to the i^th result object
95  HypoTestResult* GetResult( int index ) const ;
96 
97  double GetLastYValue( ) const { return GetYValue( fXValues.size()-1); }
98 
99  double GetLastXValue( ) const { return GetXValue( fXValues.size()-1); }
100 
101  double GetLastYError( ) const { return GetYError( fXValues.size()-1); }
102 
103  HypoTestResult * GetLastResult( ) const { return GetResult( fXValues.size()-1); }
104 
105  /// number of entries in the results array
106  int ArraySize() const { return fXValues.size(); };
107 
108  int FindIndex(double xvalue) const;
109 
110  /// set the size of the test (rate of Type I error) (eg. 0.05 for a 95% Confidence Interval)
111  virtual void SetTestSize( Double_t size ) { fConfidenceLevel = 1.-size; }
112 
113  /// set the confidence level for the interval (eg. 0.95 for a 95% Confidence Interval)
114  virtual void SetConfidenceLevel( Double_t cl ) { fConfidenceLevel = cl; }
115 
116  /// set CLs threshold for exclusion cleanup function
118 
119  /// flag to switch between using CLsb (default) or CLs as confidence level
120  void UseCLs( bool on = true ) { fUseCLs = on; }
121 
122  /// query if one sided result
123  bool IsOneSided() const { return !fIsTwoSided; }
124  /// query if two sided result
125  bool IsTwoSided() const { return fIsTwoSided; }
126 
127  /// lower and upper bound of the confidence interval (to get upper/lower limits, multiply the size( = 1-confidence level ) by 2
130 
131  /// rough estimation of the error on the computed bound of the confidence interval
132  /// Estimate of lower limit error
133  ///function evaluates only a rought error on the lower limit. Be careful when using this estimation
135 
136  /// Estimate of lower limit error
137  ///function evaluates only a rought error on the lower limit. Be careful when using this estimation
139 
140  /// return expected distribution of p-values (Cls or Clsplusb)
141 
142  SamplingDistribution * GetExpectedPValueDist(int index) const;
143 
144  SamplingDistribution * GetBackgroundTestStatDist(int index ) const;
145 
147 
148  /// same in terms of alt and null
151  }
153  return GetBackgroundTestStatDist(index);
154  }
155 
156  /// get expected lower limit distributions
157  /// implemented using interpolation
158  /// The size for the sampling distribution is given (by default is given by the average number of toy/point)
160 
161  /// get expected upper limit distributions
162  /// implemented using interpolation
164 
165  /// get Limit value correspnding at the desired nsigma level (0) is median -1 sigma is 1 sigma
166  double GetExpectedLowerLimit(double nsig = 0, const char * opt = "" ) const ;
167 
168  /// get Limit value correspnding at the desired nsigma level (0) is median -1 sigma is 1 sigma
169  double GetExpectedUpperLimit(double nsig = 0, const char * opt = "") const ;
170 
171 
172  double FindInterpolatedLimit(double target, bool lowSearch = false, double xmin=1, double xmax=0);
173 
175 
176  /// set the interpolation option, linear (kLinear ) or spline (kSpline)
178 
180 
181 private:
182 
183 
184  double CalculateEstimatedError(double target, bool lower = true, double xmin = 1, double xmax = 0);
185 
186  int FindClosestPointIndex(double target, int mode = 0, double xtarget = 0);
187 
188  SamplingDistribution* GetLimitDistribution(bool lower ) const;
189 
190  double GetExpectedLimit(double nsig, bool lower, const char * opt = "" ) const ;
191 
192  double GetGraphX(const TGraph & g, double y0, bool lowSearch, double &xmin, double &xmax) const;
193  double GetGraphX(const TGraph & g, double y0, bool lowSearch = true) const {
194  double xmin=1; double xmax = 0;
195  return GetGraphX(g,y0,lowSearch,xmin,xmax);
196  }
197 
198 
199 protected:
200 
201  bool fUseCLs;
202  bool fIsTwoSided; /// two sided scan (look for lower/upper limit)
207  InterpolOption_t fInterpolOption; /// interpolatation option (linear or spline)
208 
211 
213 
214  static double fgAsymptoticMaxSigma; /// max sigma value used to scan asymptotic expected p values
215  static int fgAsymptoticNumPoints; /// number of points used to build expected p-values
216 
217  std::vector<double> fXValues;
218 
219  TList fYObjects; /// list of HypoTestResult for each point
220  TList fExpPValues; /// list of expected sampling distribution for each point
221 
222  friend class HypoTestInverter;
223  friend class HypoTestInverterPlot;
225 
226  ClassDef(HypoTestInverterResult,5) /// HypoTestInverterResult class
227 };
228 }
229 
230 #endif
double GetExpectedLowerLimit(double nsig=0, const char *opt="") const
get Limit value correspnding at the desired nsigma level (0) is median -1 sigma is 1 sigma ...
SamplingDistribution * GetSignalAndBackgroundTestStatDist(int index) const
float xmin
Definition: THbookFile.cxx:93
SamplingDistribution * GetExpectedPValueDist(int index) const
return expected distribution of p-values (Cls or Clsplusb)
HypoTestInverter class for performing an hypothesis test inversion by scanning the hypothesis test re...
int FindClosestPointIndex(double target, int mode=0, double xtarget=0)
void SetCLsCleanupThreshold(Double_t th)
set CLs threshold for exclusion cleanup function
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) ...
double GetExpectedUpperLimit(double nsig=0, const char *opt="") const
get Limit value correspnding at the desired nsigma level (0) is median -1 sigma is 1 sigma ...
HypoTestInverterResult(const char *name=0)
default constructor
double GetXValue(int index) const
function to return the value of the parameter of interest for the i^th entry in the results ...
HypoTestResult is a base class for results from hypothesis tests.
HypoTestResult * GetResult(int index) const
return a pointer to the i^th result object
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 SetInterpolationOption(InterpolOption_t opt)
set the interpolation option, linear (kLinear ) or spline (kSpline)
Double_t x[n]
Definition: legend1.C:17
#define ClassDef(name, id)
Definition: Rtypes.h:254
double FindInterpolatedLimit(double target, bool lowSearch=false, double xmin=1, double xmax=0)
double CLsError(int index) const
return the observed CLb value for the i-th entry
SamplingDistribution * GetLowerLimitDistribution() const
get expected lower limit distributions implemented using interpolation The size for the sampling dist...
double GetExpectedLimit(double nsig, bool lower, const char *opt="") const
void UseCLs(bool on=true)
flag to switch between using CLsb (default) or CLs as confidence level
double CLsplusb(int index) const
return the observed CLsplusb value for the i-th entry
Double_t LowerLimit()
lower and upper bound of the confidence interval (to get upper/lower limits, multiply the size( = 1-c...
bool fInterpolateLowerLimit
two sided scan (look for lower/upper limit)
int ExclusionCleanup()
remove points that appear to have failed.
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:37
A doubly linked list.
Definition: TList.h:47
double CLs(int index) const
return the observed CLb value for the i-th entry
SamplingDistribution * GetBackgroundTestStatDist(int index) const
double CLbError(int index) const
return the observed CLb value for the i-th entry
double CLb(int index) const
return the observed CLb value for the i-th entry
TList fExpPValues
list of HypoTestResult for each point
SamplingDistribution * GetUpperLimitDistribution() const
get expected upper limit distributions implemented using interpolation
static int fgAsymptoticNumPoints
max sigma value used to scan asymptotic expected p values
float xmax
Definition: THbookFile.cxx:93
This class simply holds a sampling distribution of some test statistic.
bool IsTwoSided() const
query if two sided result
HypoTestResult * GetLastResult() const
InterpolOption_t GetInterpolationOption() const
Namespace for the RooStats classes.
Definition: Asimov.h:20
This class is now depratcated and to be replaced by the HypoTestInverter.
HypoTestInverterResult class: holds the array of hypothesis test results and compute a confidence int...
double Double_t
Definition: RtypesCore.h:55
Class to plot an HypoTestInverterResult, result of the HypoTestInverter calculator.
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 ...
double GetGraphX(const TGraph &g, double y0, bool lowSearch=true) const
Double_t UpperLimitEstimatedError()
Estimate of lower limit error function evaluates only a rought error on the lower limit...
HypoTestInverterResult & operator=(const HypoTestInverterResult &other)
operator =
double fLowerLimitError
interpolatation option (linear or spline)
SamplingDistribution * GetNullTestStatDist(int index) const
same in terms of alt and null
A Graph is a graphics object made of two arrays X and Y with npoints each.
Definition: TGraph.h:53
double GetGraphX(const TGraph &g, double y0, bool lowSearch, double &xmin, double &xmax) const
SamplingDistribution * GetLimitDistribution(bool lower) const
bool Add(const HypoTestInverterResult &otherResult)
merge with the content of another HypoTestInverterResult object
double result[121]
double CalculateEstimatedError(double target, bool lower=true, double xmin=1, double xmax=0)
double CLsplusbError(int index) const
return the observed CLsplusb value for the i-th entry
int ArraySize() const
number of entries in the results array
virtual void SetConfidenceLevel(Double_t cl)
set the confidence level for the interval (eg. 0.95 for a 95% Confidence Interval) ...
char name[80]
Definition: TGX11.cxx:109
bool IsOneSided() const
query if one sided result
SamplingDistribution * GetAltTestStatDist(int index) const
std::vector< double > fXValues
number of points used to build expected p-values