[R-sig-ME] cross-sex genetic correlation

Jarrod Hadfield j.hadfield at ed.ac.uk
Thu Jul 27 12:10:44 CEST 2017


Hi Simona,

I'm not sure what the problem is? The h2 have very high uncertainty 
ranging from 0 to 0.4, and as a consequence the genetic correlation 
spans a large range. Accurate estimates of h2 require 1000's of 
individuals from 100's of families, genetic correlations even more.

Cheers,

Jarrod


On 27/07/2017 10:38, Simona Kralj Fiser wrote:
>
> Dear Jarrod and Paul,
>
> Thank you for your reply. We used the suggested prior and model 
> specifications, but we also LOG transformed our weight data (L). The 
> new values are mostly negative. We ran:
>
> prior <- list(R=list(V=diag(2), nu=0.02), G=list(G1=list(V=diag(2), 
> nu=2, alpha.mu <http://alpha.mu/>=c(0,0),alpha.V=diag(2)*1000)))
>
> model14 <- MCMCglmm(L~sex, random=~us(sex):animal, 
> rcov=~idh(sex):units, prior=prior, pedigree=Ped, data=Data1, 
> nitt=100000, burnin=10000, thin=10)
>
> The resulting summary is:
>
> *Iterations = 10001:99991*
> *Thinning interval  = 10*
> *Sample size  = 9000 *
> **
> *DIC: -466.781 *
> **
> *G-structure:  ~us(sex):animal*
> **
> *                 post.mean   l-95% CI u-95% CI eff.samp*
> *sex1:sex1.animal  0.003846  4.515e-10 0.009540     4761*
> *sex2:sex1.animal  0.001122 -6.715e-04 0.003216     1436*
> *sex1:sex2.animal  0.001122 -6.715e-04 0.003216     1436*
> *sex2:sex2.animal  0.002096  1.310e-11 0.004439     5447*
> **
> *R-structure:  ~idh(sex):units*
> **
> *           post.mean l-95% CI u-95% CI eff.samp*
> *sex1.units  0.019094 0.012842 0.025643     5700*
> *sex2.units  0.007019 0.004553 0.009551     6510*
> **
> *Location effects: L ~ sex *
> **
> *            post.mean l-95% CI u-95% CI eff.samp  pMCMC *
> *(Intercept)   -0.9866  -1.0150  -0.9599     9521 <1e-04 ****
> *sex2          -0.2536  -0.2843  -0.2227     9000 <1e-04 ****
> *---*
> *Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1*
>
> With LOG values, the “subscript out of bounds” problem in gone, 
> herit() analyses run and the resulting heritability and correlation 
> values are reasonable. However, the HPD interval is extremely wide.
>
> *> 
> **herit14<-model14$VCV[,'sex1:sex1.animal']/(model14$VCV[,'sex1:sex1.animal']+model14$VCV[,'sex1.units'])*
> *> 
> **herit15<-model14$VCV[,'sex2:sex2.animal']/(model14$VCV[,'sex2:sex2.animal']+model14$VCV[,'sex2.units'])*
> *> **mean(herit14)*
> *[1] 0.1643911*
> *> **mean(herit15)*
> *[1] 0.226494*
> *> **corr.gen <- model14$VCV[, 
> 'sex1:sex2.animal']/sqrt(model14$VCV[,'sex1:sex1.animal']*model14$VCV[,'sex2:sex2.animal'])*
> *> **mean(corr.gen)*
> *[1] 0.4729393*
> *> **HPDinterval(herit14)*
> *            lower     upper*
> *var1 2.149316e-08 0.3883343*
> *attr(,"Probability")*
> *[1] 0.95*
> *> **HPDinterval(herit15)*
> *            lower     upper*
> *var1 1.539724e-09 0.4509762*
> *attr(,"Probability")*
> *[1] 0.95*
> *> **HPDinterval(corr.gen)*
> *          lower    upper*
> *var1 -0.1849416 0.999439*
> *attr(,"Probability")*
> *[1] 0.95*
>
> We are starting to run out of ideas on why this is happening or where 
> the problem lies. We’d appreciate any further advice!
>
> Eva and Simona
>
>
> On 26 July 2017 at 14:42, Jarrod Hadfield <j.hadfield at ed.ac.uk 
> <mailto:j.hadfield at ed.ac.uk>> wrote:
>
>     Hi,
>
>     The second way is a *much* better way of doing it but should give
>     the same answer. However, in both cases the residual covariance is
>     not identifiable (no individual is both male and female) and so
>     you should use idh rather than us.
>
>     The "subscript out of bounds" error is to do with your code that
>     post-processes the model output not an issue with MCMCglmm.
>     Probably you have used the wrong names for the (co)variance
>     components.
>
>     Also, you haven't passed the prior to MCMCglmm, nor is the prior a
>     valid one for the problem as it specifies scalar variances rather
>     than 2x2 covariance matrices. You could try
>
>     prior2 <- list(R=list(V=diag(2), nu=0.02),
>     G=list(G1=list(V=diag(2), nu=2, alpha.mu
>     <http://alpha.mu>=c(0,0),alpha.V=diag(2)*1000)))
>
>     Cheers,
>
>     Jarrod
>
>
>
>
>     On 26/07/2017 13:33, Simona Kralj Fiser wrote:
>
>         model <- MCMCglmm(W~sex, random=~us(sex):animal,
>         rcov=~us(sex):units,
>         prior=prior2, pedigree=Ped, data=Data1, nitt=100000,
>         burnin=10000, thin=10)
>
>
>
>     -- 
>     The University of Edinburgh is a charitable body, registered in
>     Scotland, with registration number SC005336.
>
>
>
>
> -- 
> Doc. Dr. Simona Kralj-Fišer
>
> Institute of Biology
> Scientific Research Centre of the Slovenian Academy of Sciences and Arts
> Novi trg 2, P. O. Box 306, SI-1001 Ljubljana.
> Phone: ++38614706333; fax: ++38614257797
>
> http://ezlab.zrc-sazu.si/
>
> http://bijh.zrc-sazu.si/sl/sodelavci/simona-kralj-fi%C5%A1er-sl#v

-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: not available
URL: <https://stat.ethz.ch/pipermail/r-sig-mixed-models/attachments/20170727/87c10e5d/attachment.ksh>


More information about the R-sig-mixed-models mailing list