Logo ROOT  
Reference Guide
RooProfileLL.cxx
Go to the documentation of this file.
1  /*****************************************************************************
2  * Project: RooFit *
3  * *
4  * Copyright (c) 2000-2005, Regents of the University of California *
5  * and Stanford University. All rights reserved. *
6  * *
7  * Redistribution and use in source and binary forms, *
8  * with or without modification, are permitted according to the terms *
9  * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
10  *****************************************************************************/
11 
12 /**
13 \file RooProfileLL.cxx
14 \class RooProfileLL
15 \ingroup Roofitcore
16 
17 Class RooProfileLL implements the profile likelihood estimator for
18 a given likelihood and set of parameters of interest. The value return by
19 RooProfileLL is the input likelihood nll minimized w.r.t all nuisance parameters
20 (which are all parameters except for those listed in the constructor) minus
21 the -log(L) of the best fit. Note that this function is slow to evaluate
22 as a MIGRAD minimization step is executed for each function evaluation
23 **/
24 
25 #include "Riostream.h"
26 
27 #include "RooFit.h"
28 #include "RooProfileLL.h"
29 #include "RooAbsReal.h"
30 #include "RooMinuit.h"
31 #include "RooMinimizer.h"
32 #include "RooMsgService.h"
33 #include "RooRealVar.h"
34 
35 using namespace std ;
36 
38 
39 
40 ////////////////////////////////////////////////////////////////////////////////
41 /// Default constructor
42 /// Should only be used by proof.
43 
45  RooAbsReal("RooProfileLL","RooProfileLL"),
46  _nll(),
47  _obs("paramOfInterest","Parameters of interest",this),
48  _par("nuisanceParam","Nuisance parameters",this,kFALSE,kFALSE),
49  _startFromMin(kTRUE),
50  _minimizer(0),
51  _absMinValid(kFALSE),
52  _absMin(0),
53  _neval(0)
54 {
57 }
58 
59 
60 ////////////////////////////////////////////////////////////////////////////////
61 /// Constructor of profile likelihood given input likelihood nll w.r.t
62 /// the given set of variables. The input log likelihood is minimized w.r.t
63 /// to all other variables of the likelihood at each evaluation and the
64 /// value of the global log likelihood minimum is always subtracted.
65 
66 RooProfileLL::RooProfileLL(const char *name, const char *title,
67  RooAbsReal& nllIn, const RooArgSet& observables) :
68  RooAbsReal(name,title),
69  _nll("input","-log(L) function",this,nllIn),
70  _obs("paramOfInterest","Parameters of interest",this),
71  _par("nuisanceParam","Nuisance parameters",this,kFALSE,kFALSE),
72  _startFromMin(kTRUE),
73  _minimizer(0),
74  _absMinValid(kFALSE),
75  _absMin(0),
76  _neval(0)
77 {
78  // Determine actual parameters and observables
79  RooArgSet* actualObs = nllIn.getObservables(observables) ;
80  RooArgSet* actualPars = nllIn.getParameters(observables) ;
81 
82  _obs.add(*actualObs) ;
83  _par.add(*actualPars) ;
84 
85  delete actualObs ;
86  delete actualPars ;
87 
90 }
91 
92 
93 
94 ////////////////////////////////////////////////////////////////////////////////
95 /// Copy constructor
96 
97 RooProfileLL::RooProfileLL(const RooProfileLL& other, const char* name) :
98  RooAbsReal(other,name),
99  _nll("nll",this,other._nll),
100  _obs("obs",this,other._obs),
101  _par("par",this,other._par),
102  _startFromMin(other._startFromMin),
103  _minimizer(0),
104  _absMinValid(kFALSE),
105  _absMin(0),
106  _paramFixed(other._paramFixed),
107  _neval(0)
108 {
111 
114 
115 }
116 
117 
118 
119 ////////////////////////////////////////////////////////////////////////////////
120 /// Destructor
121 
123 {
124  // Delete instance of minuit if it was ever instantiated
125  if (_minimizer) {
126  delete _minimizer ;
127  }
128 
129  delete _piter ;
130  delete _oiter ;
131 }
132 
133 
134 
135 
136 ////////////////////////////////////////////////////////////////////////////////
137 
139 {
140  validateAbsMin() ;
141  return _paramAbsMin ;
142 }
143 
144 
145 ////////////////////////////////////////////////////////////////////////////////
146 
148 {
149  validateAbsMin() ;
150  return _obsAbsMin ;
151 }
152 
153 
154 
155 
156 ////////////////////////////////////////////////////////////////////////////////
157 /// Optimized implementation of createProfile for profile likelihoods.
158 /// Return profile of original function in terms of stated parameters
159 /// of interest rather than profiling recursively.
160 
162 {
163  return nll().createProfile(paramsOfInterest) ;
164 }
165 
166 
167 
168 
169 ////////////////////////////////////////////////////////////////////////////////
170 
172 {
173  coutI(Minimization) << "RooProfileLL::evaluate(" << GetName() << ") Creating instance of MINUIT" << endl ;
174 
177  _minimizer = new MINIMIZER(const_cast<RooAbsReal&>(_nll.arg())) ;
179 
180  //_minimizer->setPrintLevel(-999) ;
181  //_minimizer->setNoWarn() ;
182  //_minimizer->setVerbose(1) ;
183 }
184 
185 
186 
187 ////////////////////////////////////////////////////////////////////////////////
188 /// Evaluate profile likelihood by minimizing likelihood w.r.t. all
189 /// parameters that are not considered observables of this profile
190 /// likelihood object.
191 
193 {
194  // Instantiate minimizer if we haven't done that already
195  if (!_minimizer) {
197  }
198 
199  // Save current value of observables
200  RooArgSet* obsSetOrig = (RooArgSet*) _obs.snapshot() ;
201 
202  validateAbsMin() ;
203 
204 
205  // Set all observables constant in the minimization
206  const_cast<RooSetProxy&>(_obs).setAttribAll("Constant",kTRUE) ;
207  ccoutP(Eval) << "." ; ccoutP(Eval).flush() ;
208 
209  // If requested set initial parameters to those corresponding to absolute minimum
210  if (_startFromMin) {
211  const_cast<RooProfileLL&>(*this)._par = _paramAbsMin ;
212  }
213 
214  _minimizer->zeroEvalCount() ;
215 
216  //TString minim=::ROOT::Math::MinimizerOptions::DefaultMinimizerType().c_str();
217  //TString algorithm = ROOT::Math::MinimizerOptions::DefaultMinimizerAlgo().c_str();
218  //if (algorithm == "Migrad") algorithm = "Minimize"; // prefer to use Minimize instead of Migrad
219  //_minimizer->minimize(minim.Data(),algorithm.Data());
220 
221  _minimizer->migrad() ;
222  _neval = _minimizer->evalCounter() ;
223 
224  // Restore original values and constant status of observables
225  TIterator* iter = obsSetOrig->createIterator() ;
226  RooRealVar* var ;
227  while((var=dynamic_cast<RooRealVar*>(iter->Next()) ) ) {
228  RooRealVar* target = (RooRealVar*) _obs.find(var->GetName()) ;
229  target->setVal(var->getVal()) ;
230  target->setConstant(var->isConstant()) ;
231  }
232  delete iter ;
233  delete obsSetOrig ;
234 
235  return _nll - _absMin ;
236 }
237 
238 
239 
240 ////////////////////////////////////////////////////////////////////////////////
241 /// Check that parameters and likelihood value for 'best fit' are still valid. If not,
242 /// because the best fit has never been calculated, or because constant parameters have
243 /// changed value or parameters have changed const/float status, the minimum is recalculated
244 
246 {
247  // Check if constant status of any of the parameters have changed
248  if (_absMinValid) {
249  _piter->Reset() ;
250  RooAbsArg* par ;
251  while((par=(RooAbsArg*)_piter->Next())) {
252  if (_paramFixed[par->GetName()] != par->isConstant()) {
253  cxcoutI(Minimization) << "RooProfileLL::evaluate(" << GetName() << ") constant status of parameter " << par->GetName() << " has changed from "
254  << (_paramFixed[par->GetName()]?"fixed":"floating") << " to " << (par->isConstant()?"fixed":"floating")
255  << ", recalculating absolute minimum" << endl ;
256  _absMinValid = kFALSE ;
257  break ;
258  }
259  }
260  }
261 
262 
263  // If we don't have the absolute minimum w.r.t all observables, calculate that first
264  if (!_absMinValid) {
265 
266  cxcoutI(Minimization) << "RooProfileLL::evaluate(" << GetName() << ") determining minimum likelihood for current configurations w.r.t all observable" << endl ;
267 
268 
269  if (!_minimizer) {
271  }
272 
273  // Save current values of non-marginalized parameters
274  RooArgSet* obsStart = (RooArgSet*) _obs.snapshot(kFALSE) ;
275 
276  // Start from previous global minimum
277  if (_paramAbsMin.getSize()>0) {
278  const_cast<RooSetProxy&>(_par).assignValueOnly(_paramAbsMin) ;
279  }
280  if (_obsAbsMin.getSize()>0) {
281  const_cast<RooSetProxy&>(_obs).assignValueOnly(_obsAbsMin) ;
282  }
283 
284  // Find minimum with all observables floating
285  const_cast<RooSetProxy&>(_obs).setAttribAll("Constant",kFALSE) ;
286 
287  //TString minim=::ROOT::Math::MinimizerOptions::DefaultMinimizerType().c_str();
288  //TString algorithm = ROOT::Math::MinimizerOptions::DefaultMinimizerAlgo().c_str();
289  //if (algorithm == "Migrad") algorithm = "Minimize"; // prefer to use Minimize instead of Migrad
290  //_minimizer->minimize(minim.Data(),algorithm.Data());
291  _minimizer->migrad() ;
292 
293  // Save value and remember
294  _absMin = _nll ;
295  _absMinValid = kTRUE ;
296 
297  // Save parameter values at abs minimum as well
299 
300  // Only store non-constant parameters here!
301  RooArgSet* tmp = (RooArgSet*) _par.selectByAttrib("Constant",kFALSE) ;
302  _paramAbsMin.addClone(*tmp) ;
303  delete tmp ;
304 
306 
307  // Save constant status of all parameters
308  _piter->Reset() ;
309  RooAbsArg* par ;
310  while((par=(RooAbsArg*)_piter->Next())) {
311  _paramFixed[par->GetName()] = par->isConstant() ;
312  }
313 
314  if (dologI(Minimization)) {
315  cxcoutI(Minimization) << "RooProfileLL::evaluate(" << GetName() << ") minimum found at (" ;
316 
317  RooAbsReal* arg ;
318  Bool_t first=kTRUE ;
319  _oiter->Reset() ;
320  while ((arg=(RooAbsReal*)_oiter->Next())) {
321  ccxcoutI(Minimization) << (first?"":", ") << arg->GetName() << "=" << arg->getVal() ;
322  first=kFALSE ;
323  }
324  ccxcoutI(Minimization) << ")" << endl ;
325  }
326 
327  // Restore original parameter values
328  const_cast<RooSetProxy&>(_obs) = *obsStart ;
329  delete obsStart ;
330 
331  }
332 }
333 
334 
335 
336 ////////////////////////////////////////////////////////////////////////////////
337 
338 Bool_t RooProfileLL::redirectServersHook(const RooAbsCollection& /*newServerList*/, Bool_t /*mustReplaceAll*/,
339  Bool_t /*nameChange*/, Bool_t /*isRecursive*/)
340 {
341  if (_minimizer) {
342  delete _minimizer ;
343  _minimizer = 0 ;
344  }
345  return kFALSE ;
346 }
347 
348 
RooSetProxy::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Overloaded RooArgSet::add() method inserts 'var' into set and registers 'var' as server to owner with...
Definition: RooSetProxy.cxx:165
first
Definition: first.py:1
RooFit::Minimization
@ Minimization
Definition: RooGlobalFunc.h:67
RooRealVar::setVal
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Definition: RooRealVar.cxx:226
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
RooAbsReal.h
RooMinuit.h
RooMsgService.h
RooProfileLL::initializeMinimizer
void initializeMinimizer() const
Definition: RooProfileLL.cxx:171
RooFit.h
RooTemplateProxy::arg
const T & arg() const
Return reference to object held in proxy.
Definition: RooTemplateProxy.h:271
RooMsgService::setSilentMode
void setSilentMode(Bool_t flag)
Definition: RooMsgService.h:178
RooProfileLL::bestFitObs
const RooArgSet & bestFitObs() const
Definition: RooProfileLL.cxx:147
RooProfileLL::_oiter
TIterator * _oiter
Iterator over profile likelihood parameters to be minimized.
Definition: RooProfileLL.h:63
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
RooSetProxy
Definition: RooSetProxy.h:23
RooProfileLL
Definition: RooProfileLL.h:26
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooAbsCollection::find
RooAbsArg * find(const char *name) const
Find object with given name in list.
Definition: RooAbsCollection.cxx:813
coutI
#define coutI(a)
Definition: RooMsgService.h:30
RooAbsReal
Definition: RooAbsReal.h:61
RooProfileLL::_paramAbsMin
RooArgSet _paramAbsMin
Definition: RooProfileLL.h:69
bool
TIterator
Definition: TIterator.h:30
RooProfileLL::bestFitParams
const RooArgSet & bestFitParams() const
Definition: RooProfileLL.cxx:138
RooProfileLL::createProfile
virtual RooAbsReal * createProfile(const RooArgSet &paramsOfInterest)
Optimized implementation of createProfile for profile likelihoods.
Definition: RooProfileLL.cxx:161
RooAbsArg::getParameters
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:546
RooArgSet::snapshot
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:134
RooProfileLL::_piter
TIterator * _piter
Definition: RooProfileLL.h:62
RooMsgService::silentMode
Bool_t silentMode() const
Definition: RooMsgService.h:177
RooProfileLL::_absMinValid
Bool_t _absMinValid
Internal minuit instance.
Definition: RooProfileLL.h:67
RooProfileLL::_obs
RooSetProxy _obs
Definition: RooProfileLL.h:58
RooProfileLL::_startFromMin
Bool_t _startFromMin
Definition: RooProfileLL.h:60
RooProfileLL::redirectServersHook
virtual Bool_t redirectServersHook(const RooAbsCollection &, Bool_t, Bool_t, Bool_t)
Definition: RooProfileLL.cxx:338
RooProfileLL::_paramFixed
std::map< std::string, bool > _paramFixed
Definition: RooProfileLL.h:71
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
RooProfileLL::validateAbsMin
void validateAbsMin() const
Check that parameters and likelihood value for 'best fit' are still valid.
Definition: RooProfileLL.cxx:245
RooProfileLL::_minimizer
MINIMIZER * _minimizer
Iterator of profile likelihood output parameter(s)
Definition: RooProfileLL.h:65
dologI
#define dologI(a)
Definition: RooMsgService.h:66
RooAbsCollection::createIterator
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
Definition: RooAbsCollection.h:118
ccoutP
#define ccoutP(a)
Definition: RooMsgService.h:39
RooAbsCollection
Definition: RooAbsCollection.h:30
RooRealVar.h
TIterator::Next
virtual TObject * Next()=0
RooAbsReal::createProfile
virtual RooAbsReal * createProfile(const RooArgSet &paramsOfInterest)
Create a RooProfileLL object that eliminates all nuisance parameters in the present function.
Definition: RooAbsReal.cxx:516
TIterator::Reset
virtual void Reset()=0
RooProfileLL::_nll
RooRealProxy _nll
Definition: RooProfileLL.h:57
RooProfileLL::_absMin
Double_t _absMin
Definition: RooProfileLL.h:68
Double_t
double Double_t
Definition: RtypesCore.h:59
RooAbsArg::getObservables
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Return the observables of this pdf given a set of observables.
Definition: RooAbsArg.h:294
RooMinimizer.h
RooProfileLL::_obsAbsMin
RooArgSet _obsAbsMin
Definition: RooProfileLL.h:70
name
char name[80]
Definition: TGX11.cxx:110
RooAbsRealLValue::setConstant
void setConstant(Bool_t value=kTRUE)
Definition: RooAbsRealLValue.h:109
MINIMIZER
#define MINIMIZER
Definition: RooProfileLL.h:24
RooProfileLL::nll
RooAbsReal & nll()
Definition: RooProfileLL.h:39
RooProfileLL::~RooProfileLL
virtual ~RooProfileLL()
Destructor.
Definition: RooProfileLL.cxx:122
RooAbsArg
Definition: RooAbsArg.h:73
RooArgSet::addClone
virtual void addClone(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling addOwned() for each element in the source...
Definition: RooArgSet.h:96
RooProfileLL::evaluate
Double_t evaluate() const
Evaluate profile likelihood by minimizing likelihood w.r.t.
Definition: RooProfileLL.cxx:192
RooProfileLL::RooProfileLL
RooProfileLL()
Default constructor Should only be used by proof.
Definition: RooProfileLL.cxx:44
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:53
RooProfileLL::_par
RooSetProxy _par
Definition: RooProfileLL.h:59
RooRealVar
Definition: RooRealVar.h:35
RooProfileLL::_neval
Int_t _neval
Definition: RooProfileLL.h:72
RooAbsCollection::removeAll
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
Definition: RooAbsCollection.cxx:643
RooFit::Eval
@ Eval
Definition: RooGlobalFunc.h:68
Riostream.h
cxcoutI
#define cxcoutI(a)
Definition: RooMsgService.h:85
RooMsgService::instance
static RooMsgService & instance()
Return reference to singleton instance.
Definition: RooMsgService.cxx:363
RooAbsCollection::selectByAttrib
RooAbsCollection * selectByAttrib(const char *name, Bool_t value) const
Create a subset of the current collection, consisting only of those elements with the specified attri...
Definition: RooAbsCollection.cxx:677
RooAbsArg::isConstant
Bool_t isConstant() const
Check if the "Constant" attribute is set.
Definition: RooAbsArg.h:360
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:171
RooProfileLL.h
RooArgSet
Definition: RooArgSet.h:28
ccxcoutI
#define ccxcoutI(a)
Definition: RooMsgService.h:86