root/RooStats/HypoTestInverter.h

0001 // @(#)root/roostats:$Id$
0002 // Author: Kyle Cranmer, Lorenzo Moneta, Gregory Schott, Wouter Verkerke
0003 /*************************************************************************
0004  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers.               *
0005  * All rights reserved.                                                  *
0006  *                                                                       *
0007  * For the licensing terms see $ROOTSYS/LICENSE.                         *
0008  * For the list of contributors see $ROOTSYS/README/CREDITS.             *
0009  *************************************************************************/
0010
0011 #ifndef ROOSTATS_HypoTestInverter
0012 #define ROOSTATS_HypoTestInverter
0013
0014
0015 #include "RooStats/IntervalCalculator.h"
0016
0017
0018 #include "RooStats/HypoTestInverterResult.h"
0019
0020 class RooRealVar;
0021 class TGraphErrors;
0022
0023 #include <memory>
0024 #include <string>
0025
0026 namespace RooStats {
0027
0028    //class HypoTestCalculator;
0029    class HybridCalculator;
0030    class FrequentistCalculator;
0031    class AsymptoticCalculator;
0032    class HypoTestCalculatorGeneric;
0033    class TestStatistic;
0034
0035
0036 class HypoTestInverter : public IntervalCalculator {
0037
0038 public:
0039
0040    enum ECalculatorType { kUndefined = 0, kHybrid = 1, kFrequentist = 2, kAsymptotic = 3};
0041
0042    /// default constructor (used only for I/O)
0043    HypoTestInverter();
0044
0045    /// constructor from generic hypotest calculator
0046    HypoTestInverter( HypoTestCalculatorGeneric & hc,
0047                      RooRealVar* scannedVariable =nullptr,
0048                      double size = 0.05) ;
0049
0050
0051    /// constructor from hybrid calculator
0052    HypoTestInverter( HybridCalculator & hc,
0053                      RooRealVar* scannedVariable = nullptr,
0054                      double size = 0.05) ;
0055
0056    /// constructor from frequentist calculator
0057    HypoTestInverter( FrequentistCalculator & hc,
0058                      RooRealVar* scannedVariable,
0059                      double size = 0.05) ;
0060
0061    /// constructor from asymptotic calculator
0062    HypoTestInverter( AsymptoticCalculator & hc,
0063                      RooRealVar* scannedVariable,
0064                      double size = 0.05) ;
0065
0066    /// constructor from two ModelConfigs (first sb (the null model) then b (the alt model)
0067    /// creating a calculator inside
0068    HypoTestInverter( RooAbsData& data, ModelConfig &sb, ModelConfig &b,
0069            RooRealVar * scannedVariable = nullptr,  ECalculatorType type = kFrequentist,
0070            double size = 0.05) ;
0071
0072
0073    HypoTestInverterResult* GetInterval() const override;
0074
0075    void Clear();
0076
0077    /// Set up to perform a fixed scan.
0078    /// \param[in] nBins Number of points to scan.
0079    /// \param[in] xMin Lower limit of range to be scanned.
0080    /// \param[in] xMax Upper limit of range to be scanned.
0081    /// \param[in] scanLog Run in logarithmic steps along x.
0082    void SetFixedScan(int nBins, double xMin = 1, double xMax = -1, bool scanLog = false ) {
0083       fNBins = nBins;
0084       fXmin = xMin; fXmax = xMax;
0085       fScanLog = scanLog;
0086    }
0087
0088    /// Use automatic scanning, i.e. adaptive
0089    void SetAutoScan() { SetFixedScan(0); }
0090
0091    /// Run a fixed scan.
0092    /// \param[in] nBins Number of points to scan.
0093    /// \param[in] xMin Lower limit of range to be scanned.
0094    /// \param[in] xMax Upper limit of range to be scanned.
0095    /// \param[in] scanLog Run in logarithmic steps along x.
0096    bool RunFixedScan( int nBins, double xMin, double xMax, bool scanLog = false ) const;
0097
0098    bool RunOnePoint( double thisX, bool adaptive = false, double clTarget = -1 ) const;
0099
0100    //bool RunAutoScan( double xMin, double xMax, double target, double epsilon=nullptr.005, unsigned int numAlgorithm=nullptr );
0101
0102    bool RunLimit(double &limit, double &limitErr, double absTol = 0, double relTol = 0, const double *hint=nullptr) const;
0103
0104    void UseCLs( bool on = true) { fUseCLs = on; if (fResults) fResults->UseCLs(on);   }
0105
0106    void  SetData(RooAbsData &) override;
0107
0108    void SetModel(const ModelConfig &) override { } // not needed
0109
0110    /// set the size of the test (rate of Type I error) ( Eg. 0.05 for a 95% Confidence Interval)
0111    void SetTestSize(double size) override {fSize = size; if (fResults) fResults->SetTestSize(size); }
0112    /// set the confidence level for the interval (eg. 0.95 for a 95% Confidence Interval)
0113    void SetConfidenceLevel(double cl) override {fSize = 1.-cl;  if (fResults) fResults->SetConfidenceLevel(cl); }
0114    /// Get the size of the test (eg. rate of Type I error)
0115    double Size() const override {return fSize;}
0116    /// Get the Confidence level for the test
0117    double ConfidenceLevel()  const override {return 1.-fSize;}
0118
0119    /// destructor
0120    ~HypoTestInverter() override ;
0121
0122    /// retrieved a reference to the internally used HypoTestCalculator
0123    /// it might be invalid when the class is deleted
0124    HypoTestCalculatorGeneric * GetHypoTestCalculator() const { return fCalculator0; }
0125
0126    /// get the upper/lower limit distribution
0127    SamplingDistribution * GetLowerLimitDistribution(bool rebuild=false, int nToys = 100);
0128    SamplingDistribution * GetUpperLimitDistribution(bool rebuild=false, int nToys = 100);
0129
0130    /// function to rebuild the distributions
0131    SamplingDistribution * RebuildDistributions(bool isUpper=true, int nToys = 100,
0132                                                TList *clsDist = nullptr, TList *clsbDist = nullptr, TList *clbDist = nullptr, const char *outputfile = "HypoTestInverterRebuiltDist.root");
0133
0134    /// get the test statistic
0135    TestStatistic * GetTestStatistic() const;
0136
0137    /// set the test statistic
0138    bool SetTestStatistic(TestStatistic& stat);
0139
0140    /// set verbose level (0,1,2)
0141    void SetVerbose(int level=1) { fVerbose = level; }
0142
0143    /// set maximum number of toys
0144    void SetMaximumToys(int ntoys) { fMaxToys = ntoys;}
0145
0146    /// set numerical error in test statistic evaluation (default is zero)
0147    void SetNumErr(double err) { fNumErr = err; }
0148
0149 protected:
0150
0151    /// copy c-tor
0152    HypoTestInverter(const HypoTestInverter & rhs);
0153
0154    /// assignment
0155    HypoTestInverter & operator=(const HypoTestInverter & rhs);
0156
0157    void CreateResults() const;
0158
0159    /// run the hybrid at a single point
0160    HypoTestResult * Eval( HypoTestCalculatorGeneric &hc, bool adaptive , double clsTarget) const;
0161
0162    /// helper functions
0163    static RooRealVar * GetVariableToScan(const HypoTestCalculatorGeneric &hc);
0164    static void CheckInputModels(const HypoTestCalculatorGeneric &hc, const RooRealVar & scanVar);
0165
0166 private:
0167
0168
0169    static unsigned int fgNToys;
0170    static double fgCLAccuracy;
0171    static double fgAbsAccuracy;
0172    static double fgRelAccuracy;
0173    static std::string fgAlgo;
0174
0175    // graph, used to compute the limit, not just for plotting!
0176    mutable std::unique_ptr<TGraphErrors> fLimitPlot;  ///<! plot of limits
0177
0178
0179    // performance counter: remember how many toys have been thrown
0180    mutable int fTotalToysRun;
0181    int fMaxToys;  ///< maximum number of toys to run
0182
0183    HypoTestCalculatorGeneric* fCalculator0;        ///< pointer to the calculator passed in the constructor
0184    std::unique_ptr<HypoTestCalculatorGeneric> fHC; ///<! pointer to the generic hypotest calculator used
0185    RooRealVar* fScannedVariable;                   ///< pointer to the constrained variable
0186    mutable HypoTestInverterResult* fResults;       ///< pointer to the result
0187
0188    bool fUseCLs;
0189    bool fScanLog;
0190    double fSize;
0191    int fVerbose;
0192    ECalculatorType fCalcType;
0193    int fNBins;
0194    double fXmin;
0195    double fXmax;
0196    double fNumErr;
0197
0198 protected:
0199
0200    ClassDefOverride(HypoTestInverter,4)  // HypoTestInverter class
0201
0202 };
0203
0204 }
0205
0206 #endif