[R] R - need more memory, or rejection sampling algorithm doesn't work?

Daniel Nordlund djnordlund at verizon.net
Tue Mar 3 07:21:10 CET 2009


> -----Original Message-----
> From: r-help-bounces at r-project.org 
> [mailto:r-help-bounces at r-project.org] On Behalf Of ekwaters
> Sent: Monday, March 02, 2009 9:39 PM
> To: r-help at r-project.org
> Subject: [R] R - need more memory, or rejection sampling 
> algorithm doesn't work?
> 
> 
> Hi all,
> 
> I am trying to run rejection sampling for the quantity z11 in 
> the function
> below.  Unfortunately I can't simplify the function further 
> so that z11 only
> appears once. 
> 
> Whenever I run the algorithm, R looks as if it is running it (no error
> messages or anything), but then nothing happens for minutes...how long
> should it take to run something like this in R? I have tried 
> in in both
> linux and windows.
> 
> I realise standard rejection sampling is pretty cumbersome, 
> so I have also
> tried to run a Gibbs sampler doing rejection sampling instead 
> of the form of
> the algorithm here, with the same result. R thinks about 
> starting to run it,
> then freezes.
> 
> Is this a memory issue, or an issue with the algorithm?
> 
> > count=0
> > k=1
> > f=matrix(NA,1000)
> 
> >while(k<1001){
> 
> z11=runif(1,min=0,max=2)
> 
> r11=(((log(z11/1-z11)^231)*((1-log(z11/1-z11))^62)*((b12)^170)
> *((1-b12)^250)*((b13)^217)*((1-b13)^^38))/2*.5)
> 
> if(r11<runif(1,min=0,max=1))
> 
> {f(k)=z11; k=k+1}
> count=count+1
> }
> 
> THe GIbbs sampler I tried looks as follows:
> 
> > p11=matrix(0.5,2000)
> 
> for(i in 2000)
> {
> z=0
> while(z==0)
> {
> u=runif(1,min=0,max=2)
> if(
> ((log(p11[i]/1-p11[i])^231)*((1-log(p11[i]/1-p11[i]))^62)*((b1
> 2)^170)*((1-b12)^250)*((b13)^217)*((1-b13)^38))>(2*runif(1,min
> =0,max=1)*.5))
> {p11[i]=u; z=1}
> }
> }
> 
> 
> Ned

Ned,

Your first rejection sampling algorithm will not run as presented.  There
are syntax errors and no indication of where variables b12 and b13 get their
values from, so it is not possible to say for sure what the problem is.
However, it is not a memory problem (unless you have a broken system).
There are also problems in your second routine.  We need reproducible code
if you want R-help to do more than guess.

Dan 

Daniel Nordlund
Bothell, WA USA




More information about the R-help mailing list