Logo ROOT   master
Reference Guide
ProfileInspector.cxx
Go to the documentation of this file.
1 // @(#)root/roostats:$Id: ProfileInspector.cxx 34109 2010-06-24 15:00:16Z moneta $
2 
3 /*************************************************************************
4  * Project: RooStats *
5  * Package: RooFit/RooStats *
6  * Authors: *
7  * Kyle Cranmer, Lorenzo Moneta, Gregory Schott, Wouter Verkerke *
8  * Akira Shibata
9  *************************************************************************
10  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
11  * All rights reserved. *
12  * *
13  * For the licensing terms see $ROOTSYS/LICENSE. *
14  * For the list of contributors see $ROOTSYS/README/CREDITS. *
15  *************************************************************************/
16 
17 
18 /** \class RooStats::ProfileInspector
19  \ingroup Roostats
20 
21 Utility class to plot conditional MLE of nuisance parameters vs. Parameters of Interest
22 
23 */
24 
25 
27 #include "RooRealVar.h"
28 #include "RooAbsReal.h"
29 #include "RooArgSet.h"
30 #include "RooAbsPdf.h"
31 #include "RooCurve.h"
32 #include "TAxis.h"
33 
35 
36 using namespace RooStats;
37 using namespace std;
38 
39 ////////////////////////////////////////////////////////////////////////////////
40 
42 {
43 }
44 
45 ////////////////////////////////////////////////////////////////////////////////
46 /// ProfileInspector destructor
47 
49 {
50 }
51 
52 ////////////////////////////////////////////////////////////////////////////////
53 ///
54 /// This tool makes a plot of the conditional maximum likelihood estimate of the nuisance parameter
55 /// vs the parameter of interest
56 ///
57 /// This enables you to discover if any of the nuisance parameters are behaving strangely
58 /// curve is the optional parameters, when used you can specify the points previously scanned
59 /// in the process of plotOn or createHistogram.
60 /// To do this, you can do the following after the plot has been made:
61 /// ~~~ {.cpp}
62 /// profile, RooRealVar * poi, RooCurve * curve ){
63 ///RooCurve * curve = 0;
64 /// ~~~
65 
67 {
68 
69  const RooArgSet* poi_set = config->GetParametersOfInterest();
70  const RooArgSet* nuis_params=config->GetNuisanceParameters();
71  RooAbsPdf* pdf = config->GetPdf();
72 
73 
74  if(!poi_set){
75  cout << "no parameters of interest" << endl;
76  return 0;
77  }
78 
79  if(poi_set->getSize()!=1){
80  cout << "only one parameter of interest is supported currently" << endl;
81  return 0;
82  }
83  RooRealVar* poi = (RooRealVar*) poi_set->first();
84 
85 
86  if(!nuis_params){
87  cout << "no nuisance parameters" << endl;
88  return 0;
89  }
90 
91  if(!pdf){
92  cout << "pdf not set" << endl;
93  return 0;
94  }
95 
96  RooAbsReal* nll = pdf->createNLL(data);
97  RooAbsReal* profile = nll->createProfile(*poi);
98 
99  TList * list = new TList;
100  Int_t curve_N=100;
101  Double_t* curve_x=0;
102 // if(curve){
103 // curve_N=curve->GetN();
104 // curve_x=curve->GetX();
105 // } else {
106  Double_t max = dynamic_cast<RooAbsRealLValue*>(poi)->getMax();
107  Double_t min = dynamic_cast<RooAbsRealLValue*>(poi)->getMin();
108  Double_t step = (max-min)/(curve_N-1);
109  curve_x=new Double_t[curve_N];
110  for(int i=0; i<curve_N; ++i){
111  curve_x[i]=min+step*i;
112  }
113 // }
114 
115  map<string, std::vector<Double_t> > name_val;
116  for(int i=0; i<curve_N; i++){
117  poi->setVal(curve_x[i]);
118  profile->getVal();
119 
120  TIterator* nuis_params_itr=nuis_params->createIterator();
121  TObject* nuis_params_obj;
122  while((nuis_params_obj=nuis_params_itr->Next())){
123  RooRealVar* nuis_param = dynamic_cast<RooRealVar*>(nuis_params_obj);
124  if(nuis_param) {
125  string name = nuis_param->GetName();
126  if(nuis_params->getSize()==0) continue;
127  if(nuis_param && (! nuis_param->isConstant())){
128  if(name_val.find(name)==name_val.end()) name_val[name]=std::vector<Double_t>(curve_N);
129  name_val[name][i]=nuis_param->getVal();
130 
131  if(i==curve_N-1){
132  TGraph* g = new TGraph(curve_N, curve_x, &(name_val[name].front()));
133  g->SetName((name+"_"+string(poi->GetName())+"_profile").c_str());
134  g->GetXaxis()->SetTitle(poi->GetName());
135  g->GetYaxis()->SetTitle(nuis_param->GetName());
136  g->SetTitle("");
137  list->Add(g);
138  }
139  }
140  }
141  }
142  }
143 
144  delete [] curve_x;
145 
146 
147  delete nll;
148  delete profile;
149  return list;
150 }
virtual RooAbsReal * createNLL(RooAbsData &data, const RooLinkedList &cmdList)
Construct representation of -log(L) of PDFwith given dataset.
Definition: RooAbsPdf.cxx:916
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Definition: ModelConfig.h:30
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
#define g(i)
Definition: RSha256.hxx:105
STL namespace.
Iterator abstract base class.
Definition: TIterator.h:30
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:35
A doubly linked list.
Definition: TList.h:44
virtual void setVal(Double_t value)
Set value of variable to &#39;value&#39;.
Definition: RooRealVar.cxx:225
Int_t getSize() const
RooAbsArg * first() const
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:90
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:44
Namespace for the RooStats classes.
Definition: Asimov.h:19
RooAbsPdf * GetPdf() const
get model PDF (return NULL if pdf has not been specified or does not exist)
Definition: ModelConfig.h:234
const RooArgSet * GetParametersOfInterest() const
get RooArgSet containing the parameter of interest (return NULL if not existing)
Definition: ModelConfig.h:237
#define ClassImp(name)
Definition: Rtypes.h:361
double Double_t
Definition: RtypesCore.h:57
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:60
Utility class to plot conditional MLE of nuisance parameters vs.
virtual ~ProfileInspector()
ProfileInspector destructor.
Mother of all ROOT objects.
Definition: TObject.h:37
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
RooAbsPdf, the base class of all PDFs
Definition: RooAbsPdf.h:40
virtual void Add(TObject *obj)
Definition: TList.h:87
const RooArgSet * GetNuisanceParameters() const
get RooArgSet containing the nuisance parameters (return NULL if not existing)
Definition: ModelConfig.h:240
TList * GetListOfProfilePlots(RooAbsData &data, RooStats::ModelConfig *config)
This tool makes a plot of the conditional maximum likelihood estimate of the nuisance parameter vs th...
virtual RooAbsReal * createProfile(const RooArgSet &paramsOfInterest)
Create a RooProfileLL object that eliminates all nuisance parameters in the present function...
Definition: RooAbsReal.cxx:515
virtual TObject * Next()=0
A TGraph is an object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
char name[80]
Definition: TGX11.cxx:109
Bool_t isConstant() const
Check if the "Constant" attribute is set.
Definition: RooAbsArg.h:360