Re: generating non-IID random sequences
- To: mathgroup at smc.vnet.net
- Subject: [mg79169] Re: generating non-IID random sequences
- From: Mark Fisher <particlefilter at gmail.com>
- Date: Fri, 20 Jul 2007 03:24:41 -0400 (EDT)
- References: <f7kegd$58u$1@smc.vnet.net>
On Jul 18, 3:12 am, Yaroslav Bulatov <yarosla... at gmail.com> wrote:
> I'm looking for a fast way to sample from a Markov-1 sequence of
> random bits. The method below is 3600 times slower than built-in
> RandomInteger function, can it be made much faster?
>
> p = 0.9; (* the probability of encountering 00 or 11 *)
> f = RandomChoice[{p^# (1 - p)^(1 - #), p^(1 - #) (1 - p)^#} -> {1, 0}]
> &;
> NestList[f, RandomChoice[{0, 1}], 100000]
Here are 3 functions. The first is your version; the second is your
version compiled; the third is a slightly different version also
compiled.
fun = Function[{n, p},
NestList[
RandomChoice[{p^# (1 - p)^(1 - #), p^(1 - #) (1 - p)^#} -> {1,
0}] &, RandomChoice[{0, 1}], n]];
fun1 = Compile[{{n, _Integer}, {p, _Real}},
NestList[
RandomChoice[{p^# (1 - p)^(1 - #), p^(1 - #) (1 - p)^#} -> {1,
0}] &, RandomChoice[{0, 1}], n]];
fun2 = Compile[{{n, _Integer}, {p, _Real}},
NestList[
If[# == 0, If[RandomReal[] < p, 0, 1],
If[RandomReal[] < 1 - p, 0, 1]] &, RandomChoice[{0, 1}], n]];
This gives the relative timings:
First[#]/#& @ Table[Timing[i[10^5, .9]][[1]], {i, {fun, fun1, fun2}}]
On my machine I get a factor of 5 speed up for fun1 and a factor of 25
speed up for fun2 relative to fun.
--Mark.