# BlockRandom and NestList: how get random maps with the same noise?

I need to compare solutions in random systems with the same noise and different initial conditions.

Here I run twice a random map with the same initial condition, inside a BlockRandom. The output expected should be exactly the same orbit but they are different.

What have I done wrong?

map[x_, α_] := α/2 ArcTan[x] + 1/2 x

ClearAll[α]
noisymap[x_, α_, σ_] := map[x, α] + σ RandomReal[{-1, 1}]

ClearAll[noisy1, noisy2]
BlockRandom[
noisy1 = NestList[With [{α = 2, σ = 0.1}, noisymap[#, α, σ] &], -0.32, 1000];
noisy2 = NestList[With [{α = 2, σ = 0.1}, noisymap[#, α, σ] &], -0.32, 1000];
]

ListLinePlot[{noisy1, noisy2}]

ListLinePlot[Abs[noisy1 - noisy2]]

-
Add a SeedRandom[<some same value>] before each, should be peachy. –  ciao Mar 2 '14 at 23:07
Simple mistake: use two separate BlockRandoms, one for each definition. (And without producing any random numbers or reseeding the RNGs in between.) –  Oleksandr R. Mar 3 '14 at 1:55

BlockRandom simply isolates anything affecting random number generation, etc. that is changed (like seeding). If you do nothing, there's nothing to isolate, and if you do nothing, the results will be the same as outside such a block: the second list gets subsequent random variates from the stream.

Use SeedRandom like so (snippet of your code modified):

map[x_, \[Alpha]_] := \[Alpha]/2 ArcTan[x] + 1/2 x

ClearAll[\[Alpha]]
noisymap[x_, \[Alpha]_, \[Sigma]_] :=
map[x, \[Alpha]] + \[Sigma] RandomReal[{-1, 1}]

ClearAll[noisy1, noisy2]
BlockRandom[

SeedRandom[1];
noisy1 =
NestList[With[{\[Alpha] = 2, \[Sigma] = 0.1},
noisymap[#, \[Alpha], \[Sigma]] &], -0.32, 1000];
SeedRandom[1];
noisy2 =
NestList[With[{\[Alpha] = 2, \[Sigma] = 0.1},
noisymap[#, \[Alpha], \[Sigma]] &], -0.32, 1000];]

noisy1 == noisy2

(* True *)


If you want both to just have same noise that differs each run, substitute some random value for the SeedRandom of both, or generate noise once and add to both.

Here's an example of BlockRandom behavior:

(* set some seed*)
SeedRandom[1];
list1 = RandomInteger[10, 10];
Print["First ten with seed 1 is: ", list1];

(*reset seed, so stream will be *the same*)
SeedRandom[1];

(* get first five *)
list2 = RandomInteger[10, 5];
Print["First five with seed 1 is: ", list2]

(* make a block *)
BlockRandom[SeedRandom[1]; list3 = RandomInteger[10, 10];
Print["First ten with seed 1 in block: ", list3]]

(* we're out of block...*)
list4 = RandomInteger[10, 5];
Print["Next five from original stream: ", list4];

(*
First ten with seed 1 is: {1,4,0,7,0,0,8,6,0,4}
First five with seed 1 is: {1,4,0,7,0}
First ten with seed 1 in block: {1,4,0,7,0,0,8,6,0,4}
Next five from original stream: {0,8,6,0,4}
*)


Note how the variates within the block are completely independent from what's happening "outside": when we leave the block, the behavior and state pick up where we left off. This includes seeding, setting of generator type, etc.

-
@anna: of course you get the same results. That's expected. BlockRandom is for when you want to isolate some changes you make to the state of MMA's "random objects" that you don't want to be seen outside of the blocked code, like say you needed super-quality numbers via mersenne twister (slower generation) for some operation: you could BlockRandom, tell MMA to use mersenne, when the block is left, MMA reverts to whatever was in use before (Incuding all prior state(s)). I'll add an example to my answer, might help you... –  ciao Mar 2 '14 at 23:47
thanks, I'll wait for the example. Just to clarify: fixing the seed, do I generate the same succession of random reals for both maps (which is what I want to do) or the same random real number at each step? (and why doesn't stackexcange alert me when an answer is posted?) –  anna Mar 3 '14 at 0:06
@anna: Yes, when you use SeedRandom with some value, the calls to random generators will return the same stream, so you'd either seed with the same value before generating the noise for each, or more efficiently, generate the noise once and add to each list. No idea why you don't get notified... –  ciao Mar 3 '14 at 0:09
thanks for the example –  anna Mar 3 '14 at 0:29
FWIW, I think your explanation is wrong. @anna's surprise that BlockRandom didn't do what she expected seems warranted, except that she ought to have wrapped the definitions of noisy1 and noisy2 in separate BlockRandoms. –  Oleksandr R. Mar 3 '14 at 1:52