if (FALSE) {
data(rhizoctonia)
### Define the model
corrf <- "spherical"
kappa <- 0
ssqdf <- 1
ssqsc <- 1
betm0 <- 0
betQ0 <- .01
family <- "binomial.probit"
### Skeleton points
philist <- c(100, 140, 180)
omglist <- c(.5, 1)
parlist <- expand.grid(linkp=0, phi=philist, omg=omglist, kappa = kappa)
### MCMC sizes
Nout <- 100
Nthin <- 1
Nbi <- 0
### Take MCMC samples
runs <- list()
for (i in 1:NROW(parlist)) {
runs[[i]] <- mcsglmm(Infected ~ 1, family, rhizoctonia, weights = Total,
atsample = ~ Xcoord + Ycoord,
Nout = Nout, Nthin = Nthin, Nbi = Nbi,
betm0 = betm0, betQ0 = betQ0,
ssqdf = ssqdf, ssqsc = ssqsc,
phi = parlist$phi[i], omg = parlist$omg[i],
linkp = parlist$linkp[i], kappa = parlist$kappa[i],
corrfcn = corrf,
corrtuning=list(phi = 0, omg = 0, kappa = 0))
}
bf <- bf1skel(runs)
bfall <- bf2new(bf, phi = seq(100, 200, 10), omg = seq(0, 2, .2))
plotbf2(bfall, c("phi", "omg"))
plotbf2(bfall, c("phi", "omg"), profile = TRUE, type = "b", ylab="log(BF)")
}
Run the code above in your browser using DataLab