# NOT RUN {
# Create AR(0)-SV model for population growth rate
r <- ts(ew[2:167]/ew[1:166]-1, start=1841)
y <- diff(r)
plot(y, type = "l")
rv0 <- rv.bugs(y)
print(rv0)
# AR(2)-RV model with alternative priors
rv2 <- rv.bugs(y, rv.eps.prior = "dbeta(1,20)")
print(rv2)
# AR(0)-RV model with posterior simulations
rv0 <- rv.bugs(y, sim = TRUE)
print(rv0)
# }
# NOT RUN {
# Run in OpenBUGS
writeLines(rv0$bug, "rv0.txt")
library("R2OpenBUGS")
# Run model (can take some time, depending on data length)
rv0.bug <- bugs(data = rv0$data,
inits = list(inits(rv0)),
param = c(nodes(rv0, "prior")$name, "y.sim", "y.new"),
model = "rv0.txt",
n.iter = 11000, n.burnin = 1000, n.chains = 1)
# Plot the parameters posteriors and traces
library("coda")
param.mcmc <- as.mcmc(rv0.bug$sims.matrix[, nodes(rv0, "prior")$name])
plot(param.mcmc)
# Plot posterior simulations using fanplot
library("fanplot")
y.mcmc <- rv0.bug$sims.list$y.sim
y.pn <- pn(y.mcmc, st = start(y)[1])
plot(y, type = "n")
fan(y.pn)
lines(y)
# Plot volatility
h.mcmc <- rv0.bug$sims.list$h
h.pn <- pn(h.mcmc, st = start(y)[1])
sigma.pn <- pn(sims = sqrt(exp(h.mcmc)), st = start(y)[1])
par(mfrow = c(2, 1), mar = rep(2, 4))
plot(NULL, type = "n", xlim = tsp(h.pn)[1:2], ylim = range(h.pn[, 5:95]), main = "h_t")
fan(h.pn)
plot(NULL, type = "n", xlim = tsp(sigma.pn)[1:2], ylim = range(sigma.pn[, 1:95]), main = "sigma_t")
fan(sigma.pn)
# }
Run the code above in your browser using DataLab