AsymptoticCalculator.h

Go to the documentation of this file.

// @(#)root/roostats:$Id$
// Author: Sven Kreiss, Kyle Cranmer   Nov 2010
/*************************************************************************
 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers.               *
 * All rights reserved.                                                  *
 *                                                                       *
 * For the licensing terms see $ROOTSYS/LICENSE.                         *
 * For the list of contributors see $ROOTSYS/README/CREDITS.             *
 *************************************************************************/
 
#ifndef ROOSTATS_AsymptoticCalculator
#define ROOSTATS_AsymptoticCalculator
 
#include "RooStats/HypoTestCalculatorGeneric.h"
#include "RooArgSet.h"
#include "Rtypes.h"
 
class RooArgList;
class RooCategory;
class RooRealVar;
class RooPoisson;
class RooProdPdf;
 
 
namespace RooStats {
 
   class AsymptoticCalculator : public HypoTestCalculatorGeneric {
 
   public:
      AsymptoticCalculator(
         RooAbsData &data,  // need to pass non-const since RooAbsPdf::fitTo takes a non-const data set
         const ModelConfig &altModel,
         const ModelConfig &nullModel,
         bool nominalAsimov = false
         );
      //    HypoTestCalculatorGeneric(data, altModel, nullModel, 0)
      // {
      // }
 
      ~AsymptoticCalculator() {
      }
 
      /// initialize the calculator by performing a global fit and make the Asimov data set
      bool Initialize() const;
 
      /// re-implement HypoTest computation using the asymptotic
      virtual HypoTestResult *GetHypoTest() const;
 
      /// make the asimov data from the ModelConfig and list of poi - return data set and snapshot of global obs
      /// poiValues is the snapshot of POI used for finding the best nuisance parameter values (conditioned at these values)
      /// genPoiValues is optionally a different set of POI values used for generating. By default the same POI are used for generating and for finding the nuisance parameters
      static RooAbsData * MakeAsimovData( RooAbsData & data, const ModelConfig & model,  const RooArgSet & poiValues, RooArgSet & globObs, const RooArgSet * genPoiValues = 0);
 
 
      /// make a nominal asimov data from the ModelConfig and parameter values
      /// The parameter values (including the nuisance) could be given from a fit to data or be at the nominal values
      static RooAbsData * MakeAsimovData( const ModelConfig & model,  const RooArgSet & allParamValues, RooArgSet & globObs);
 
 
 
      static RooAbsData * GenerateAsimovData(const RooAbsPdf & pdf, const RooArgSet & observables );
 
      /// function given the null and the alt p value - return the expected one given the N - sigma value
      static double GetExpectedPValues(double pnull, double palt, double nsigma, bool usecls, bool oneSided = true );
 
      /// set test statistic for one sided (upper limits)
      void SetOneSided(bool on) { fOneSided = on; }
 
      /// set the test statistics for two sided (in case of upper limits
      /// for discovery does not make really sense)
      void SetTwoSided() { fOneSided = false; fOneSidedDiscovery = false;}
 
      /// set the test statistics for one-sided discovery
      void SetOneSidedDiscovery(bool on) { fOneSidedDiscovery = on; }
 
      /// re-implementation of  setters since they are needed to re-initialize the calculator
      virtual void SetNullModel(const ModelConfig &nullModel) {
         HypoTestCalculatorGeneric::SetNullModel(nullModel);
         fIsInitialized = false;
      }
      virtual void SetAlternateModel(const ModelConfig &altModel) {
         HypoTestCalculatorGeneric::SetAlternateModel(altModel);
         fIsInitialized = false;
      }
      virtual void SetData(RooAbsData &data) {
         HypoTestCalculatorGeneric::SetData(data);
         fIsInitialized = false;
      }
 
 
      bool IsTwoSided() const { return (!fOneSided && !fOneSidedDiscovery); }
      bool IsOneSidedDiscovery() const { return fOneSidedDiscovery; }
 
 
      /// set using of qtilde, by default is controlled if RoORealVar is limited or not
      void SetQTilde(bool on) { fUseQTilde = on; }
 
      /// return snapshot of the best fit parameter
      const RooArgSet & GetBestFitPoi() const { return fBestFitPoi; }
      /// return best fit parameter (firs of poi)
      const RooRealVar * GetMuHat() const { return dynamic_cast<RooRealVar*>(fBestFitPoi.first()); }
      /// return best fit value for all parameters
      const RooArgSet & GetBestFitParams() const { return fBestFitPoi; }
 
      static void SetPrintLevel(int level);
 
   protected:
      // // configure TestStatSampler for the Null run
      // int PreNullHook(RooArgSet *parameterPoint, double obsTestStat) const;
 
      // // configure TestStatSampler for the Alt run
      // int PreAltHook(RooArgSet *parameterPoint, double obsTestStat) const;
 
 
      static RooAbsData * GenerateAsimovDataSinglePdf(const RooAbsPdf & pdf, const RooArgSet & obs,  const RooRealVar & weightVar,
                                                      RooCategory * channelCat = 0);
 
      static RooAbsData * GenerateCountingAsimovData(RooAbsPdf & pdf, const RooArgSet & obs,  const RooRealVar & weightVar,
                                                      RooCategory * channelCat = 0);
 
 
      static void FillBins(const RooAbsPdf & pdf, const RooArgList &obs, RooAbsData & data, int &index,  double
                           &binVolume, int &ibin);
 
      static double EvaluateNLL(RooAbsPdf & pdf, RooAbsData& data, const RooArgSet * condObs, const RooArgSet * globObs, const RooArgSet *poiSet = 0 );
 
      static bool SetObsToExpected(RooAbsPdf &pdf, const RooArgSet &obs);
      static bool SetObsToExpected(RooProdPdf &prod, const RooArgSet &obs);
 
   protected:
      ClassDef(AsymptoticCalculator,2)
 
   private:
 
      bool fOneSided;                     // for one sided PL test statistic (upper limits)
      mutable bool fOneSidedDiscovery;    // for one sided PL test statistic (for discovery)
      bool fNominalAsimov;                // make Asimov at nominal parameter values
      mutable bool fIsInitialized;        //! flag to check if calculator is initialized
      mutable int fUseQTilde;             // flag to indicate if using qtilde or not (-1 (default based on RooRealVar)), 0 false, 1 (true)
      static int fgPrintLevel;            // control print level  (0 minimal, 1 normal, 2 debug)
      mutable double fNLLObs;
      mutable double fNLLAsimov;
 
      mutable RooAbsData * fAsimovData;   // asimov data set
      mutable RooArgSet  fAsimovGlobObs;  // snapshot of Asimov global observables
      mutable RooArgSet  fBestFitPoi;     // snapshot of best fitted POI values
      mutable RooArgSet  fBestFitParams;  // snapshot of all best fitted Parameter values
 
 
   };
}
 
#endif