splitmix-distributions 0.1.0.0 → 0.3.0.0

raw patch · 3 files changed

+85/−20 lines, 3 filesPVP ok

version bump matches the API change (PVP)

API changes (from Hackage documentation)

Files

• ChangeLog.md +7/−2
• splitmix-distributions.cabal +2/−2
• src/System/Random/SplitMix/Distributions.hs +76/−16

@@ -1,3 +1,8 @@-# Changelog for splitmix-distributions+0.3 : -## Unreleased changes+- type signatures of all generators are now parametrized over some Monad m rather than Identity. This allows for more flexibility on the use site.++0.2 :++- add Pareto, Dirichlet, multinomial distributions+

@@ -1,5 +1,5 @@ name:           splitmix-distributions-version:        0.1.0.0+version:        0.3.0.0 description:    Random samplers for some common distributions, as well as a convenient interface for composing them, based on splitmix. Please see the README on GitHub at <https://github.com/ocramz/splitmix-distributions#readme> homepage:       https://github.com/ocramz/splitmix-distributions#readme bug-reports:    https://github.com/ocramz/splitmix-distributions/issues@@ -34,7 +34,7 @@     , transformers   default-language: Haskell2010 -test-suite splitmix-distributions-test+test-suite test   type: exitcode-stdio-1.0   main-is: Spec.hs   hs-source-dirs:

@@ -40,8 +40,11 @@   stdNormal, normal,   beta,   gamma,+  pareto,+  dirichlet,   -- ** Discrete-  bernoulli,+  bernoulli, fairCoin,+  multinomial,   -- * PRNG   -- ** Pure   Gen, sample,@@ -52,7 +55,9 @@  import Control.Monad (replicateM) import Control.Monad.IO.Class (MonadIO(..))+import Data.Foldable (toList) import Data.Functor.Identity (Identity(..))+import Data.List (findIndex) import GHC.Word (Word64)  -- erf@@ -89,30 +94,58 @@   -- | Bernoulli trial-bernoulli :: Double -- ^ bias parameter \( 0 \lt p \lt 1 \)-          -> Gen Bool+bernoulli :: Monad m =>+             Double -- ^ bias parameter \( 0 \lt p \lt 1 \)+          -> GenT m Bool bernoulli p = withGen (bernoulliF p) +-- | A fair coin toss returns either value with probability 0.5+fairCoin :: Monad m => GenT m Bool+fairCoin = bernoulli 0.5++-- | Multinomial distribution+--+-- NB : returns @Nothing@ if any of the input probabilities is negative+multinomial :: (Monad m, Foldable t) =>+               Int -- ^ number of Bernoulli trials \( n \gt 0 \)+            -> t Double -- ^ probability vector \( p_i \gt 0 , \forall i \) (does not need to be normalized)+            -> GenT m (Maybe [Int])+multinomial n ps = do+    let (cumulative, total) = runningTotals (toList ps)+    ms <- replicateM n $ do+      z <- uniformR 0 total+      pure $ findIndex (> z) cumulative+        -- Just g  -> return g+        -- Nothing -> error "splitmix-distributions: invalid probability vector"+    pure $ sequence ms+  where+    runningTotals :: Num a => [a] -> ([a], a)+    runningTotals xs = let adds = scanl1 (+) xs in (adds, sum xs)+{-# INLINABLE multinomial #-}++ -- | Uniform between two values-uniformR :: Double -- ^ low+uniformR :: Monad m =>+            Double -- ^ low          -> Double -- ^ high-         -> Gen Double+         -> GenT m Double uniformR lo hi = scale <$> stdUniform   where     scale x = x * (hi - lo) + lo --- | Standard normal-stdNormal :: Gen Double+-- | Standard normal distribution+stdNormal :: Monad m => GenT m Double stdNormal = normal 0 1  -- | Uniform in [0, 1)-stdUniform :: Gen Double+stdUniform :: Monad m => GenT m Double stdUniform = withGen nextDouble  -- | Beta distribution, from two standard uniform samples-beta :: Double -- ^ shape parameter \( \alpha \gt 0 \) +beta :: Monad m =>+        Double -- ^ shape parameter \( \alpha \gt 0 \)       -> Double -- ^ shape parameter \( \beta \gt 0 \)-     -> Gen Double+     -> GenT m Double beta a b = go   where     go = do@@ -128,9 +161,10 @@ -- | Gamma distribution, using Ahrens-Dieter accept-reject (algorithm GD): -- -- Ahrens, J. H.; Dieter, U (January 1982). "Generating gamma variates by a modified rejection technique". Communications of the ACM. 25 (1): 47–54-gamma :: Double -- ^ shape parameter \( k \gt 0 \)+gamma :: Monad m =>+         Double -- ^ shape parameter \( k \gt 0 \)       -> Double -- ^ scale parameter \( \theta \gt 0 \)-      -> Gen Double+      -> GenT m Double gamma k th = do   xi <- sampleXi   us <- replicateM n (log <$> stdUniform)@@ -153,16 +187,42 @@             let xi = 1 - log v             in (xi, w * exp (- xi)) +-- | Pareto distribution+pareto :: Monad m =>+          Double -- ^ shape parameter \( \alpha \gt 0 \)+       -> Double -- ^ scale parameter \( x_{min} \gt 0 \)+       -> GenT m Double+pareto a xmin = do+  y <- exponential a+  return $ xmin * exp y+{-# INLINABLE pareto #-} +-- | The Dirichlet distribution with the provided concentration parameters.+--   The dimension of the distribution is determined by the number of+--   concentration parameters supplied.+--+--   >>> sample 1234 (dirichlet [0.1, 1, 10])+--   [2.3781130220132788e-11,6.646079701567026e-2,0.9335392029605486]+dirichlet :: (Monad m, Traversable f) =>+             f Double -- ^ concentration parameters \( \gamma_i \gt 0 , \forall i \)+          -> GenT m (f Double)+dirichlet as = do+  zs <- traverse (`gamma` 1) as+  return $ fmap (/ sum zs) zs+{-# INLINABLE dirichlet #-}++ -- | Normal distribution-normal :: Double -- ^ mean+normal :: Monad m =>+          Double -- ^ mean        -> Double -- ^ standard deviation \( \sigma \gt 0 \)-       -> Gen Double+       -> GenT m Double normal mu sig = withGen (normalF mu sig)  -- | Exponential distribution-exponential :: Double -- ^ rate parameter \( \lambda > 0 \)-            -> Gen Double+exponential :: Monad m =>+               Double -- ^ rate parameter \( \lambda > 0 \)+            -> GenT m Double exponential l = withGen (exponentialF l)  -- | Wrap a 'splitmix' PRNG function