rpns
rpns

Reputation: 13

Iterate and getStdRandom

I have the following code (trying to simulate a Gaussian random walk):

import System.Random
import Control.Applicative

getStdNormal :: StdGen -> (Float,StdGen)
getStdNormal g = let l = drop 1 . take 49 $ iterate (\(_,g') -> random g') (0,g) in
                ((/ 2) . (+ (-24)) . sum . map fst $ l, snd . last $ l)

mphi = getStdRandom getStdNormal

l1 = setStdGen (mkStdGen 20180916) >> (sequence . take 10 $
                iterate (\s -> (\a -> \b ->  b) <$> s <*> mphi) (return 0))

l2 = setStdGen (mkStdGen 20180916) >> (sequence . take 10 $
                iterate (\s -> (\a -> \b -> a + b) <$> s <*> mphi) (return 0))

l1 is evaluated to the list of (apparently normal) random values, just as I can expect:

[0.0,1.1155739,-0.24667645,0.7793722,-0.18391132,0.23517609,-0.80208874,-1.5305595,-0.28670216,0.53894806]

However, l2 behaves far from my expectations. The list obviously is not the series of l1 partial sums, and is very unlikely to have standard normal differences:

[0.0,1.1155739,0.55951977,0.22015285,-2.0858078,0.6170025,-5.20298,0.3877325,1.410594,0.8647003]

(Yet note that the first two members are correct.) In fact, I can't explain where do the numbers in l2 come from, although I tend to blame laziness in mphi.

Why doesn't the code work as intended? Many thanks!

Upvotes: 1

Views: 86

Answers (1)

Daniel Wagner
Daniel Wagner

Reputation: 153182

Here are the first few elements of iterate (\s -> (\a b -> b) <$> s <*> mphi) (return 0):

return 0
(\a b -> b) <$> return 0 <*> mphi
(\a b -> b) <$> ((\a b -> b) <$> return 0 <*> mphi) <*> mphi
(\a b -> b) <$> ((\a b -> b) <$> ((\a b -> b) <$> return 0 <*> mphi) <*> mphi) <*> mphi

When you sequence these, here's what you get:

do
    v0 <- return 0
    v1 <- (\a b -> b) <$> return 0 <*> mphi
    v2 <- (\a b -> b) <$> ((\a b -> b) <$> return 0 <*> mphi) <*> mphi
    v3 <- (\a b -> b) <$> ((\a b -> b) <$> ((\a b -> b) <$> return 0 <*> mphi) <*> mphi) <*> mphi
    return [v0, v1, v2, v3]

Notice in particular that the random seed is not reset between these lines, and the things being thrown away by \a b -> b are result values, not IO actions. This means that we are getting the constant 0 in v0; the first element of the random sequence in v1 (not v0!); the third element of the random sequence (completely skipping the second element) in v2 (this one is probably what you wanted, but just by accident...); the sixth element of the random sequence in v3 (not v5! and skipping the fourth and fifth elements of the sequence entirely); and so on up the triangle numbers, throwing away most of the random numbers we generate.

What you almost certainly wanted instead was

do
    v0 <- mphi
    v1 <- mphi
    v2 <- mphi
    v3 <- mphi
    return [v0, v1, v2, v3]

which can be achieved by using repeat mphi instead of iterate (\s -> ...) (return 0). This lands us at:

sequence . take 10 $ repeat mphi

Idiomatically, one would then change take 10 (repeat mphi) to replicate 10 mphi, then change sequence (replicate 10 mphi) to replicateM 10 mphi.

Similar comments apply to your summing version: element n of the resulting list has generated n random numbers (distinct from any of the previously used ones for sums earlier in the result list) and summed them. Maybe this is what you want, but I suspect what you really wanted was just:

scanl (+) 0 <$> replicateM 10 mphi

Upvotes: 1

Related Questions