Re: [R] Random seed problem in MCMC coupling of chains

From: Gabor Grothendieck <ggrothendieck_at_gmail.com>
Date: Thu 09 Jun 2005 - 02:55:07 EST

That could be addressed like this (where changing the offset changes the experiment).

offset <- 123

niter <- 3
nchain <- 2
for (i in 1:niter) { # iterations
 for (j in 1:nchain) { # chains
   set.seed(i+offset)
   a <- runif(1)
   cat("iter:", i, "chain:", j, "runif:", a, "\n")  }

}

On 6/8/05, Paul Gilbert <pgilbert@bank-banque-canada.ca> wrote:
> Beware that your easy trick will give you the same result every time you
> run it. You need a better scheme if you actually intend to get a new
> experiment each time you run it.
>
> Paul
>
> Gorjanc Gregor wrote:
>
> > Thanks to Duncan, Dimitris as well as James for answers. I'll provide
> > here also example from James, which seems to be the easiest of them
> > all and was not posted to the list:
> >
> > niter <- 3
> > nchain <- 2
> > for (i in 1:niter) { # iterations
> > for (j in 1:nchain) { # chains
> > set.seed(i)
> > a <- runif(1)
> > cat("iter:", i, "chain:", j, "runif:", a, "\n")
> > }
> > }
> >
> > Note that seed is set with iteration counter. This is really neat and
> > simple. I am just wondering if this is OK from "RNG point of view". Can
> > someone comment on that?

> >
> > Lep pozdrav / With regards,
> > Gregor Gorjanc
> >
> > ----------------------------------------------------------------------
> > University of Ljubljana
> > Biotechnical Faculty URI: http://www.bfro.uni-lj.si/MR/ggorjan
> > Zootechnical Department mail: gregor.gorjanc <at> bfro.uni-lj.si
> > Groblje 3 tel: +386 (0)1 72 17 861
> > SI-1230 Domzale fax: +386 (0)1 72 17 888
> > Slovenia, Europe
> > ----------------------------------------------------------------------
> > "One must learn by doing the thing; for though you think you know it,
> > you have no certainty until you try." Sophocles ~ 450 B.C.
> > ----------------------------------------------------------------------
> >
> > -----Original Message-----
> > From: Duncan Murdoch [mailto:murdoch@stats.uwo.ca]
> > Sent: sre 2005-06-08 15:53
> > To: Gorjanc Gregor
> > Cc: r-help@stat.math.ethz.ch
> > Subject: Re: [R] Random seed problem in MCMC coupling of chains
> >
> > On 6/8/2005 9:27 AM, Gorjanc Gregor wrote:
> >
> >>Hello!
> >>
> >>I am performing coupling of chains in MCMC and I need the same value
> >>of seed for two chains. I will show demo of what I want:
> >>
> >>R code, which might show my example is:
> >>niter <- 3
> >>nchain <- 2
> >>tmpSeed <- 123
> >>for (i in 1:niter) { # iterations
> >> for (j in 1:nchain) { # chains
> >> set.seed(tmpSeed)
> >> a <- runif(1)
> >> cat("iter:", i, "chain:", j, "runif:", a, "\n")
> >> tmpSeed <- .Random.seed
> >> }
> >>}
> >>
> >>I get this:
> >>
> >>iter: 1 chain: 1 runif: 0.43588
> >>iter: 1 chain: 2 runif: 0.43588
> >>iter: 2 chain: 1 runif: 0.43588
> >>iter: 2 chain: 2 runif: 0.43588
> >>iter: 3 chain: 1 runif: 0.43588
> >>iter: 3 chain: 2 runif: 0.43588
> >>
> >>but I would like to get:
> >>
> >>iter: 1 chain: 1 runif: 0.43588
> >>iter: 1 chain: 2 runif: 0.43588
> >>iter: 2 chain: 1 runif: 0.67676
> >>iter: 2 chain: 2 runif: 0.67676
> >>iter: 3 chain: 1 runif: 0.12368
> >>iter: 3 chain: 2 runif: 0.12368
> >>
> >>Note that seed value is of course changing, but it is parallel
> >>between chains.
> >
> >
> > set.seed takes an integer, but .Random.seed is a complicated vector.
> > You need to play with .Random.seed directly, and move your setting of
> > tmpSeed out of the inner loop, i.e.
> >
> > > niter <- 3
> > > nchain <- 2
> > > set.seed(123)
> > > tmpSeed <- .Random.seed
> > > for (i in 1:niter) { # iterations
> > + for (j in 1:nchain) { # chains
> > + .Random.seed <- tmpSeed
> > + a <- runif(1)
> > + cat("iter:", i, "chain:", j, "runif:", a, "\n")
> > + }
> > + tmpSeed <- .Random.seed
> > + }
> > iter: 1 chain: 1 runif: 0.2875775
> > iter: 1 chain: 2 runif: 0.2875775
> > iter: 2 chain: 1 runif: 0.7883051
> > iter: 2 chain: 2 runif: 0.7883051
> > iter: 3 chain: 1 runif: 0.4089769
> > iter: 3 chain: 2 runif: 0.4089769
> >
> > However, heed the warnings in ?set.seed: with some generators
> > .Random.seed *does not* contain the full state of the RNG. As far as I
> > know there is no way to obtain the full state.
> >
> > Duncan Murdoch
> >
> > ______________________________________________
> > R-help@stat.math.ethz.ch mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
>



R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Thu Jun 09 03:04:01 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:32:27 EST