if (FALSE) {
## We conduct CV on the classic "dyes" BUGS model.
dyesCode <- nimbleCode({
for (i in 1:BATCHES) {
for (j in 1:SAMPLES) {
y[i,j] ~ dnorm(mu[i], tau.within);
}
mu[i] ~ dnorm(theta, tau.between);
}
theta ~ dnorm(0.0, 1.0E-10);
tau.within ~ dgamma(0.001, 0.001); sigma2.within <- 1/tau.within;
tau.between ~ dgamma(0.001, 0.001); sigma2.between <- 1/tau.between;
})
dyesData <- list(y = matrix(c(1545, 1540, 1595, 1445, 1595,
1520, 1440, 1555, 1550, 1440,
1630, 1455, 1440, 1490, 1605,
1595, 1515, 1450, 1520, 1560,
1510, 1465, 1635, 1480, 1580,
1495, 1560, 1545, 1625, 1445),
nrow = 6, ncol = 5))
dyesConsts <- list(BATCHES = 6,
SAMPLES = 5)
dyesInits <- list(theta = 1500, tau.within = 1, tau.between = 1)
dyesModel <- nimbleModel(code = dyesCode,
constants = dyesConsts,
data = dyesData,
inits = dyesInits)
# Define the fold function.
# This function defines the data to leave out for the i'th fold
# as the i'th row of the data matrix y. This implies we will have
# 6 folds.
dyesFoldFunction <- function(i){
foldNodes_i <- paste0('y[', i, ', ]') # will return 'y[1,]' for i = 1 e.g.
return(foldNodes_i)
}
# We define our own loss function as well.
# The function below will compute the root mean squared error.
RMSElossFunction <- function(simulatedDataValues, actualDataValues){
dataLength <- length(simulatedDataValues) # simulatedDataValues is a vector
SSE <- 0
for(i in 1:dataLength){
SSE <- SSE + (simulatedDataValues[i] - actualDataValues[i])^2
}
MSE <- SSE / dataLength
RMSE <- sqrt(MSE)
return(RMSE)
}
dyesMCMCconfiguration <- configureMCMC(dyesModel)
crossValOutput <- runCrossValidate(MCMCconfiguration = dyesMCMCconfiguration,
k = 6,
foldFunction = dyesFoldFunction,
lossFunction = RMSElossFunction,
MCMCcontrol = list(niter = 5000,
nburnin = 500))
}
Run the code above in your browser using DataLab