MCMC gradually slows down
-----Original Message-----
From: r-help-bounces at r-project.org
[mailto:r-help-bounces at r-project.org] On Behalf Of Jens Malmros
Sent: Sunday, November 08, 2009 11:11 AM
To: r-help at r-project.org
Subject: [R] MCMC gradually slows down
Hello,
I have written a simple Metropolis-Hastings MCMC algorithm for a
binomial parameter:
MHastings = function(n,p0,d){
theta = c()
Change that theta <- c() to
theta <- numeric(n)
and it will go faster. Growing datasets
can result in a lot of unnecessary memory
management. The original had an obvious
quadratic quality in the plot of n vs. MHastings(n,.2,2)
and preallocating theta to its final length
made it look linear up to n=100000 and at n=100000
the time for the original was 35 seconds vs 3.75
for the preallocated version.
Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com
theta[1] = p0
t =1
while(t<=n){
phi = log(theta[t]/(1-theta[t]))
phisim = phi + rnorm(1,0,d)
thetasim = exp(phisim)/(1+exp(phisim))
r =
(thetasim)^4*(1-thetasim)^8/(theta[t]^4*(1-theta[t])^8)
if(runif(1,0,1)<r){
theta[t+1] = thetasim
} else {
theta[t+1] = theta[t]
}
t = t+1
if(t%%1000==0) print(t) # diagnostic
}
data.frame(theta)
}
The problem is that it gradually slows down. It is very fast in the
beginning, but slows down and gets very slow as you reach about 50000
iterations and I need do to plenty more.
I know there are more fancy MCMC routines available, but I am really
just interested in this to work.
Thank you for your help,
Jens Malmros
______________________________________________ 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.