Dear Stephane,
When you say crash do you mean crash in the sense of a segfault or
in the sense that it stops with the errors:
-Mixed model equations singular: use a (stronger) prior
-ill-conditioned G/R structure: use proper priors if you haven't or
rescale data if you have
If the latter, it may just require a rescaling of your continuous
covariates by using something like scale(). If the former, it would
be good for me to have a reproducible example as it means there is a
bug.
Cheers,
Jarrod
uoting chantepie at mnhn.fr on Fri, 09 Mar 2012 10:33:42 +0100:
Dear all,
I am using MCMCglmm function to construct bivariate animal models
of bustard sperm production according to age-classes.
The problem is that the models can stochastically crash before the
end of the run (at 2000 iterations or 120000 or other) or can
finish. For the model which does not finish, R returns different
errors as:
-Mixed model equations singular: use a (stronger) prior
-ill-conditioned G/R structure: use proper priors if you haven't or
rescale data if you have
For the models which reach the end, the estimations of genetic
additive variance appear quite good (nice bell shaped posterior
disctribution).
The problem still remains when I remove the animal term.
When I run univariate models, it works fine and the posteriors
distributions look very good.
Strangely, the more data I have, the more models crash (the largest
amount of data I have is 65000 data for 2400 individuals for one
model).
The model looks like:
priorExp<-list(G=list(G1=list(V=diag(2), nu=2,
alpha.mu=rep(0,2),alpha.V=diag(2)*100000),
G2=list(V=diag(2), nu=2, alpha.mu=rep(0,2),alpha.V=diag(2)*100000),
G3=list(V=diag(2), nu=2, alpha.mu=rep(0,2),alpha.V=diag(2)*100000)),
R=list(V=diag(2), nu=2))
spz<-MCMCglmm(cbind(age1_2,age5_6)~trait-1 + trait:tse+
trait:I(tse^2)+ trait:joe + trait:I(joe^2),
random=~us(trait):animal+us(trait):ID+us(trait):annee ,
rcov=~us(trait):units,nitt=150000, thin=1000, burnin=100000,
prior=priorExp, verbose=TRUE, pedigree=ped,
family=c("gaussian","gaussian"), data=dat)
For the fixed effects : I use 4 continuous parameters as correction
for each trait
For the random effects: I use, individuals, years and animal parameters
I have also tried more informative prior (as described in WAMWIKI)
but the problem was the same.
To give you an example :
Because of computing limitation, I use multi-chain process. I run
several times the same model (as above) and concatenate results
(same prior,same burning, same thin and random seed) to obtain at
least 1000 estimates (50 estimates by model). In this context, I
ran 50 bivariable models with the age-class age1_2 and the
age-class age5_6 but only 9 models of the 50 models reached the end.
When we look fixed parameters estimates (estimate are binded for
the nine models : http://ubuntuone.com/3Gi8GwjcRk3P01MxJp2qLe ), we
can see that the estimates are really close to 0. Could it be the
problem?
When we look ramdom parameters estimates (estimate are binded for
the nine models : http://ubuntuone.com/42oaP9euG1m2LNipMawcHX ),
the residual estimates do not look very good. Could it be the
problem?
Last thing, if I try to add a cubic effect, all the models crash
(same error than before or memory mapped error).
I really do not know where the problem comes from. Do you have an idea?
Thanks
--
Stephane Chantepie
CNRS Phd candidate
Mus?um national d'Histoire naturelle
55 rue Buffon
75005 paris
E-mail : chantepie_at_mnhn.fr
--