Logo ROOT   6.16/01
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
17Class RooProfileLL implements the profile likelihood estimator for
18a given likelihood and set of parameters of interest. The value return by
19RooProfileLL 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
21the -log(L) of the best fit. Note that this function is slow to evaluate
22as 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#include "RooMsgService.h"
35
36using namespace std ;
37
39
40
41////////////////////////////////////////////////////////////////////////////////
42/// Default constructor
43/// Should only be used by proof.
44
46 RooAbsReal("RooProfileLL","RooProfileLL"),
47 _nll(),
48 _obs("paramOfInterest","Parameters of interest",this),
49 _par("nuisanceParam","Nuisance parameters",this,kFALSE,kFALSE),
50 _startFromMin(kTRUE),
51 _minimizer(0),
52 _absMinValid(kFALSE),
53 _absMin(0),
54 _neval(0)
55{
58}
59
60
61////////////////////////////////////////////////////////////////////////////////
62/// Constructor of profile likelihood given input likelihood nll w.r.t
63/// the given set of variables. The input log likelihood is minimized w.r.t
64/// to all other variables of the likelihood at each evaluation and the
65/// value of the global log likelihood minimum is always subtracted.
66
67RooProfileLL::RooProfileLL(const char *name, const char *title,
68 RooAbsReal& nllIn, const RooArgSet& observables) :
69 RooAbsReal(name,title),
70 _nll("input","-log(L) function",this,nllIn),
71 _obs("paramOfInterest","Parameters of interest",this),
72 _par("nuisanceParam","Nuisance parameters",this,kFALSE,kFALSE),
73 _startFromMin(kTRUE),
74 _minimizer(0),
75 _absMinValid(kFALSE),
76 _absMin(0),
77 _neval(0)
78{
79 // Determine actual parameters and observables
80 RooArgSet* actualObs = nllIn.getObservables(observables) ;
81 RooArgSet* actualPars = nllIn.getParameters(observables) ;
82
83 _obs.add(*actualObs) ;
84 _par.add(*actualPars) ;
85
86 delete actualObs ;
87 delete actualPars ;
88
91}
92
93
94
95////////////////////////////////////////////////////////////////////////////////
96/// Copy constructor
97
98RooProfileLL::RooProfileLL(const RooProfileLL& other, const char* name) :
99 RooAbsReal(other,name),
100 _nll("nll",this,other._nll),
101 _obs("obs",this,other._obs),
102 _par("par",this,other._par),
103 _startFromMin(other._startFromMin),
104 _minimizer(0),
105 _absMinValid(kFALSE),
106 _absMin(0),
107 _paramFixed(other._paramFixed),
108 _neval(0)
109{
112
115
116}
117
118
119
120////////////////////////////////////////////////////////////////////////////////
121/// Destructor
122
124{
125 // Delete instance of minuit if it was ever instantiated
126 if (_minimizer) {
127 delete _minimizer ;
128 }
129
130 delete _piter ;
131 delete _oiter ;
132}
133
134
135
136
137////////////////////////////////////////////////////////////////////////////////
138
140{
142 return _paramAbsMin ;
143}
144
145
146////////////////////////////////////////////////////////////////////////////////
147
149{
151 return _obsAbsMin ;
152}
153
154
155
156
157////////////////////////////////////////////////////////////////////////////////
158/// Optimized implementation of createProfile for profile likelihoods.
159/// Return profile of original function in terms of stated parameters
160/// of interest rather than profiling recursively.
161
163{
164 return nll().createProfile(paramsOfInterest) ;
165}
166
167
168
169
170////////////////////////////////////////////////////////////////////////////////
171
173{
174 coutI(Minimization) << "RooProfileLL::evaluate(" << GetName() << ") Creating instance of MINUIT" << endl ;
175
178 _minimizer = new MINIMIZER(const_cast<RooAbsReal&>(_nll.arg())) ;
180
181 //_minimizer->setPrintLevel(-999) ;
182 //_minimizer->setNoWarn() ;
183 //_minimizer->setVerbose(1) ;
184}
185
186
187
188////////////////////////////////////////////////////////////////////////////////
189/// Evaluate profile likelihood by minimizing likelihood w.r.t. all
190/// parameters that are not considered observables of this profile
191/// likelihood object.
192
194{
195 // Instantiate minimizer if we haven't done that already
196 if (!_minimizer) {
198 }
199
200 // Save current value of observables
201 RooArgSet* obsSetOrig = (RooArgSet*) _obs.snapshot() ;
202
204
205
206 // Set all observables constant in the minimization
207 const_cast<RooSetProxy&>(_obs).setAttribAll("Constant",kTRUE) ;
208 ccoutP(Eval) << "." ; ccoutP(Eval).flush() ;
209
210 // If requested set initial parameters to those corresponding to absolute minimum
211 if (_startFromMin) {
212 const_cast<RooProfileLL&>(*this)._par = _paramAbsMin ;
213 }
214
215 _minimizer->zeroEvalCount() ;
216
217 //TString minim=::ROOT::Math::MinimizerOptions::DefaultMinimizerType().c_str();
218 //TString algorithm = ROOT::Math::MinimizerOptions::DefaultMinimizerAlgo().c_str();
219 //if (algorithm == "Migrad") algorithm = "Minimize"; // prefer to use Minimize instead of Migrad
220 //_minimizer->minimize(minim.Data(),algorithm.Data());
221
222 _minimizer->migrad() ;
223 _neval = _minimizer->evalCounter() ;
224
225 // Restore original values and constant status of observables
226 TIterator* iter = obsSetOrig->createIterator() ;
227 RooRealVar* var ;
228 while((var=dynamic_cast<RooRealVar*>(iter->Next()) ) ) {
229 RooRealVar* target = (RooRealVar*) _obs.find(var->GetName()) ;
230 target->setVal(var->getVal()) ;
231 target->setConstant(var->isConstant()) ;
232 }
233 delete iter ;
234 delete obsSetOrig ;
235
236 return _nll - _absMin ;
237}
238
239
240
241////////////////////////////////////////////////////////////////////////////////
242/// Check that parameters and likelihood value for 'best fit' are still valid. If not,
243/// because the best fit has never been calculated, or because constant parameters have
244/// changed value or parameters have changed const/float status, the minimum is recalculated
245
247{
248 // Check if constant status of any of the parameters have changed
249 if (_absMinValid) {
250 _piter->Reset() ;
251 RooAbsArg* par ;
252 while((par=(RooAbsArg*)_piter->Next())) {
253 if (_paramFixed[par->GetName()] != par->isConstant()) {
254 cxcoutI(Minimization) << "RooProfileLL::evaluate(" << GetName() << ") constant status of parameter " << par->GetName() << " has changed from "
255 << (_paramFixed[par->GetName()]?"fixed":"floating") << " to " << (par->isConstant()?"fixed":"floating")
256 << ", recalculating absolute minimum" << endl ;
258 break ;
259 }
260 }
261 }
262
263
264 // If we don't have the absolute minimum w.r.t all observables, calculate that first
265 if (!_absMinValid) {
266
267 cxcoutI(Minimization) << "RooProfileLL::evaluate(" << GetName() << ") determining minimum likelihood for current configurations w.r.t all observable" << endl ;
268
269
270 if (!_minimizer) {
272 }
273
274 // Save current values of non-marginalized parameters
275 RooArgSet* obsStart = (RooArgSet*) _obs.snapshot(kFALSE) ;
276
277 // Start from previous global minimum
278 if (_paramAbsMin.getSize()>0) {
279 const_cast<RooSetProxy&>(_par).assignValueOnly(_paramAbsMin) ;
280 }
281 if (_obsAbsMin.getSize()>0) {
282 const_cast<RooSetProxy&>(_obs).assignValueOnly(_obsAbsMin) ;
283 }
284
285 // Find minimum with all observables floating
286 const_cast<RooSetProxy&>(_obs).setAttribAll("Constant",kFALSE) ;
287
288 //TString minim=::ROOT::Math::MinimizerOptions::DefaultMinimizerType().c_str();
289 //TString algorithm = ROOT::Math::MinimizerOptions::DefaultMinimizerAlgo().c_str();
290 //if (algorithm == "Migrad") algorithm = "Minimize"; // prefer to use Minimize instead of Migrad
291 //_minimizer->minimize(minim.Data(),algorithm.Data());
292 _minimizer->migrad() ;
293
294 // Save value and remember
295 _absMin = _nll ;
297
298 // Save parameter values at abs minimum as well
300
301 // Only store non-constant parameters here!
302 RooArgSet* tmp = (RooArgSet*) _par.selectByAttrib("Constant",kFALSE) ;
303 _paramAbsMin.addClone(*tmp) ;
304 delete tmp ;
305
307
308 // Save constant status of all parameters
309 _piter->Reset() ;
310 RooAbsArg* par ;
311 while((par=(RooAbsArg*)_piter->Next())) {
312 _paramFixed[par->GetName()] = par->isConstant() ;
313 }
314
315 if (dologI(Minimization)) {
316 cxcoutI(Minimization) << "RooProfileLL::evaluate(" << GetName() << ") minimum found at (" ;
317
318 RooAbsReal* arg ;
320 _oiter->Reset() ;
321 while ((arg=(RooAbsReal*)_oiter->Next())) {
322 ccxcoutI(Minimization) << (first?"":", ") << arg->GetName() << "=" << arg->getVal() ;
323 first=kFALSE ;
324 }
325 ccxcoutI(Minimization) << ")" << endl ;
326 }
327
328 // Restore original parameter values
329 const_cast<RooSetProxy&>(_obs) = *obsStart ;
330 delete obsStart ;
331
332 }
333}
334
335
336
337////////////////////////////////////////////////////////////////////////////////
338
339Bool_t RooProfileLL::redirectServersHook(const RooAbsCollection& /*newServerList*/, Bool_t /*mustReplaceAll*/,
340 Bool_t /*nameChange*/, Bool_t /*isRecursive*/)
341{
342 if (_minimizer) {
343 delete _minimizer ;
344 _minimizer = 0 ;
345 }
346 return kFALSE ;
347}
348
349
#define coutI(a)
Definition: RooMsgService.h:31
#define cxcoutI(a)
Definition: RooMsgService.h:83
#define ccoutP(a)
Definition: RooMsgService.h:39
#define dologI(a)
Definition: RooMsgService.h:64
#define ccxcoutI(a)
Definition: RooMsgService.h:84
#define MINIMIZER
Definition: RooProfileLL.h:24
const Bool_t kFALSE
Definition: RtypesCore.h:88
bool Bool_t
Definition: RtypesCore.h:59
double Double_t
Definition: RtypesCore.h:55
const Bool_t kTRUE
Definition: RtypesCore.h:87
#define ClassImp(name)
Definition: Rtypes.h:363
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Definition: RooAbsArg.h:194
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:533
Bool_t isConstant() const
Definition: RooAbsArg.h:266
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
Int_t getSize() const
RooAbsCollection * snapshot(Bool_t deepCopy=kTRUE) const
Take a snap shot of current collection contents: An owning collection is returned containing clones o...
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...
TIterator * createIterator(Bool_t dir=kIterForward) const
RooAbsArg * find(const char *name) const
Find object with given name in list.
void setConstant(Bool_t value=kTRUE)
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
virtual RooAbsReal * createProfile(const RooArgSet &paramsOfInterest)
Create a RooProfileLL object that eliminates all nuisance parameters in the present function.
Definition: RooAbsReal.cxx:462
Double_t getVal(const RooArgSet *set=0) const
Evaluate object. Returns either cached value or triggers a recalculation.
Definition: RooAbsReal.h:64
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
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
static RooMsgService & instance()
Return reference to singleton instance.
Bool_t silentMode() const
void setSilentMode(Bool_t flag)
Class RooProfileLL implements the profile likelihood estimator for a given likelihood and set of para...
Definition: RooProfileLL.h:26
virtual ~RooProfileLL()
Destructor.
const RooArgSet & bestFitObs() const
RooArgSet _paramAbsMin
Definition: RooProfileLL.h:69
RooSetProxy _obs
Definition: RooProfileLL.h:58
RooProfileLL()
Default constructor Should only be used by proof.
void initializeMinimizer() const
TIterator * _piter
Definition: RooProfileLL.h:62
MINIMIZER * _minimizer
Iterator of profile likelihood output parameter(s)
Definition: RooProfileLL.h:65
RooAbsReal & nll()
Definition: RooProfileLL.h:39
std::map< std::string, bool > _paramFixed
Definition: RooProfileLL.h:71
Double_t evaluate() const
Evaluate profile likelihood by minimizing likelihood w.r.t.
RooArgSet _obsAbsMin
Definition: RooProfileLL.h:70
Bool_t _startFromMin
Definition: RooProfileLL.h:60
TIterator * _oiter
Iterator over profile likelihood parameters to be minimized.
Definition: RooProfileLL.h:63
void validateAbsMin() const
Check that parameters and likelihood value for 'best fit' are still valid.
RooSetProxy _par
Definition: RooProfileLL.h:59
virtual RooAbsReal * createProfile(const RooArgSet &paramsOfInterest)
Optimized implementation of createProfile for profile likelihoods.
virtual Bool_t redirectServersHook(const RooAbsCollection &, Bool_t, Bool_t, Bool_t)
RooRealProxy _nll
Definition: RooProfileLL.h:57
Double_t _absMin
Definition: RooProfileLL.h:68
const RooArgSet & bestFitParams() const
Bool_t _absMinValid
Internal minuit instance.
Definition: RooProfileLL.h:67
const RooAbsReal & arg() const
Definition: RooRealProxy.h:43
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Definition: RooRealVar.cxx:204
RooSetProxy is the concrete proxy for RooArgSet objects.
Definition: RooSetProxy.h:24
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...
Iterator abstract base class.
Definition: TIterator.h:30
virtual void Reset()=0
virtual TObject * Next()=0
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
@ Minimization
Definition: RooGlobalFunc.h:57
Definition: first.py:1
STL namespace.