Repeating a monadic computation multiple times and printing out results? (MonadRandom)

Question

Right now I'm using the MonadRandom library. I have a computation:

metroChain :: (RandomGen g) => Rand g Double

I'd like to perform it multiple times, and sequentially print out the results. Or rather, I'd like to create some kind of list of multiple computations.

To do it once, I would use

main = do
  result <- evalRandIO metroChain
  print result

or

main = evalRandIO metroChain >>= (\result -> print result)

However, I'm having a lot of trouble being able to print out an arbitrary (n) amount of metroChain results.

Each result should use the RandomGen given by the end of the last result...that's how MonadRandom is supposed to work, right?

I've looked into replicateM, fmap, and a bit into transformers (although i admit I can't seem to understand them enough to grasp their application to my problem).

Can anyone help me achieve the functionality I'm looking for? I feel like I'm missing something really simple. But I'm pretty new to Haskell.

Answer 1

replicateM is what you want in building up the random computation. (Unless tel's guess is correct.)

foo :: Int -> IO ()
foo n = do
    results <- evalRandIO (replicateM n metroStep)
    mapM_ print results

Then you want mapM_ to aid in actually printing out the results.

Does this do what you want? Is there any of this you'd like me to expand on?

Answer 2

I'm going to make a leap and assume that metroStep is a MCMC Metropolis-Hastings iteration.

The problem you have is that you want the MH steps to be Markovian, but simply sharing RandomGen state, which is exactly what replicateM n metroStep does, is insufficient. That only makes it so that each step is capable of being based on independent random variables. To compare, if the RandomGen state weren't shared then immutability would guarantee that every metroStep is identical.

So what you really need is something that has both RandomGen state in order to provide an chain of psuedorandom numbers for generating independent variable samples and a fixed state so that at each step you can have P(x_i | theta, x_(i-1)). We build a transformer stack to do this—I'll use the mtl library and random-fu because I just wrote an MCMC using those libraries a few days ago.

metroStep :: (MonadRandom m, MonadState StateSpace m) => m StateSpace

where StateSpace is a point in state space including both observed an unobserved variables—it's every parameter on the right side of your likelihood function. Now, replicateM n metroStep :: (MonadRandom m, MonadState StateSpace m) => m [StateSpace] is a list of Markov-sequential StateSpace points.

Then we "run" a concrete version of this monad stack like this

do steps <- (`runRVar` StdRandom) . (`evalStateT` ss0) $ (replicateM n metroStep)
   mapM_ print steps