Logo ROOT   6.18/05
Reference Guide
HypoTestCalculatorGeneric.cxx
Go to the documentation of this file.
1// @(#)root/roostats:$Id$
2// Author: Kyle Cranmer, Sven Kreiss 23/05/10
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/** \class RooStats::HypoTestCalculatorGeneric
12 \ingroup Roostats
13
14Common base class for the Hypothesis Test Calculators.
15It is not designed to use directly but via its derived classes
16
17Same purpose as HybridCalculatorOriginal, but different implementation.
18
19This is the "generic" version that works with any TestStatSampler. The
20HybridCalculator derives from this class but explicitly uses the
21ToyMCSampler as its TestStatSampler.
22
23*/
24
30
31#include "RooAddPdf.h"
32
33#include "RooRandom.h"
34
35
37
38using namespace RooStats;
39using namespace std;
40
41////////////////////////////////////////////////////////////////////////////////
42/// Constructor. When test stat sampler is not provided
43/// uses ToyMCSampler and RatioOfProfiledLikelihoodsTestStat
44/// and nToys = 1000.
45/// User can : GetTestStatSampler()->SetNToys( # )
46
47HypoTestCalculatorGeneric::HypoTestCalculatorGeneric(
48 const RooAbsData &data,
49 const ModelConfig &altModel,
50 const ModelConfig &nullModel,
51 TestStatSampler *sampler
52 ) :
53 fAltModel(&altModel),
54 fNullModel(&nullModel),
55 fData(&data),
56 fTestStatSampler(sampler),
57 fDefaultSampler(0),
58 fDefaultTestStat(0),
59 fAltToysSeed(0)
60{
61 if(!sampler){
64 *altModel.GetPdf(),
65 altModel.GetSnapshot());
66
69 }
70
71
72}
73
74////////////////////////////////////////////////////////////////////////////////
75/// common setup for both models
76
81
82 // for this model
83 model.LoadSnapshot();
87 // global observables or nuisance pdf will be set by the derived classes
88 // (e.g. Frequentist or HybridCalculator)
89}
90
91////////////////////////////////////////////////////////////////////////////////
92
96}
97
98////////////////////////////////////////////////////////////////////////////////
99/// several possibilities:
100/// no prior nuisance given and no nuisance parameters: ok
101/// no prior nuisance given but nuisance parameters: error
102/// prior nuisance given for some nuisance parameters:
103/// - nuisance parameters are constant, so they don't float in test statistic
104/// - nuisance parameters are floating, so they do float in test statistic
105
107
108 // initial setup
109 PreHook();
110 const_cast<ModelConfig*>(fNullModel)->GuessObsAndNuisance(*fData);
111 const_cast<ModelConfig*>(fAltModel)->GuessObsAndNuisance(*fData);
112
113 const RooArgSet * nullSnapshot = fNullModel->GetSnapshot();
114 if(nullSnapshot == NULL) {
115 oocoutE((TObject*)0,Generation) << "Null model needs a snapshot. Set using modelconfig->SetSnapshot(poi)." << endl;
116 return 0;
117 }
118
119 // CheckHook
120 if(CheckHook() != 0) {
121 oocoutE((TObject*)0,Generation) << "There was an error in CheckHook(). Stop." << endl;
122 return 0;
123 }
124
126 oocoutE((TObject*)0,InputArguments) << "Test Statistic Sampler or Test Statistics not defined. Stop." << endl;
127 return 0;
128 }
129
130 // get a big list of all variables for convenient switching
131 RooArgSet *nullParams = fNullModel->GetPdf()->getParameters(*fData);
132 RooArgSet *altParams = fAltModel->GetPdf()->getParameters(*fData);
133 // save all parameters so we can set them back to what they were
134 RooArgSet *bothParams = fNullModel->GetPdf()->getParameters(*fData);
135 bothParams->add(*altParams,false);
136 RooArgSet *saveAll = (RooArgSet*) bothParams->snapshot();
137
138 // check whether we have a ToyMCSampler and if so, keep a pointer to it
139 ToyMCSampler* toymcs = dynamic_cast<ToyMCSampler*>( fTestStatSampler );
140
141
142 // evaluate test statistic on data
143 RooArgSet nullP(*nullSnapshot);
144 double obsTestStat;
145
146 RooArgList* allTS = NULL;
147 if( toymcs ) {
148 allTS = toymcs->EvaluateAllTestStatistics(*const_cast<RooAbsData*>(fData), nullP);
149 if (!allTS) return 0;
150 //oocoutP((TObject*)0,Generation) << "All Test Statistics on data: " << endl;
151 //allTS->Print("v");
152 RooRealVar* firstTS = (RooRealVar*)allTS->at(0);
153 obsTestStat = firstTS->getVal();
154 if (allTS->getSize()<=1) {
155 delete allTS;
156 allTS= 0; // don't save
157 }
158 }else{
159 obsTestStat = fTestStatSampler->EvaluateTestStatistic(*const_cast<RooAbsData*>(fData), nullP);
160 }
161 oocoutP((TObject*)0,Generation) << "Test Statistic on data: " << obsTestStat << endl;
162
163 // set parameters back ... in case the evaluation of the test statistic
164 // modified something (e.g. a nuisance parameter that is not randomized
165 // must be set here)
166 *bothParams = *saveAll;
167
168
169
170 // Generate sampling distribution for null
172 RooArgSet paramPointNull(*fNullModel->GetParametersOfInterest());
173 if(PreNullHook(&paramPointNull, obsTestStat) != 0) {
174 oocoutE((TObject*)0,Generation) << "PreNullHook did not return 0." << endl;
175 }
176 SamplingDistribution* samp_null = NULL;
177 RooDataSet* detOut_null = NULL;
178 if(toymcs) {
179 detOut_null = toymcs->GetSamplingDistributions(paramPointNull);
180 if( detOut_null ) {
181 samp_null = new SamplingDistribution( detOut_null->GetName(), detOut_null->GetTitle(), *detOut_null );
182 if (detOut_null->get()->getSize()<=1) {
183 delete detOut_null;
184 detOut_null= 0;
185 }
186 }
187 }else samp_null = fTestStatSampler->GetSamplingDistribution(paramPointNull);
188
189 // set parameters back
190 *bothParams = *saveAll;
191
192 // Generate sampling distribution for alternate
195 if(PreAltHook(&paramPointAlt, obsTestStat) != 0) {
196 oocoutE((TObject*)0,Generation) << "PreAltHook did not return 0." << endl;
197 }
198 SamplingDistribution* samp_alt = NULL;
199 RooDataSet* detOut_alt = NULL;
200 if(toymcs) {
201
202 // case of re-using same toys for every points
203 // set a given seed
204 unsigned int prevSeed = 0;
205 if (fAltToysSeed > 0) {
206 prevSeed = RooRandom::integer(std::numeric_limits<unsigned int>::max()-1)+1; // want to avoid zero value
208 }
209
210 detOut_alt = toymcs->GetSamplingDistributions(paramPointAlt);
211 if( detOut_alt ) {
212 samp_alt = new SamplingDistribution( detOut_alt->GetName(), detOut_alt->GetTitle(), *detOut_alt );
213 if (detOut_alt->get()->getSize()<=1) {
214 delete detOut_alt;
215 detOut_alt= 0;
216 }
217 }
218
219 // restore the seed
220 if (prevSeed > 0) {
222 }
223
224 }else samp_alt = fTestStatSampler->GetSamplingDistribution(paramPointAlt);
225
226
227 // create result
228 string resultname = "HypoTestCalculator_result";
229 HypoTestResult* res = new HypoTestResult(resultname.c_str());
231 res->SetTestStatisticData(obsTestStat);
232 res->SetAltDistribution(samp_alt);
233 res->SetNullDistribution(samp_null);
234 res->SetAltDetailedOutput( detOut_alt );
235 res->SetNullDetailedOutput( detOut_null );
236 res->SetAllTestStatisticsData( allTS );
237
238 const RooArgSet *aset = GetFitInfo();
239 if (aset != NULL) {
240 RooDataSet *dset = new RooDataSet(NULL, NULL, *aset);
241 dset->add(*aset);
242 res->SetFitInfo( dset );
243 }
244
245 *bothParams = *saveAll;
246 delete allTS;
247 delete bothParams;
248 delete saveAll;
249 delete altParams;
250 delete nullParams;
251 delete nullSnapshot;
252 PostHook();
253 return res;
254}
255
256////////////////////////////////////////////////////////////////////////////////
257/// to re-use same toys for alternate hypothesis
258
260 fAltToysSeed = RooRandom::integer(std::numeric_limits<unsigned int>::max()-1)+1;
261}
#define oocoutE(o, a)
Definition: RooMsgService.h:47
#define oocoutP(o, a)
Definition: RooMsgService.h:45
#define ClassImp(name)
Definition: Rtypes.h:365
RooArgSet * getParameters(const RooAbsData *data, Bool_t stripDisconnected=kTRUE) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
Definition: RooAbsArg.cxx:543
Int_t getSize() const
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:37
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:81
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:88
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:134
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
Definition: RooArgSet.h:88
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:31
virtual const RooArgSet * get(Int_t index) const
Return RooArgSet with coordinates of event 'index'.
virtual void add(const RooArgSet &row, Double_t weight=1.0, Double_t weightError=0)
Add a data point, with its coordinates specified in the 'data' argset, to the data set.
static UInt_t integer(UInt_t max, TRandom *generator=randomGenerator())
Return an integer uniformly distributed from [0,n-1].
Definition: RooRandom.cxx:102
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
Definition: RooRandom.cxx:54
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
Common base class for the Hypothesis Test Calculators.
virtual int PreNullHook(RooArgSet *, double) const
void SetupSampler(const ModelConfig &model) const
common setup for both models
virtual int PreAltHook(RooArgSet *, double) const
virtual HypoTestResult * GetHypoTest() const
inherited methods from HypoTestCalculator interface
virtual const RooArgSet * GetFitInfo() const
void UseSameAltToys()
Set this for re-using always the same toys for alternate hypothesis in case of calls at different nul...
HypoTestResult is a base class for results from hypothesis tests.
void SetAltDetailedOutput(RooDataSet *d)
void SetNullDetailedOutput(RooDataSet *d)
void SetAllTestStatisticsData(const RooArgList *tsd)
void SetPValueIsRightTail(Bool_t pr)
void SetTestStatisticData(const Double_t tsd)
void SetNullDistribution(SamplingDistribution *null)
void SetFitInfo(RooDataSet *d)
void SetAltDistribution(SamplingDistribution *alt)
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Definition: ModelConfig.h:30
const RooArgSet * GetParametersOfInterest() const
get RooArgSet containing the parameter of interest (return NULL if not existing)
Definition: ModelConfig.h:231
const RooArgSet * GetNuisanceParameters() const
get RooArgSet containing the nuisance parameters (return NULL if not existing)
Definition: ModelConfig.h:234
void LoadSnapshot() const
load the snapshot from ws if it exists
const RooArgSet * GetObservables() const
get RooArgSet for observables (return NULL if not existing)
Definition: ModelConfig.h:243
const RooArgSet * GetSnapshot() const
get RooArgSet for parameters for a particular hypothesis (return NULL if not existing)
RooAbsPdf * GetPdf() const
get model PDF (return NULL if pdf has not been specified or does not exist)
Definition: ModelConfig.h:228
TestStatistic that returns the ratio of profiled likelihoods.
This class simply holds a sampling distribution of some test statistic.
TestStatSampler is an interface class for a tools which produce RooStats SamplingDistributions.
virtual Double_t EvaluateTestStatistic(RooAbsData &data, RooArgSet &paramsOfInterest)=0
virtual void SetObservables(const RooArgSet &)=0
virtual TestStatistic * GetTestStatistic() const =0
virtual void SetParametersForTestStat(const RooArgSet &)=0
virtual void SetNuisanceParameters(const RooArgSet &)=0
virtual SamplingDistribution * GetSamplingDistribution(RooArgSet &paramsOfInterest)=0
virtual void SetSamplingDistName(const char *name)=0
virtual void SetPdf(RooAbsPdf &)=0
virtual bool PValueIsRightTail(void) const
Defines the sign convention of the test statistic. Overwrite function if necessary.
Definition: TestStatistic.h:45
ToyMCSampler is an implementation of the TestStatSampler interface.
Definition: ToyMCSampler.h:72
virtual RooArgList * EvaluateAllTestStatistics(RooAbsData &data, const RooArgSet &poi)
Evaluate all test statistics, returning result and any detailed output.
virtual RooDataSet * GetSamplingDistributions(RooArgSet &paramPoint)
Use for serial and parallel runs.
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
Mother of all ROOT objects.
Definition: TObject.h:37
virtual void SetSeed(ULong_t seed=0)
Set the random generator seed.
Definition: TRandom.cxx:597
@ Generation
Definition: RooGlobalFunc.h:57
@ InputArguments
Definition: RooGlobalFunc.h:58
Namespace for the RooStats classes.
Definition: Asimov.h:20