[R-sig-ME] Error in variance estimation using bivariate animal model

[Stephane Chantepie]

Dear all,

I am using MCMCglmm function to construct bivariate animal models of bustard 
sperm production according to age. Until six years old, each age are 
considered as a trait. After 6 years old, age-class "7-8" and "9-11" are 
computed to have enought information.

In a first part, additive genetic variance have been estimated by using simple 
animal model. Results are not sensitive to prior and give quite realistic 
estimations

In a second part, I have tried to estimate covariance between these ages. And 
I have some problems in estimating variance-covariance for the bivariate 
models which use the age class "7-8" or "9-11". For example, when I compute 
the model cbind(age_1,age_9_11), the estimation of variance for age_1 and 
age_9_11 are similar to estimations of the  simple animal models. The 
estimation of covariance also appears consistent with the others covariance 
results, so it pretty good. But when I do cbind(age_2,age_9_11), the variance 
estimation of age_9_11 become really hight (out of the range I have found in 
the simple model). So the covariance estimation appears really hight too. I do 
not really understand why it works when I use age_1 against my age-class and 
not when I use the age 2 3 4 5 6. The chains seem to converge and there is no 
autocorrelation 
All variance estimation given by bivariate model wich do not use age-class are 
consistent with the simple models.

result for cbind(age_1,age_9_11) : http://ubuntuone.com/2xewilItvXmqvZZhMW3sTN
result for cbind(age_2,age_9_11) : http://ubuntuone.com/7hVy0DZ2z0qNHkfiZUZXis

I run five times the same model (as below) and concatenate results (same  
prior,same burning, same thin and random seed)
I have just added a ID random parameter when I use a age-class as one of 
bivariate trait (because I have repeated measures)

** spz_var=as.numeric(cbind(var(na.omit(age_2)/2),var(na.omit(age_7_8)/2)))
** prior<-list(R=list(V=diag(2)*spz_var,nu=2), 		
G=list(G1=list(V=diag(2)*spz_var, nu=2),
G2=list(V=diag(2)*spz_var, nu=2),
G3=list(V=diag(2)*spz_var, nu=2)))
** spz<-MCMCglmm(cbind(age_2,age_7_8)~ trait-1, random=~us(trait):animal + 
us(trait):birth + us(trait):ID, rcov=~us(trait):units,nitt=1300000, thin=1000, 
burnin=300000, prior=prior, verbose=TRUE, pedigree=ped, 
family=c("gaussian","gaussian"), data=dat)

Should I modify  "us(trait):ID" to fit random parameter only on the second 
trait? How can I do that? Is my results arenormal?

If you have any idea

Many thanks for your reply

--
Stephane Chantepie
CNRS Phd  candidate
Muséum national d'Histoire naturelle
55 rue Buffon
75005 paris
E-mail : chantepie_at_mnhn_dot_fr
--