include/root/TRolke.h

0001
0002 //////////////////////////////////////////////////////////////////////////////
0003 //
0004 //  TRolke
0005 //
0006 //  This class computes confidence intervals for the rate of a Poisson
0007 //  in the presence of background and efficiency with a fully frequentist
0008 //  treatment of the uncertainties in the efficiency and background estimate
0009 //  using the profile likelihood method.
0010 //
0011 //      Author: Jan Conrad (CERN) 2004
0012 //      Updated: Johan Lundberg (CERN) 2009
0013 //
0014 //      Copyright CERN 2004,2009           Jan.Conrad@cern.ch,
0015 //                                     Johan.Lundberg@cern.ch
0016 //
0017 //  For information about the statistical meaning of the parameters
0018 //  and the syntax, consult TRolke.cxx
0019 //                  ------------------
0020 //
0021 //  Examples are found in the file Rolke.C
0022 //  --------------------------------------
0023 //
0024 //////////////////////////////////////////////////////////////////////////////
0025
0026 #ifndef ROOT_TRolke
0027 #define ROOT_TRolke
0028
0029 #include "TObject.h"
0030 #include "TMath.h"
0031
0032 // Class definition. This class is not intended to be used as a base class.
0033 class TRolke : public TObject
0034 {
0035
0036 private:
0037    Double_t fCL;         // confidence level as a fraction [0.9 for 90% ]
0038    Double_t fUpperLimit; // the calculated upper limit
0039    Double_t fLowerLimit; // the calculated lower limit
0040    bool fBounding;       // false for unbounded likelihood
0041                          // true for bounded likelihood
0042    Int_t fNumWarningsDeprecated1;
0043    Int_t fNumWarningsDeprecated2;
0044
0045    /* ----------------------------------------------------------------- */
0046    /* These variables are set by the Set methods for the various models */
0047    Int_t f_x;
0048    Int_t f_y;
0049    Int_t f_z;
0050    Double_t f_bm;
0051    Double_t f_em;
0052    Double_t f_e;
0053    Int_t f_mid;
0054    Double_t f_sde;
0055    Double_t f_sdb;
0056    Double_t f_tau;
0057    Double_t f_b;
0058    Int_t f_m;
0059
0060    /* ----------------------------------------------------------------- */
0061    /* Internal helper functions and methods */
0062    // The Calculator
0063    Double_t Interval(Int_t x, Int_t y, Int_t z, Double_t bm, Double_t em, Double_t e, Int_t mid, Double_t sde, Double_t sdb, Double_t tau, Double_t b, Int_t m);
0064
0065    // LIKELIHOOD ROUTINE
0066    Double_t Likelihood(Double_t mu, Int_t x, Int_t y, Int_t z, Double_t bm, Double_t em, Int_t mid, Double_t sde, Double_t sdb, Double_t tau, Double_t b, Int_t m, Int_t what);
0067
0068    //MODEL 1
0069    Double_t EvalLikeMod1(Double_t mu, Int_t x, Int_t y, Int_t z, Double_t tau, Int_t m, Int_t what);
0070    Double_t LikeMod1(Double_t mu, Double_t b, Double_t e, Int_t x, Int_t y, Int_t z, Double_t tau, Int_t m);
0071    void     ProfLikeMod1(Double_t mu, Double_t &b, Double_t &e, Int_t x, Int_t y, Int_t z, Double_t tau, Int_t m);
0072    Double_t LikeGradMod1(Double_t e, Double_t mu, Int_t x, Int_t y, Int_t z, Double_t tau, Int_t m);
0073
0074    //MODEL 2
0075    Double_t EvalLikeMod2(Double_t mu, Int_t x, Int_t y, Double_t em, Double_t sde, Double_t tau, Int_t what);
0076
0077    Double_t LikeMod2(Double_t mu, Double_t b, Double_t e, Int_t x, Int_t y, Double_t em, Double_t tau, Double_t v);
0078
0079    //MODEL 3
0080    Double_t EvalLikeMod3(Double_t mu, Int_t x, Double_t bm, Double_t em, Double_t sde, Double_t sdb, Int_t what);
0081    Double_t LikeMod3(Double_t mu, Double_t b, Double_t e, Int_t x, Double_t bm, Double_t em, Double_t u, Double_t v);
0082
0083    //MODEL 4
0084    Double_t EvalLikeMod4(Double_t mu, Int_t x, Int_t y, Double_t tau, Int_t what);
0085    Double_t LikeMod4(Double_t mu, Double_t b, Int_t x, Int_t y, Double_t tau);
0086
0087    //MODEL 5
0088    Double_t EvalLikeMod5(Double_t mu, Int_t x, Double_t bm, Double_t sdb, Int_t what);
0089    Double_t LikeMod5(Double_t mu, Double_t b, Int_t x, Double_t bm, Double_t u);
0090
0091    //MODEL 6
0092    Double_t EvalLikeMod6(Double_t mu, Int_t x, Int_t z, Double_t b, Int_t m, Int_t what);
0093    Double_t LikeMod6(Double_t mu, Double_t b, Double_t e, Int_t x, Int_t z, Int_t m);
0094
0095    //MODEL 7
0096    Double_t EvalLikeMod7(Double_t mu, Int_t x, Double_t em, Double_t sde, Double_t b, Int_t what);
0097    Double_t LikeMod7(Double_t mu, Double_t b, Double_t e, Int_t x, Double_t em, Double_t v);
0098
0099    //MISC
0100    static Double_t EvalPolynomial(Double_t x, const Int_t coef[], Int_t N);
0101    static Double_t EvalMonomial(Double_t x, const Int_t coef[], Int_t N);
0102    Double_t LogFactorial(Int_t n);
0103
0104    Double_t ComputeInterval(Int_t x, Int_t y, Int_t z, Double_t bm, Double_t em, Double_t e, Int_t mid, Double_t sde, Double_t sdb, Double_t tau, Double_t b, Int_t m);
0105
0106    void SetModelParameters(Int_t x, Int_t y, Int_t z, Double_t bm, Double_t em, Double_t e, Int_t mid, Double_t sde, Double_t sdb, Double_t tau, Double_t b, Int_t m);
0107
0108    void SetModelParameters();
0109
0110    Double_t GetBackground();
0111
0112 public:
0113
0114    /* Constructor */
0115    TRolke(Double_t CL = 0.9, Option_t *option = "");
0116
0117    /* Destructor */
0118    ~TRolke() override;
0119
0120    /* Get and set the Confidence Level */
0121    Double_t GetCL() const         {
0122       return fCL;
0123    }
0124    void     SetCL(Double_t CL)  {
0125       fCL = CL;
0126    }
0127
0128    /* Set the Confidence Level in terms of Sigmas. */
0129    void SetCLSigmas(Double_t CLsigmas) {
0130       fCL = TMath::Erf(CLsigmas / TMath::Sqrt(2.0)) ;
0131    }
0132
0133    // The Set methods for the different models are described in Rolke.cxx
0134    // model 1
0135    void SetPoissonBkgBinomEff(Int_t x, Int_t y, Int_t z, Double_t tau, Int_t m);
0136
0137    // model 2
0138    void SetPoissonBkgGaussEff(Int_t x, Int_t y, Double_t em, Double_t tau, Double_t sde);
0139
0140    // model 3
0141    void SetGaussBkgGaussEff(Int_t x, Double_t bm, Double_t em, Double_t sde, Double_t sdb);
0142
0143    // model 4
0144    void SetPoissonBkgKnownEff(Int_t x, Int_t y, Double_t tau, Double_t e);
0145
0146    // model 5
0147    void SetGaussBkgKnownEff(Int_t x, Double_t bm, Double_t sdb, Double_t e);
0148
0149    // model 6
0150    void SetKnownBkgBinomEff(Int_t x, Int_t z, Int_t m, Double_t b);
0151
0152    // model 7
0153    void SetKnownBkgGaussEff(Int_t x, Double_t em, Double_t sde, Double_t b);
0154
0155    /* Deprecated interface method (read Rolke.cxx). May be removed from future releases */
0156    Double_t CalculateInterval(Int_t x, Int_t y, Int_t z, Double_t bm, Double_t em, Double_t e, Int_t mid, Double_t sde, Double_t sdb, Double_t tau, Double_t b, Int_t m);
0157
0158    // get the upper and lower limits based on the specified model
0159    bool GetLimits(Double_t& low, Double_t& high);
0160    Double_t GetUpperLimit();
0161    Double_t GetLowerLimit();
0162
0163    // get the upper and lower average limits
0164    bool GetSensitivity(Double_t& low, Double_t& high, Double_t pPrecision = 0.00001);
0165
0166    // get the upper and lower limits for the outcome corresponding to
0167    // a given quantile.
0168    bool GetLimitsQuantile(Double_t& low, Double_t& high, Int_t& out_x, Double_t integral = 0.5);
0169
0170    // get the upper and lower limits for the most likely outcome.
0171    bool GetLimitsML(Double_t& low, Double_t& high, Int_t& out_x);
0172
0173    // get the value of x corresponding to rejection of the null hypothesis.
0174    bool GetCriticalNumber(Int_t& ncrit,Int_t maxtry=-1);
0175
0176    /* Get the bounding mode flag. True activates bounded mode. Read
0177       TRolke.cxx and the references therein for details. */
0178    bool GetBounding() const {
0179       return fBounding;
0180    }
0181
0182    /* Get the bounding mode flag. True activates bounded mode. Read
0183       TRolke.cxx and the references therein for details. */
0184    void SetBounding(const bool bnd) {
0185       fBounding = bnd;
0186    }
0187
0188    /* Deprecated name for SetBounding. */
0189    void SetSwitch(bool bnd) ;
0190
0191    /* Dump internals. Option is not used */
0192    void Print(Option_t*) const override;
0193
0194    ClassDefOverride(TRolke, 2)
0195 };
0196
0197 //calculate confidence limits using the Rolke method
0198 #endif
0199