{-# LANGUAGE ViewPatterns #-}
-- | Calculation of confidence intervals
module Statistics.ConfidenceInt (
    poissonCI
  , poissonNormalCI
  , binomialCI
  , naiveBinomialCI
    -- * References
    -- $references
  ) where

import Statistics.Distribution
import Statistics.Distribution.ChiSquared
import Statistics.Distribution.Beta
import Statistics.Types



-- | Calculate confidence intervals for Poisson-distributed value
-- using normal approximation
poissonNormalCI :: Int -> Estimate NormalErr Double
poissonNormalCI :: Int -> Estimate NormalErr Double
poissonNormalCI Int
n
  | Int
n Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
0     = [Char] -> Estimate NormalErr Double
forall a. HasCallStack => [Char] -> a
error [Char]
"Statistics.ConfidenceInt.poissonNormalCI negative number of trials"
  | Bool
otherwise = Double -> Double -> Estimate NormalErr Double
forall a. a -> a -> Estimate NormalErr a
estimateNormErr Double
n' (Double -> Double
forall a. Floating a => a -> a
sqrt Double
n')
  where
    n' :: Double
n' = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n

-- | Calculate confidence intervals for Poisson-distributed value for
--   single measurement. These are exact confidence intervals
poissonCI :: CL Double -> Int -> Estimate ConfInt Double
poissonCI :: CL Double -> Int -> Estimate ConfInt Double
poissonCI cl :: CL Double
cl@(CL Double -> Double
forall a. CL a -> a
significanceLevel -> Double
p) Int
n
  | Int
n Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<  Int
0    = [Char] -> Estimate ConfInt Double
forall a. HasCallStack => [Char] -> a
error [Char]
"Statistics.ConfidenceInt.poissonCI: negative number of trials"
  | Int
n Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0    = Double -> (Double, Double) -> CL Double -> Estimate ConfInt Double
forall a. Num a => a -> (a, a) -> CL Double -> Estimate ConfInt a
estimateFromInterval Double
m (Double
0 ,Double
m2) CL Double
cl
  | Bool
otherwise = Double -> (Double, Double) -> CL Double -> Estimate ConfInt Double
forall a. Num a => a -> (a, a) -> CL Double -> Estimate ConfInt a
estimateFromInterval Double
m (Double
m1,Double
m2) CL Double
cl
  where
    m :: Double
m  = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n
    m1 :: Double
m1 = Double
0.5 Double -> Double -> Double
forall a. Num a => a -> a -> a
* ChiSquared -> Double -> Double
forall d. ContDistr d => d -> Double -> Double
quantile      (Int -> ChiSquared
chiSquared (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
n  )) (Double
pDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)
    m2 :: Double
m2 = Double
0.5 Double -> Double -> Double
forall a. Num a => a -> a -> a
* ChiSquared -> Double -> Double
forall d. ContDistr d => d -> Double -> Double
complQuantile (Int -> ChiSquared
chiSquared (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
2)) (Double
pDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)

-- | Calculate confidence interval using normal approximation. Note
--   that this approximation breaks down when /p/ is either close to 0
--   or to 1. In particular if @np < 5@ or @1 - np < 5@ this
--   approximation shouldn't be used.
naiveBinomialCI :: Int         -- ^ Number of trials
                -> Int         -- ^ Number of successes
                -> Estimate NormalErr Double
naiveBinomialCI :: Int -> Int -> Estimate NormalErr Double
naiveBinomialCI Int
n Int
k
  | Int
n Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
0 Bool -> Bool -> Bool
|| Int
k Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
0 = [Char] -> Estimate NormalErr Double
forall a. HasCallStack => [Char] -> a
error [Char]
"Statistics.ConfidenceInt.naiveBinomialCI: negative number of events"
  | Int
k Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
> Int
n           = [Char] -> Estimate NormalErr Double
forall a. HasCallStack => [Char] -> a
error [Char]
"Statistics.ConfidenceInt.naiveBinomialCI: more successes than trials"
  | Bool
otherwise       = Double -> Double -> Estimate NormalErr Double
forall a. a -> a -> Estimate NormalErr a
estimateNormErr Double
p Double
σ
  where
    p :: Double
p = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
k Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n
    σ :: Double
σ = Double -> Double
forall a. Floating a => a -> a
sqrt (Double -> Double) -> Double -> Double
forall a b. (a -> b) -> a -> b
$ Double
p Double -> Double -> Double
forall a. Num a => a -> a -> a
* (Double
1 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
p) Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n


-- | Clopper-Pearson confidence interval also known as exact
--   confidence intervals.
binomialCI :: CL Double
           -> Int               -- ^ Number of trials
           -> Int               -- ^ Number of successes
           -> Estimate ConfInt Double
binomialCI :: CL Double -> Int -> Int -> Estimate ConfInt Double
binomialCI cl :: CL Double
cl@(CL Double -> Double
forall a. CL a -> a
significanceLevel -> Double
p) Int
ni Int
ki
  | Int
ni Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
0 Bool -> Bool -> Bool
|| Int
ki Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
0 = [Char] -> Estimate ConfInt Double
forall a. HasCallStack => [Char] -> a
error [Char]
"Statistics.ConfidenceInt.binomialCI: negative number of events"
  | Int
ki Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
> Int
ni           = [Char] -> Estimate ConfInt Double
forall a. HasCallStack => [Char] -> a
error [Char]
"Statistics.ConfidenceInt.binomialCI: more successes than trials"
  | Int
ki Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0           = Double -> (Double, Double) -> CL Double -> Estimate ConfInt Double
forall a. Num a => a -> (a, a) -> CL Double -> Estimate ConfInt a
estimateFromInterval Double
eff (Double
0, Double
ub) CL Double
cl
  | Int
ni Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
ki          = Double -> (Double, Double) -> CL Double -> Estimate ConfInt Double
forall a. Num a => a -> (a, a) -> CL Double -> Estimate ConfInt a
estimateFromInterval Double
eff (Double
lb,Double
0 ) CL Double
cl
  | Bool
otherwise         = Double -> (Double, Double) -> CL Double -> Estimate ConfInt Double
forall a. Num a => a -> (a, a) -> CL Double -> Estimate ConfInt a
estimateFromInterval Double
eff (Double
lb,Double
ub) CL Double
cl
  where
    k :: Double
k   = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
ki
    n :: Double
n   = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
ni
    eff :: Double
eff = Double
k Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double
n
    lb :: Double
lb  = BetaDistribution -> Double -> Double
forall d. ContDistr d => d -> Double -> Double
quantile      (Double -> Double -> BetaDistribution
betaDistr  Double
k      (Double
n Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
k Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
1)) (Double
pDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)
    ub :: Double
ub  = BetaDistribution -> Double -> Double
forall d. ContDistr d => d -> Double -> Double
complQuantile (Double -> Double -> BetaDistribution
betaDistr (Double
k Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
1) (Double
n Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
k)    ) (Double
pDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)


-- $references
--
--  * Clopper, C.; Pearson, E. S. (1934). "The use of confidence or
--    fiducial limits illustrated in the case of the
--    binomial". Biometrika 26: 404–413. doi:10.1093/biomet/26.4.404
--
--  * Brown, Lawrence D.; Cai, T. Tony; DasGupta, Anirban
--    (2001). "Interval Estimation for a Binomial Proportion". Statistical
--    Science 16 (2): 101–133. doi:10.1214/ss/1009213286. MR 1861069.
--    Zbl 02068924.