[R-sig-ME] Priors for and estimating heritability from an ordinal and a Poisson animal model

[Pierre B. de Villemereuil]

(I thought I responded to this e-mail a while ago... but my response 
might have not leave my computer !)

Concerning the priors, two things :
- For the first model, it is not very clear, but I assume you use 
"deformed toes" as a binary variable. In that case, a better prior would 
be the chi2 distribution. You can have it by using the following prior :

prior2=list(R=list(V=1, fix=1), G=list(G1=list(V=1, nu=1000, alpha.mu=0,alpha.V=1)))

This prior is a bit skewed toward little heritabilities, but the one you 
used is actually skewed toward strong heritabilities.
- For the second model, what is "p.var" ? If it is the phenotypic 
variance calculated from the data, you should not use it. It is 
ill-advised to use the data under investigation to set up the prior. I 
think you try there to use your first prior or the following one :

prior1.1<-list(G=list(G1=list(V=1,n=0.002)),R=list(V=1,n=0.002))

I hope this is answering some of your questions. Don't hesitate to check 
different priors to see if your posterior distribution is influenced or 
not by the prior.

Also, normally, the heritability being an intra-class coefficient, it 
should be independant from the fact that the residual variance is fixed. 
You should check the MCMCglmm Course Notes from Jarrod Hadfield 
concerning this point.

Finally, to add the random factor, using the chi2 distribution, you 
should just use :

prior2=list(R=list(V=1, fix=1), G=list(G1=list(V=1, nu=1000, alpha.mu=0,alpha.V=1),
G2=list(V=1, nu=1000, alpha.mu=0,alpha.V=1)))

This would not alter much the shape of the prior on the heritability.

Cheers,
Pierre.


Le 06/09/2012 17:13, Helen Ward a écrit :
> Hello everybody,
>
> I've just come back to some models I was working on a while ago and was
> hoping someone could help me clear up a couple of issues. My goal is to
> calculate the heritability of two independent traits: number of deformed
> toes, and age of 1st reproduction.
>
> Number of deformed toes is a data set of integers ranging from 0 to 6:
> it looks most like a zero-inflated Poisson distribution.
> Age of 1st reproduction is a data set of integers ranging from 0
> upwards: it has a clear Poisson distribution.
>
> I have built animal models for these traits in ASReml specifying that
> both data sets are Poisson distributed. However, I cannot compare my
> heritability estimates with those estimated for other Gaussian traits as
> in a Poisson model in ASReml the residual variance is fixed at 1.
>
> I am now trying to model these traits in MCMCglmm, deformed toes in an
> ordinal model and age of 1st repro in a Poisson model.
>
> At the moment for deformed toes I am using
>
> prior2=list(R=list(V=1, fix=1), G=list(G1=list(V=1, nu=1,
> alpha.mu=0,alpha.V=100)))
>
> model1<-MCMCglmm(Toes~1,random=~animal,family="ordinal",pedigree=ped,data=data,prior=prior2,nitt=500000,thin=300,burnin=300000,verbose=FALSE,pl=TRUE)
> then
>
> posterior.heritability1<-model1$VCV[,"animal"]/(model1$VCV[,"animal"]+model1$VCV[,"units"]+1)to
> calculate heritabilty.
>
> For age of 1st repro I am using
>
> prior1.1<-list(G=list(G1=list(V=matrix(p.var/2),n=1)),
> R=list(V=matrix(p.var/2),n=1))
>
> model1.1<-MCMCglmm(Poisson1strepro~1,random=~animal,family="poisson",pedigree=ped,data=data,prior=prior1.1,nitt=100000,thin=75,burnin=25000,verbose=FALSE)
>
> then
>
> posterior.heritability1.1<-model1.1$VCV[,"animal"]/(model1.1$VCV[,"animal"]+model1.1$VCV[,"units"])
> to calculate heritability
>
>
> Please can someone comment on whether they think these priors and
> methods of estimating heritability are OK?
>
>
> Please can someone also tell me whether the heritability estimate from
> the ordinal model means the same thing as a heritability estimated from
> a model in which the residual variance is not fixed?
>
>
> Finally, I would like to add a random factor (Year of birth) to both
> models. Please can someone suggest how I might alter the prior in the
> ordinal toes model to do so?
>
>
> Many (many) thanks,
> Helen
>
>
> 	[[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models