module Data.Random.Distribution.Poisson where
import Data.Random.Source
import Data.Random.Distribution
import Data.Random.RVar
import Data.Random.Distribution.Uniform
import Data.Random.Distribution.Gamma
import Data.Random.Distribution.Binomial
import Data.Int
import Data.Word
import Control.Monad
integralPoisson :: (Integral a, RealFloat b, Distribution StdUniform b, Distribution (Erlang a) b, Distribution (Binomial b) a) => b -> RVar a
integralPoisson mu = psn 0 mu
where
psn :: (Integral a, RealFloat b, Distribution StdUniform b, Distribution (Erlang a) b, Distribution (Binomial b) a) => a -> b -> RVar a
psn k mu
| mu > 10 = do
let m = floor (mu * (7/8))
x <- erlang m
if x >= mu
then do
b <- binomial (m 1) (mu / x)
return (k + b)
else psn (k + m) (mu x)
| otherwise = prod 1 k
where
emu = exp (mu)
prod p k = do
u <- stdUniform
if p * u > emu
then prod (p * u) (k + 1)
else return k
fractionalPoisson :: (Num a, Distribution (Poisson b) Integer) => b -> RVar a
fractionalPoisson mu = liftM fromInteger (poisson mu)
poisson :: (Distribution (Poisson b) a) => b -> RVar a
poisson mu = rvar (Poisson mu)
data Poisson b a = Poisson b
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Int ) b, Distribution (Binomial b) Int ) => Distribution (Poisson b) Int where rvar (Poisson mu) = integralPoisson mu
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Int8 ) b, Distribution (Binomial b) Int8 ) => Distribution (Poisson b) Int8 where rvar (Poisson mu) = integralPoisson mu
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Int16 ) b, Distribution (Binomial b) Int16 ) => Distribution (Poisson b) Int16 where rvar (Poisson mu) = integralPoisson mu
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Int32 ) b, Distribution (Binomial b) Int32 ) => Distribution (Poisson b) Int32 where rvar (Poisson mu) = integralPoisson mu
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Int64 ) b, Distribution (Binomial b) Int64 ) => Distribution (Poisson b) Int64 where rvar (Poisson mu) = integralPoisson mu
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Word8 ) b, Distribution (Binomial b) Word8 ) => Distribution (Poisson b) Word8 where rvar (Poisson mu) = integralPoisson mu
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Word16 ) b, Distribution (Binomial b) Word16 ) => Distribution (Poisson b) Word16 where rvar (Poisson mu) = integralPoisson mu
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Word32 ) b, Distribution (Binomial b) Word32 ) => Distribution (Poisson b) Word32 where rvar (Poisson mu) = integralPoisson mu
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Word64 ) b, Distribution (Binomial b) Word64 ) => Distribution (Poisson b) Word64 where rvar (Poisson mu) = integralPoisson mu
instance (RealFloat b, Distribution StdUniform b, Distribution (Erlang Integer ) b, Distribution (Binomial b) Integer ) => Distribution (Poisson b) Integer where rvar (Poisson mu) = integralPoisson mu
instance (Distribution (Poisson b) Integer) => Distribution (Poisson b) Float where rvar (Poisson mu) = fractionalPoisson mu
instance (Distribution (Poisson b) Integer) => Distribution (Poisson b) Double where rvar (Poisson mu) = fractionalPoisson mu