From: ekwaters <ekwaters_at_unimelb.edu.au>

Date: Mon, 02 Mar 2009 21:39:27 -0800 (PST)

Date: Mon, 02 Mar 2009 21:39:27 -0800 (PST)

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:

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)*((b12)^170)*((1-b12)^250)*((b13)^217)*((1-b13)^38))>(2*runif(1,min=0,max=1)*.5))
{p11[i]=u; z=1}

}

}

Ned

-- View this message in context: http://www.nabble.com/R---need-more-memory%2C-or-rejection-sampling-algorithm-doesn%27t-work--tp22302800p22302800.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.Received on Tue 03 Mar 2009 - 04:42:38 GMT

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.2.0, at Tue 03 Mar 2009 - 06:30:23 GMT.

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help.
Please read the posting
guide before posting to the list.
*