Haskell软件包，用于从标准概率分布中采样

Question

我想在Haskell中进行一些蒙特卡洛分析。我希望能够编写如下代码：

do n <- poisson lambda
   xs <- replicateM n $ normal mu sigma
   return $ maximum xs

对应于随机模型

n ~ Poisson(lambda)
for (i in 1:n)
  x[i] ~ Normal(mu, sigma)
y = max_{i=1}^n x[i]

我可以很轻松地看到如何创建必要的随机采样单子。但是，我希望不必为所有标准概率分布实施采样器。是否有已实施这些的Haskell软件包？

我看过package random-fu，它在0.2.7版本中停滞了三年，但我无法理解。它取决于typeclass MonadRandom和RandomSource，并没有得到很好的解释。

我也看过package mwc-probability，但我也听不懂它-看来您已经了解了PrimMonadand PrimState类。

这两个软件包都令我感到惊讶，因为它们具有过于复杂的API，并且似乎完全放弃了Haskell的标准随机数生成框架（如中的所示）System.Random。

任何意见，将不胜感激。

Answer 1

好吧，如果您希望能够编写如下代码：

do n <- poisson lambda
   xs <- replicateM n $ normal mu sigma
   return $ maximum xs

那么您大概想使用random-fu：

import Control.Monad
import Data.Random
import Data.Random.Distribution.Poisson
import Data.Random.Distribution.Normal

foo :: RVar Double
foo = do
  n <- poisson lambda
  xs <- replicateM (n+1) $ normal mu sigma
  return $ maximum xs

  where lambda = 10 :: Double
        mu = 0
        sigma = 6

main :: IO ()
main = print =<< replicateM 10 (sample foo)

我不确定过去三年缺乏更新是否应该成为决定因素。是否真的已经是伽玛分布的世界上许多令人兴奋的进展？

实际上，看起来mwc-probability效果差不多：

import Control.Monad
import System.Random.MWC.Probability

foo :: Prob IO Double
foo = do
  n <- poisson lambda
  xs <- replicateM (n+1) $ normal mu sigma
  return $ maximum xs

  where lambda = 10 :: Double
        mu = 0
        sigma = 6

main :: IO ()
main = do
  gen <- createSystemRandom
  print =<< replicateM 10 (sample foo gen)