using TheRandomMonad.hs

TheRandomMonad.hs

import Control.Monad.Random
import             Data.List.Split as LS

bounds = [(10.0, 20.0), (5.0, 6.0), (3.0, 4.0), (100.0, 200.0)]

-- INPUT:
-- RandomGen          random number generator for uniforms
-- [(Double, Double)] list of upper and lower bounds
-- Int                number of lists of random numbers
--
-- OUTPUT:
-- List of list of random numbers, where the first element is sampled
--  uniformly from the bounds given in the first tuple, and so on
--
--
-- E.g. for input [(0,1), (3,5), (100,200)], I would expect output:
-- [ [0.5, 4.0, 150.0], [.23, 4.23, 122.123], ...]
make_randoms :: (RandomGen g, Monad m) => g -> [(Double, Double)] -> Int -> RandT g m [[Double]]
make_randoms gen bounds n = do 
  y  <-  unifn n
  let dim         = length bounds
  let n_samples   = n * dim
  let y_finite    = take n_samples y :: [Double]
  let xis         = LS.chunksOf dim y_finite
  let btrans xi   = zipWith (\a (lb, ub) -> lb + (a * ub)) xi bounds 
  return $ map btrans xis
    

make_10x2 :: (RandomGen g, Monad m) =>  RandT g m [[Double]]
make_10x2 = do
  l1 <- unifn 10
  l2 <- unifn 10
  return $ [l1, l2]

die :: (RandomGen g) => Rand g Int
die = getRandomR (1,6)

dice :: (RandomGen g) => Int -> Rand g [Int]
dice n = sequence (replicate n die)

dieT :: (RandomGen g, Monad m) => RandT g m Int
dieT = getRandomR (1,6)

diceT :: (RandomGen g, Monad m) => Int -> RandT g m [Int]
diceT n = sequence (replicate n dieT)

unif :: (RandomGen g, Monad m) => RandT g m Double
unif = getRandomR (0,1)

unifn :: (RandomGen g, Monad m) => Int -> RandT g m [Double]
unifn n = sequence (replicate n unif)

average :: [Double] -> Double
average x = (1.0/(fromIntegral (length x))) * (sum  x)

average_sim :: (RandomGen g, Monad m) => Int -> RandT g m Double
average_sim n = do
  y <- unifn n
  return $ average y

run_many_sims :: (RandomGen g, Monad m) => Int -> RandT g m [Double]
run_many_sims nsim = sequence ( replicate nsim  (average_sim 100))
  
run_many_sims_2 :: (RandomGen g, Monad m) => Int -> RandT g m [Double]
run_many_sims_2 nsim = do
  y <- sequence $ replicate nsim  $ average_sim 100
  return y


-- run_many_sims_3 :: (RandomGen g, Monad m) => Int -> RandT g m [Double]
-- run_many_sims_3 nsim = res
--   where                        
--     output_multiplier = 10.0
--     do
--       y <- sequence $ replicate nsim  $ average_sim 100
--       return y
--     res = output_multiplier * y



main = do

  print "Rand in IO"
  run1  <- evalRandIO (dice 2)
  run2  <- evalRandIO (dice 3)
  print run1
  print run2

  print "RandT Dice"
  gen <- newStdGen
  (y1, g')  <- runRandT (diceT 10) gen
  (y2, g2)  <- runRandT (diceT 10) g'
  print y1
  print y2

  print "RandT Uniforms"
  (u1, g3) <- runRandT (unifn 10) g2
  (u2, g4) <- runRandT (unifn 10) g3
  print u1
  print u2

  print "RandT Bounded Lists"
  (b1, g5) <- runRandT (make_randoms gen bounds 10 ) g4
  print b1