computePredictedValues - Error in Eta[[r]][as.character(dfPiNew[, r]), ] : no 'dimnames' attribute for array

[wilkesma]

I fitted several models and saved them externally through a nohup call. But when I read them back in and try to computePredictedValues I get one of two errors depending on the random effects structure.

With just one random effect - "Error in Eta[[r]][as.character(dfPiNew[, r]), ] : no 'dimnames' attribute for array".

With nested random effects - "Error in PiNew[, r] <- sapply(dfPiNew[, r], function(s) which(rowNames == : incorrect number of subscripts on matrix"

In testing, when I fitted model and computePredictedValues in the same R session, there were no issues. But this is an impractical and risky approach for running the full models, and I would like to be able to work flexibly after fitting the model.

See below for an example of how I am specifying the models. I tried defining the random effects as both levels() and unique() but I get the same issue.

studyDesign <- data.frame(season=meta$season, site=as.character(meta$SITE_ID))
season <- HmscRandomLevel(units=levels(studyDesign$season))
site <- HmscRandomLevel(units=levels(studyDesign$site))
X.data <- data.frame(intercept=1, time=meta$days)

mod <- Hmsc(Y=taxa, XData=X.data, XFormula=~time, studyDesign=as.matrix(studyDesign),
              TrData=traits, TrFormula=TrFormula, phyloTree=tree,
              ranLevels=list("site"=site, "season"=season), distr=link)
mod <- sampleMcmc(mod, samples = samples, thin = thin, transient = transient,
                         nChains = nChains, nParallel = nChains)

Any advice would be greatly appreciated.

[wilkesma]

I finally got to the bottom of this. The solution had to do with the the class of variables in the study design and random levels. The fact that I was passing the study design to the Hmsc function as a matrix was then also an issue. I notice there are a few corners of the internet where people encountered this problem, so the following script might help

'''
X.data <- data.frame(intercept=1, time=meta$days)
studyDesign <- data.frame(season=meta$season, site=as.character(meta$SITE_ID), stringsAsFactors=T)
rL1 <- HmscRandomLevel(units=levels(studyDesign$season))
rL2 <- HmscRandomLevel(units=unique(studyDesign$site))
mod <- Hmsc(Y=taxa, XData=X.data, XFormula=~time, studyDesign=studyDesign,
TrData=traits, TrFormula=TrFormula, phyloTree=tree,
ranLevels=list("site"=rL2, "season"=rL1), distr=link)
mod <- sampleMcmc(mod, samples = samples, thin = thin, transient = transient,
nChains = nChains, nParallel = nChains)
'''

Defining SITE_ID as a string, setting stringsAsFactors=T, and using unique() for rL2 was the key!