[R] foreach loop - rejection method

[Steve Weston]

Hi Ben,

It looks like you're getting into trouble because you've
converted the body of a for-loop into a function.
You don't need to do that with foreach, and some of
your problems will be solved by reverting to something
more like your original for-loop, I suspect.

However, it also looks like you need to get away
from assigning to the "f2" vector.  That won't work
in parallel, using %dopar%.  It's possible that something
like this will work:

f2 <- foreach(l=2:(n-1), .combine=c, .packages='MCMCpack',
                       .options.nws=list(chunkSize=250)) %dopar% {
  repeat {
    cand2[l-1] <- rinvgamma(1, phi, lambda[l-1])
    q2 <- dinvgamma(cand2[l-1], phi, lambda[l-1])
    p2 <- cand2[l-1]^-1.5*exp(-y[l]^2/(2*cand2[l-1]))*exp(-((log(cand2[l-1])-mu_t_cand[l-1])^2)/2*sigmasq)
    ratio<-p2/(c[l-1]*q2)
    if (runif(1) < ratio) {break}
  }
  min(p2,c[l-1]*q2)
}

This is pretty complicated code, and I haven't
analyzed it enough to feel at all confident that this
is correct.  But I'm guessing that this is pretty much
like your original for-loop, with the final assignment
removed.

I also recommend trying to get this to work with the
%do% operator first, and then try %dopar%.

- Steve Weston
REvolution Computing


> From: dunno87 <benparker1987 at hotmail.com>
> Date: Wed, Oct 7, 2009 at 5:56 PM
> Subject: [R] foreach loop - rejection method
> To: r-help at r-project.org
>
>
>
> Hi Everybody,
>
> Thanks in advance for your help.
> This is my first time using the foreach statement and I cant get it to work
> properly so here is what i have
>
> test<-function(){
> repeat {
> cand2[l-1]<-rinvgamma(1,phi,lambda[l-1])
> q2<-dinvgamma(cand2[l-1],phi,lambda[l-1])
> p2<-cand2[l-1]^-1.5*exp(-y[l]^2/(2*cand2[l-1]))*exp(-((log(cand2[l-1])-mu_t_cand[l-1])^2)/2*sigmasq)
> ratio<-p2/(c[l-1]*q2)
> if (runif(1)< ratio) {break}
> }
> f2[l-1]<-min(p2,c[l-1]*q2)
> }}
> foreach(l=2:(n-1),.combine=c,.packages='MCMCpack',.options.nws=li
> st(chunkSize=250)) %dopar% test()
>
> even though i've created the array earlier called cand2 when i run this code
> i get the error object "cand2" not found.
>
> What I'm trying to do is the rejection method for my cand2. i need cand2 and
> f2 as vectors as output both of length l-2 in the correct order.
>
> if i do a for loop as normal it works perfectly. Any ideas?
>
> Ben
> --
> View this message in context:
> http://www.nabble.com/foreach-loop---rejection-method-tp25796927p25796927.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.