Haskell 伯努利分布布尔值的无限列表
Haskell infinite list of Bernoulli distributed booleans
我需要一个有偏见的随机布尔值列表。每个布尔值都需要具有相同的为真概率(伯努利分布)。这些布尔值被传递给一个函数,该函数为每个输入布尔值生成零个或多个输出布尔值。我需要一个无限列表,因为我事先不知道需要多少布尔值才能提供足够的输出。请参阅以下(简化的)代码:
import System.Random.MWC
import System.Random.MWC.Distributions
foo :: [Bool] -> [Bool] -- foo outputs zero or more Bools per input Bool
main = do
gen <- create
bits <- sequence . repeat $ bernoulli 0.25 gen
print . take 32 . foo $ bits
不幸的是,这段代码挂在 main
的第二行。我猜 Control.Monad.ST
?
某处发生了一些非懒惰的事情
(我可以用 System.Random.randoms
做这样的事情,但结果值没有所需的分布。)
我可以在继续使用 System.Random.MWC
库的同时解决这个问题吗?还是这需要我切换到其他实施方式?
罪魁祸首是 sequence . repeat
- 这将(几乎?)每个 monad 挂起,因为您必须执行可能无限数量的效果。
最简单的解决方案是使用不同的库 - 如果您依赖于 mwc-random 生成的数字的质量,这可能是不可能的。下一个最简单的解决方案是重写 foo
以具有类型 [IO Bool] -> IO [Bool]
并传递给它 repeat (bernoulli 0.25 gen)
- 这将允许 foo
选择何时停止执行由无限列表。但是把你的逻辑放在 IO
里面不是很好。
当您需要无限的随机数列表时,标准技巧是使用纯函数 f :: StdGen -> (Result, StdGen)
。然后unfoldr (Just . f) :: StdGen -> [Result]
,输出是一个无限列表。乍一看,mwc-random 似乎只有 monadic 函数,没有纯接口。然而,事实并非如此,因为 ST s
是 PrimMonad
的一个实例。您还可以使用将 Gen
转换为 Seed
的函数。使用这些,您可以获得任何单子函数的纯 RNG 函数:
{-# LANGUAGE RankNTypes #-}
import System.Random.MWC
import System.Random.MWC.Distributions
import Control.Monad.ST
import Data.List
pureRand :: (forall s . GenST s -> ST s t) -> Seed -> (t, Seed)
pureRand f s = runST $ do
s' <- restore s
r <- f s'
s'' <- save s'
return (r, s'')
pureBernoulli :: Double -> Seed -> (Bool, Seed)
pureBernoulli a = pureRand (bernoulli a)
foo :: [Bool] -> [Bool]
foo = id
main = do
gen <- create >>= save
let bits = unfoldr (Just . pureBernoulli 0.25) gen
print . take 32 . foo $ bits
不幸的是,mwc-random 默认情况下不公开此类接口,但很容易获得。
另一个选项稍微有点可怕——使用不安全的函数。
import System.IO.Unsafe
repeatM rand = go where
go = do
x <- rand
xs <- unsafeInterleaveIO go
return (x : xs)
main2 = do
gen <- create
bits <- repeatM (bernoulli 0.25 gen)
print . take 32 . foo $ bits
自然地,这伴随着围绕 unsafe
的常见警告 - 仅当纯函数给您带来极大不便时才使用它。 unsafeInterleaveIO
可能会重新排序或永远不会执行效果 - 例如,如果 foo
忽略一个元素,它将永远不会被计算并且更新存储在 gen
中的状态的相应效果可能不会发生。例如,以下将不打印任何内容:
snd <$> ((,) <$> unsafeInterleaveIO (putStrLn "Hello") <*> return ())
mwc-random
包提供了两个 PrimMonad
实例,一个用于 IO
,另一个用于 ST s
。只要 ST
计算在所有状态标签 s
上参数化,我们就可以 运行 计算并使用 runST :: (forall s. ST s a) -> a
. By itself this wouldn't be very useful since we'd lose the state: the seed of the random generator, but mwc-random also provides explicit ways to handle the seeds:
提取值
save :: PrimMonad m => Gen (PrimState m) -> m Seed
restore :: PrimMonad m => Seed -> m (Gen (PrimState m))
只要生成器在 forall s. ST s
.
中,我们就可以使用它们进行计算,从生成单个值的计算中生成值流
{-# LANGUAGE RankNTypes #-}
{-# LANGUAGE ScopedTypeVariables #-}
import System.Random.MWC
import Control.Monad.ST
import System.Random.MWC.Distributions
randomStream :: forall s a. (forall s. GenST s -> ST s a) -> GenST s -> ST s [a]
randomStream item = go
where
go :: forall s. GenST s -> ST s [a]
go gen = do
x <- item gen
seed <- save gen
return (x:runST (restore seed >>= go))
有了这个我们可以把你的例子写成
main = do
bits <- withSystemRandom (randomStream (bernoulli 0.25))
print . take 32 $ bits
我们实际上可以构建比为流中的每个项目使用相同生成器更复杂的生成器。我们可以沿着流线程化一个状态,这样每个值都可以依赖于前一个值的结果。
unfoldStream :: forall s a b. (forall s. b -> GenST s -> ST s (a, b)) -> b -> GenST s -> ST s [a]
unfoldStream item = go
where
go :: forall s. b -> GenST s -> ST s [a]
go b gen = do
(x,b') <- item b gen
seed <- save gen
return (x:runST (restore seed >>= go b'))
以下示例流的结果每次出现的可能性都会增加 False
。
import Control.Monad.Primitive
interesting :: (PrimMonad m) => Double -> Gen (PrimState m) -> m (Bool, Double)
interesting p gen = do
result <- bernoulli p gen
let p' = if result then p else p + (1-p)*0.25
return (result, p')
main = do
bits <- withSystemRandom (unfoldStream interesting 0)
print . take 32 $ bits
我需要一个有偏见的随机布尔值列表。每个布尔值都需要具有相同的为真概率(伯努利分布)。这些布尔值被传递给一个函数,该函数为每个输入布尔值生成零个或多个输出布尔值。我需要一个无限列表,因为我事先不知道需要多少布尔值才能提供足够的输出。请参阅以下(简化的)代码:
import System.Random.MWC
import System.Random.MWC.Distributions
foo :: [Bool] -> [Bool] -- foo outputs zero or more Bools per input Bool
main = do
gen <- create
bits <- sequence . repeat $ bernoulli 0.25 gen
print . take 32 . foo $ bits
不幸的是,这段代码挂在 main
的第二行。我猜 Control.Monad.ST
?
(我可以用 System.Random.randoms
做这样的事情,但结果值没有所需的分布。)
我可以在继续使用 System.Random.MWC
库的同时解决这个问题吗?还是这需要我切换到其他实施方式?
罪魁祸首是 sequence . repeat
- 这将(几乎?)每个 monad 挂起,因为您必须执行可能无限数量的效果。
最简单的解决方案是使用不同的库 - 如果您依赖于 mwc-random 生成的数字的质量,这可能是不可能的。下一个最简单的解决方案是重写 foo
以具有类型 [IO Bool] -> IO [Bool]
并传递给它 repeat (bernoulli 0.25 gen)
- 这将允许 foo
选择何时停止执行由无限列表。但是把你的逻辑放在 IO
里面不是很好。
当您需要无限的随机数列表时,标准技巧是使用纯函数 f :: StdGen -> (Result, StdGen)
。然后unfoldr (Just . f) :: StdGen -> [Result]
,输出是一个无限列表。乍一看,mwc-random 似乎只有 monadic 函数,没有纯接口。然而,事实并非如此,因为 ST s
是 PrimMonad
的一个实例。您还可以使用将 Gen
转换为 Seed
的函数。使用这些,您可以获得任何单子函数的纯 RNG 函数:
{-# LANGUAGE RankNTypes #-}
import System.Random.MWC
import System.Random.MWC.Distributions
import Control.Monad.ST
import Data.List
pureRand :: (forall s . GenST s -> ST s t) -> Seed -> (t, Seed)
pureRand f s = runST $ do
s' <- restore s
r <- f s'
s'' <- save s'
return (r, s'')
pureBernoulli :: Double -> Seed -> (Bool, Seed)
pureBernoulli a = pureRand (bernoulli a)
foo :: [Bool] -> [Bool]
foo = id
main = do
gen <- create >>= save
let bits = unfoldr (Just . pureBernoulli 0.25) gen
print . take 32 . foo $ bits
不幸的是,mwc-random 默认情况下不公开此类接口,但很容易获得。
另一个选项稍微有点可怕——使用不安全的函数。
import System.IO.Unsafe
repeatM rand = go where
go = do
x <- rand
xs <- unsafeInterleaveIO go
return (x : xs)
main2 = do
gen <- create
bits <- repeatM (bernoulli 0.25 gen)
print . take 32 . foo $ bits
自然地,这伴随着围绕 unsafe
的常见警告 - 仅当纯函数给您带来极大不便时才使用它。 unsafeInterleaveIO
可能会重新排序或永远不会执行效果 - 例如,如果 foo
忽略一个元素,它将永远不会被计算并且更新存储在 gen
中的状态的相应效果可能不会发生。例如,以下将不打印任何内容:
snd <$> ((,) <$> unsafeInterleaveIO (putStrLn "Hello") <*> return ())
mwc-random
包提供了两个 PrimMonad
实例,一个用于 IO
,另一个用于 ST s
。只要 ST
计算在所有状态标签 s
上参数化,我们就可以 运行 计算并使用 runST :: (forall s. ST s a) -> a
. By itself this wouldn't be very useful since we'd lose the state: the seed of the random generator, but mwc-random also provides explicit ways to handle the seeds:
save :: PrimMonad m => Gen (PrimState m) -> m Seed
restore :: PrimMonad m => Seed -> m (Gen (PrimState m))
只要生成器在 forall s. ST s
.
{-# LANGUAGE RankNTypes #-}
{-# LANGUAGE ScopedTypeVariables #-}
import System.Random.MWC
import Control.Monad.ST
import System.Random.MWC.Distributions
randomStream :: forall s a. (forall s. GenST s -> ST s a) -> GenST s -> ST s [a]
randomStream item = go
where
go :: forall s. GenST s -> ST s [a]
go gen = do
x <- item gen
seed <- save gen
return (x:runST (restore seed >>= go))
有了这个我们可以把你的例子写成
main = do
bits <- withSystemRandom (randomStream (bernoulli 0.25))
print . take 32 $ bits
我们实际上可以构建比为流中的每个项目使用相同生成器更复杂的生成器。我们可以沿着流线程化一个状态,这样每个值都可以依赖于前一个值的结果。
unfoldStream :: forall s a b. (forall s. b -> GenST s -> ST s (a, b)) -> b -> GenST s -> ST s [a]
unfoldStream item = go
where
go :: forall s. b -> GenST s -> ST s [a]
go b gen = do
(x,b') <- item b gen
seed <- save gen
return (x:runST (restore seed >>= go b'))
以下示例流的结果每次出现的可能性都会增加 False
。
import Control.Monad.Primitive
interesting :: (PrimMonad m) => Double -> Gen (PrimState m) -> m (Bool, Double)
interesting p gen = do
result <- bernoulli p gen
let p' = if result then p else p + (1-p)*0.25
return (result, p')
main = do
bits <- withSystemRandom (unfoldStream interesting 0)
print . take 32 $ bits