I am trying to run a bayesian multivariate generalized linear model with OFT1
, MIS1
, wt
, and g.rate
as the response variables, grid
, collar
, sex
, fieldBirthDate
, and trialnumber
as the predictor variables, and sq_id
as the random effect.
My data looks like this:
sq_id sex grid trialnumber collar OFT1 MIS1 fieldBirthDate wt g.rate
22640 F KL 2 Y -2.1947851 0.08686934 -2.036220 2.04349949 0.202092846
22640 F KL 1 Y 0.7661517 2.65544077 -2.036220 -0.09300674 0.546493570
22641 F KL 1 Y 0.9689955 1.38944543 -2.036220 0.04942701 0.582793646
22653 F SU 1 Y -2.1524967 1.03831742 -1.665209 0.44824150 0.691384500
22657 M KL 1 N 1.0918323 -2.03883227 -1.788879 -0.40636099 0.008592439
22657 M KL 2 N 3.1173521 -2.34449199 -1.788879 2.44231398 0.249185968
And, is available here (this link is no longer active as of 06/22/19). The R package used was MCMCglmm
.
I start by stating my prior:
#inverse whishart prior
prior.miw<-list(R=list(V=diag(4), nu=0.002),
G=list(G1=list(V=diag(4),
nu=0.002,
alpha.mu=c(0,0,0,0),
alpha.V=diag(4)*1000)))
Then my model:
mod.1 <- MCMCglmm(
cbind(OFT1, MIS1, wt, g.rate) ~ (trait-1):grid + (trait-1):collar + (trait-1):sex + (trait-1):fieldBirthDate + (trait-1):trialnumber,
random = ~us(trait):sq_id,
rcov = ~us(trait):units, #allows for trait to have different residuals
family = c("gaussian", "gaussian", "gaussian", "gaussian"), #state response variables distributions
data=multi_data,
prior = prior.miw,
verbose = FALSE,
nitt=103000, #number of iterations
thin=100, #the interval at which the Markov chain is stored
burnin=3000) #number of iterations before samples are stored
But, I then get the following two errors (which change each time I run my code):
G-structure 1 is ill-conditioned: use proper priors if you haven't or rescale data if you have
and
ill-conditioned cross-product: can't form Cholesky factor
I made sure to scale my variables properly (mean centered and standardized) and I have tried other priors such as the inverse gamma prior (that works). However, I want to use the inverse whishart prior. If I modify my inverse whishart prior code by removing the alpha.mu
and alpha.V
functions from it, like so:
prior.miw<-list(R=list(V=diag(4), nu=0.002),
G=list(G1=list(V=diag(4),
nu=0.002)))
my multivariate model runs. But, I would like to keep alpha.mu
and alpha.V
in my prior.
I have two questions:
Why am I getting these errors? (i.e. why is my inverse whishart prior causing me this issue as I've currently written it?)
What is the proper inverse whishart prior with
alpha.mu
andalpha.V
functions in it?
Any suggestions or ideas would be much appreciated!