mwc-probability 1.3.0 → 2.0.0

raw patch · 4 files changed

+164/−35 lines, 4 filesdep ~basenew-uploader

Dependency ranges changed: base

Files

CHANGELOG +11/−0
README.md +40/−0
mwc-probability.cabal +10/−7
src/System/Random/MWC/Probability.hs +103/−28

+ CHANGELOG view

@@ -0,0 +1,11 @@+# Changelog++	- 2.0.0 (2018-01-29)+	* Add Laplace and Zipf-Mandelbrot distribution+	* Rename `isoGauss` to `isoNormal` and `standard` to `standardNormal` to uniform naming scheme	+	* Divide Haddock in sections+	+	- 1.3.0 (2016-12-04)+	* Generalize a couple of samplers to use Traversable rather than lists.++

+ README.md view

@@ -0,0 +1,40 @@+# mwc-probability++[![Build Status](https://secure.travis-ci.org/jtobin/mwc-probability.png)](http://travis-ci.org/jtobin/mwc-probability)+[![Hackage Version](https://img.shields.io/hackage/v/mwc-probability.svg)](http://hackage.haskell.org/package/mwc-probability)+[![MIT License](https://img.shields.io/badge/license-MIT-blue.svg)](https://github.com/jtobin/mwc-probability/blob/master/LICENSE)++Sampling function-based probability distributions.++A simple probability distribution type, where distributions are characterized+by sampling functions.++This implementation is a thin layer over `mwc-random`, which handles RNG+state-passing automatically by using a `PrimMonad` like `IO` or `ST s` under+the hood.++Examples+--------++Transform a distribution's support while leaving its density structure+invariant:++    -- uniform over [0, 1] to uniform over [1, 2]+    succ <$> uniform++Sequence distributions together using bind:++    -- a beta-binomial composite distribution+    beta 1 10 >>= binomial 10++Use do-notation to build complex joint distributions from composable,+local conditionals:++    hierarchicalModel = do+      [c, d, e, f] <- replicateM 4 $ uniformR (1, 10)+      a <- gamma c d+      b <- gamma e f+      p <- beta a b+      n <- uniformR (5, 10)+      binomial n p+

mwc-probability.cabal view

@@ -1,13 +1,14 @@ name:                mwc-probability-version:             1.3.0+version:             2.0.0 homepage:            http://github.com/jtobin/mwc-probability license:             MIT license-file:        LICENSE-author:              Jared Tobin-maintainer:          jared@jtobin.ca+author:              Jared Tobin, Marco Zocca+maintainer:          jared@jtobin.ca, zocca.marco gmail category:            Math build-type:          Simple cabal-version:       >= 1.10+tested-with:         GHC == 8.0.2, GHC == 8.2.2                      synopsis:            Sampling function-based probability distributions. description: @@ -41,6 +42,8 @@   >   p <- beta a b   >   n <- uniformR (5, 10)   >   binomial n p+extra-source-files:  README.md+                     CHANGELOG  Source-repository head   Type:     git@@ -51,8 +54,8 @@   default-language:    Haskell2010   hs-source-dirs:      src   build-depends:-      base          >  4 && < 6-    , mwc-random    >  0.13 && < 0.14-    , primitive     >= 0.6 && < 1.0-    , transformers  >= 0.5 && < 1.0+      base          >=  4.8 && < 6+    , mwc-random    >   0.13 && < 0.14+    , primitive     >=  0.6 && < 1.0+    , transformers  >=  0.5 && < 1.0

src/System/Random/MWC/Probability.hs view

@@ -4,10 +4,10 @@  -- | -- Module: System.Random.MWC.Probability--- Copyright: (c) 2015 Jared Tobin+-- Copyright: (c) 2015-2017 Jared Tobin, Marco Zocca -- License: MIT ----- Maintainer: Jared Tobin <jared@jtobin.ca>+-- Maintainer: Jared Tobin <jared@jtobin.ca>, Marco Zocca <zocca.marco gmail> -- Stability: unstable -- Portability: ghc --@@ -42,34 +42,54 @@ -- For example, @'uniform'@ is a distribution over the 0-1 interval, but @fmap -- (+ 1) uniform@ is the translated distribution over the 1-2 interval. ----- >>> sample (fmap (+ 1) uniform) gen+-- >>> create >>= sample (fmap (+ 1) uniform) -- 1.5480073474340754+--+-- == Running the examples+--+-- In the following we will assume an interactive GHCi session; the user should first declare a random number generator: +--+-- >>> gen <- create+--+-- which will be reused throughout all examples.+-- Note: creating a random generator is an expensive operation, so it should be only performed once in the code (usually in the top-level IO action, e.g `main`). + module System.Random.MWC.Probability (     module MWC   , Prob(..)   , samples +  -- * Distributions+  -- ** Continuous-valued   , uniform   , uniformR-  , discreteUniform-  , categorical-  , standard   , normal+  , standardNormal+  , isoNormal       , logNormal   , exponential+  , laplace   , gamma   , inverseGamma+  , normalGamma+  , weibull   , chiSquare   , beta+  , student+  -- *** Dirichlet process   , dirichlet-  , symmetricDirichlet+  , symmetricDirichlet  +  -- ** Discrete-valued+  , discreteUniform+  , zipf+  , categorical   , bernoulli   , binomial   , multinomial-  , student-  , isoGauss-  , poisson+  , poisson  ++   ) where  import Control.Applicative@@ -89,7 +109,6 @@  -- | A probability distribution characterized by a sampling function. ----- >>> gen <- create -- >>> sample uniform gen -- 0.4208881170464097 newtype Prob m a = Prob { sample :: Gen (PrimState m) -> m a }@@ -102,13 +121,20 @@ samples n model gen = replicateM n (sample model gen) {-# INLINABLE samples #-} -instance Monad m => Functor (Prob m) where-  fmap h (Prob f) = Prob $ fmap h . f+instance Functor m => Functor (Prob m) where+  fmap h (Prob f) = Prob (fmap h . f)  instance Monad m => Applicative (Prob m) where-  pure  = return+  pure  = Prob . const . pure   (<*>) = ap +instance Monad m => Monad (Prob m) where+  return = pure+  m >>= h = Prob $ \g -> do+    z <- sample m g+    sample (h z) g+  {-# INLINABLE (>>=) #-}+ instance (Monad m, Num a) => Num (Prob m a) where   (+)         = liftA2 (+)   (-)         = liftA2 (-)@@ -117,13 +143,8 @@   signum      = fmap signum   fromInteger = pure . fromInteger -instance Monad m => Monad (Prob m) where-  return  = Prob . const . return-  m >>= h = Prob $ \g -> do-    z <- sample m g-    sample (h z) g-  {-# INLINABLE (>>=) #-} + instance MonadTrans Prob where   lift m = Prob $ const m @@ -137,7 +158,6 @@  -- | The uniform distribution over a type. -----   >>> gen <- create --   >>> sample uniform gen :: IO Double --   0.29308497534914946 --   >>> sample uniform gen :: IO Bool@@ -168,9 +188,9 @@  -- | The standard normal or Gaussian distribution (with mean 0 and standard --   deviation 1).-standard :: PrimMonad m => Prob m Double-standard = Prob MWC.Dist.standard-{-# INLINABLE standard #-}+standardNormal :: PrimMonad m => Prob m Double+standardNormal = Prob MWC.Dist.standard+{-# INLINABLE standardNormal #-}  -- | The normal or Gaussian distribution with a specified mean and standard --   deviation.@@ -188,6 +208,23 @@ exponential r = Prob $ MWC.Dist.exponential r {-# INLINABLE exponential #-} +-- | The Laplace distribution with provided location and scale parameters.+laplace :: (Floating a, Variate a, PrimMonad m) => a -> a -> Prob m a+laplace mu sigma = do+  u <- uniformR (-0.5, 0.5)+  let b = sigma / sqrt 2+  return $ mu - b * signum u * log (1 - 2 * abs u)+{-# INLINABLE laplace #-}  +++-- | The Weibull distribution with provided shape and scale parameters.+weibull :: (Floating a, Variate a, PrimMonad m) => a -> a -> Prob m a+weibull a b = do+  x <- uniform+  return $ (- 1/a * log (1 - x)) ** 1/b+{-# INLINABLE weibull #-}++ -- | The gamma distribution with shape parameter a and scale parameter b. -- --   This is the parameterization used more traditionally in frequentist@@ -204,6 +241,14 @@ inverseGamma a b = recip <$> gamma a b {-# INLINABLE inverseGamma #-} +-- | The Normal-Gamma distribution of parameters mu, lambda, a, b+normalGamma :: PrimMonad m => Double -> Double -> Double -> Double -> Prob m Double+normalGamma mu lambda a b = do+  tau <- gamma a b+  let xsd = sqrt $ 1 / (lambda * tau)+  normal mu xsd+{-# INLINABLE normalGamma #-}+ -- | The chi-square distribution. chiSquare :: PrimMonad m => Int -> Prob m Double chiSquare k = Prob $ MWC.Dist.chiSquare k@@ -257,11 +302,12 @@   normal m sd {-# INLINABLE student #-} --- | An isotropic or spherical Gaussian distribution.-isoGauss+-- | An isotropic or spherical Gaussian distribution with specified mean+-- vector and scalar standard deviation parameter.+isoNormal   :: (Traversable f, PrimMonad m) => f Double -> Double -> Prob m (f Double)-isoGauss ms sd = traverse (`normal` sd) ms-{-# INLINABLE isoGauss #-}+isoNormal ms sd = traverse (`normal` sd) ms+{-# INLINABLE isoNormal #-}  -- | The Poisson distribution. poisson :: PrimMonad m => Double -> Prob m Int@@ -278,3 +324,32 @@     _   -> error "categorical: invalid return value" {-# INLINABLE categorical #-} ++-- | The Zipf-Mandelbrot distribution, generated with the rejection+-- sampling algorithm X.6.1 shown in+-- L.Devroye, Non-Uniform Random Variate Generation.+--+-- The parameter should be positive, but values close to 1 should be+-- avoided as they are very computationally intensive. The following+-- code illustrates this behaviour.+-- +-- >>> samples 10 (zipf 1.1) gen+-- [11315371987423520,2746946,653,609,2,13,85,4,256184577853,50]+-- +-- >>> samples 10 (zipf 1.5) gen+-- [19,3,3,1,1,2,1,191,2,1]+zipf :: (PrimMonad m, Integral b) => Double -> Prob m b+zipf a = do+  let+    b = 2 ** (a - 1)+    go = do+        u <- uniform+        v <- uniform+        let xInt = floor (u ** (- 1 / (a - 1))) +            x = fromIntegral xInt+            t = (1 + 1 / x) ** (a - 1)+        if v * x * (t - 1) / (b - 1) <= t / b+          then return xInt+          else go+  go+{-# INLINABLE zipf #-}