File: Laplace.hs

package info (click to toggle)
haskell-statistics 0.16.2.1-2
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 640 kB
  • sloc: haskell: 6,819; ansic: 35; python: 33; makefile: 9
file content (163 lines) | stat: -rw-r--r-- 5,234 bytes parent folder | download
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
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
{-# LANGUAGE MultiParamTypeClasses #-}
{-# LANGUAGE OverloadedStrings #-}
{-# LANGUAGE DeriveDataTypeable, DeriveGeneric #-}
-- |
-- Module    : Statistics.Distribution.Laplace
-- Copyright : (c) 2015 Mihai Maruseac
-- License   : BSD3
--
-- Maintainer  : mihai.maruseac@maruseac.com
-- Stability   : experimental
-- Portability : portable
--
-- The Laplace distribution.  This is the continuous probability
-- defined as the difference of two iid exponential random variables
-- or a Brownian motion evaluated as exponentially distributed times.
-- It is used in differential privacy (Laplace Method), speech
-- recognition and least absolute deviations method (Laplace's first
-- law of errors, giving a robust regression method)
--
module Statistics.Distribution.Laplace
    (
      LaplaceDistribution
    -- * Constructors
    , laplace
    , laplaceE
    -- * Accessors
    , ldLocation
    , ldScale
    ) where

import Control.Applicative
import Data.Aeson           (FromJSON(..), ToJSON, Value(..), (.:))
import Data.Binary          (Binary(..))
import Data.Data            (Data, Typeable)
import GHC.Generics         (Generic)
import qualified Data.Vector.Generic             as G
import qualified Statistics.Distribution         as D
import qualified Statistics.Quantile             as Q
import qualified Statistics.Sample               as S
import Statistics.Internal


data LaplaceDistribution = LD {
      ldLocation :: {-# UNPACK #-} !Double
    -- ^ Location.
    , ldScale    :: {-# UNPACK #-} !Double
    -- ^ Scale.
    } deriving (Eq, Typeable, Data, Generic)

instance Show LaplaceDistribution where
  showsPrec i (LD l s) = defaultShow2 "laplace" l s i
instance Read LaplaceDistribution where
  readPrec = defaultReadPrecM2 "laplace" laplaceE

instance ToJSON LaplaceDistribution
instance FromJSON LaplaceDistribution where
  parseJSON (Object v) = do
    l <- v .: "ldLocation"
    s <- v .: "ldScale"
    maybe (fail $ errMsg l s) return $ laplaceE l s
  parseJSON _ = empty

instance Binary LaplaceDistribution where
  put (LD l s) = put l >> put s
  get = do
    l <- get
    s <- get
    maybe (fail $ errMsg l s) return $ laplaceE l s

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

instance D.ContDistr LaplaceDistribution where
    density    (LD l s) x = exp (- abs (x - l) / s) / (2 * s)
    logDensity (LD l s) x = - abs (x - l) / s - log 2 - log s
    quantile      = quantile
    complQuantile = complQuantile

instance D.Mean LaplaceDistribution where
    mean (LD l _) = l

instance D.Variance LaplaceDistribution where
    variance (LD _ s) = 2 * s * s

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

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

instance D.Entropy LaplaceDistribution where
  entropy (LD _ s) = 1 + log (2 * s)

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

instance D.ContGen LaplaceDistribution where
  genContVar = D.genContinuous

cumulative :: LaplaceDistribution -> Double -> Double
cumulative (LD l s) x
  | x <= l    = 0.5 * exp ( (x - l) / s)
  | otherwise = 1 - 0.5 * exp ( - (x - l) / s )

complCumulative :: LaplaceDistribution -> Double -> Double
complCumulative (LD l s) x
  | x <= l    = 1 - 0.5 * exp ( (x - l) / s)
  | otherwise = 0.5 * exp ( - (x - l) / s )

quantile :: LaplaceDistribution -> Double -> Double
quantile (LD l s) p
  | p == 0             = -inf
  | p == 1             = inf
  | p == 0.5           = l
  | p > 0   && p < 0.5 = l + s * log (2 * p)
  | p > 0.5 && p < 1   = l - s * log (2 - 2 * p)
  | otherwise          =
    error $ "Statistics.Distribution.Laplace.quantile: p must be in [0,1] range. Got: "++show p
  where
    inf = 1 / 0

complQuantile :: LaplaceDistribution -> Double -> Double
complQuantile (LD l s) p
  | p == 0             = inf
  | p == 1             = -inf
  | p == 0.5           = l
  | p > 0   && p < 0.5 = l - s * log (2 * p)
  | p > 0.5 && p < 1   = l + s * log (2 - 2 * p)
  | otherwise          =
    error $ "Statistics.Distribution.Laplace.quantile: p must be in [0,1] range. Got: "++show p
  where
    inf = 1 / 0

-- | Create an Laplace distribution.
laplace :: Double         -- ^ Location
        -> Double        -- ^ Scale
        -> LaplaceDistribution
laplace l s = maybe (error $ errMsg l s) id $ laplaceE l s

-- | Create an Laplace distribution.
laplaceE :: Double         -- ^ Location
         -> Double        -- ^ Scale
         -> Maybe LaplaceDistribution
laplaceE l s
  | s >= 0    = Just (LD l s)
  | otherwise = Nothing

errMsg :: Double -> Double -> String
errMsg _ s = "Statistics.Distribution.Laplace.laplace: scale parameter must be positive. Got " ++ show s


-- | Create Laplace distribution from sample.  The location is estimated
--   as the median of the sample, and the scale as the mean absolute
--   deviation of the median.
instance D.FromSample LaplaceDistribution Double where
  fromSample xs
    | G.null xs = Nothing
    | otherwise = Just $! LD s l
    where
      s = Q.median Q.medianUnbiased xs
      l = S.mean $ G.map (\x -> abs $ x - s) xs