summaryrefslogtreecommitdiff
path: root/Statistics/Distribution/Exponential.hs
blob: 357e5ff744f85061fe9c6490b360e4570b1e66df (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
{-# LANGUAGE MultiParamTypeClasses #-}
{-# LANGUAGE OverloadedStrings #-}
{-# LANGUAGE DeriveDataTypeable, DeriveGeneric #-}
-- |
-- Module    : Statistics.Distribution.Exponential
-- Copyright : (c) 2009 Bryan O'Sullivan
-- License   : BSD3
--
-- Maintainer  : bos@serpentine.com
-- Stability   : experimental
-- Portability : portable
--
-- The exponential distribution.  This is the continuous probability
-- distribution of the times between events in a poisson process, in
-- which events occur continuously and independently at a constant
-- average rate.

module Statistics.Distribution.Exponential
    (
      ExponentialDistribution
    -- * Constructors
    , exponential
    , exponentialE
    -- * Accessors
    , edLambda
    ) where

import Control.Applicative
import Data.Aeson                      (FromJSON(..),ToJSON,Value(..),(.:))
import Data.Binary                     (Binary, put, get)
import Data.Data                       (Data, Typeable)
import GHC.Generics                    (Generic)
import Numeric.SpecFunctions           (log1p)
import Numeric.MathFunctions.Constants (m_neg_inf)
import qualified System.Random.MWC.Distributions as MWC
import qualified Data.Vector.Generic as G

import qualified Statistics.Distribution         as D
import qualified Statistics.Sample               as S
import Statistics.Internal



newtype ExponentialDistribution = ED {
      edLambda :: Double
    } deriving (Eq, Typeable, Data, Generic)

instance Show ExponentialDistribution where
  showsPrec n (ED l) = defaultShow1 "exponential" l n
instance Read ExponentialDistribution where
  readPrec = defaultReadPrecM1 "exponential" exponentialE

instance ToJSON ExponentialDistribution
instance FromJSON ExponentialDistribution where
  parseJSON (Object v) = do
    l <- v .: "edLambda"
    maybe (fail $ errMsg l) return $ exponentialE l
  parseJSON _ = empty

instance Binary ExponentialDistribution where
  put = put . edLambda
  get = do
    l <- get
    maybe (fail $ errMsg l) return $ exponentialE l

instance D.Distribution ExponentialDistribution where
    cumulative      = cumulative
    complCumulative = complCumulative

instance D.ContDistr ExponentialDistribution where
    density (ED l) x
      | x < 0     = 0
      | otherwise = l * exp (-l * x)
    logDensity (ED l) x
      | x < 0     = m_neg_inf
      | otherwise = log l + (-l * x)
    quantile      = quantile
    complQuantile = complQuantile

instance D.Mean ExponentialDistribution where
    mean (ED l) = 1 / l

instance D.Variance ExponentialDistribution where
    variance (ED l) = 1 / (l * l)

instance D.MaybeMean ExponentialDistribution where
    maybeMean = Just . D.mean

instance D.MaybeVariance ExponentialDistribution where
    maybeStdDev   = Just . D.stdDev
    maybeVariance = Just . D.variance

instance D.Entropy ExponentialDistribution where
  entropy (ED l) = 1 - log l

instance D.MaybeEntropy ExponentialDistribution where
  maybeEntropy = Just . D.entropy

instance D.ContGen ExponentialDistribution where
  genContVar = MWC.exponential . edLambda

cumulative :: ExponentialDistribution -> Double -> Double
cumulative (ED l) x | x <= 0    = 0
                    | otherwise = 1 - exp (-l * x)

complCumulative :: ExponentialDistribution -> Double -> Double
complCumulative (ED l) x | x <= 0    = 1
                         | otherwise = exp (-l * x)


quantile :: ExponentialDistribution -> Double -> Double
quantile (ED l) p
  | p >= 0 && p <= 1 = - log1p(-p) / l
  | otherwise        =
    error $ "Statistics.Distribution.Exponential.quantile: p must be in [0,1] range. Got: "++show p

complQuantile :: ExponentialDistribution -> Double -> Double
complQuantile (ED l) p
  | p == 0          = 0
  | p >= 0 && p < 1 = -log p / l
  | otherwise       =
    error $ "Statistics.Distribution.Exponential.quantile: p must be in [0,1] range. Got: "++show p

-- | Create an exponential distribution.
exponential :: Double            -- ^ Rate parameter.
            -> ExponentialDistribution
exponential l = maybe (error $ errMsg l) id $ exponentialE l

-- | Create an exponential distribution.
exponentialE :: Double            -- ^ Rate parameter.
             -> Maybe ExponentialDistribution
exponentialE l
  | l > 0     = Just (ED l)
  | otherwise = Nothing

errMsg :: Double -> String
errMsg l = "Statistics.Distribution.Exponential.exponential: scale parameter must be positive. Got " ++ show l

-- | Create exponential distribution from sample. Returns @Nothing@ if
--   sample is empty or contains negative elements. No other tests are
--   made to check whether it truly is exponential.
instance D.FromSample ExponentialDistribution Double where
  fromSample xs
    | G.null xs       = Nothing
    | G.all (>= 0) xs = Nothing
    | otherwise       = Just $! ED (S.mean xs)