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

Simona Kralj Fiser simonakf at gmail.com
Thu Jul 27 13:06:36 CEST 2017


I see :-(

Thanks a lot for your help.

Simona

On 27 July 2017 at 12:10, Jarrod Hadfield <j.hadfield at ed.ac.uk> wrote:

> 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=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> 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=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 <(01)%20470%2063%2033>; fax: ++38614257797
> <(01)%20425%2077%2097>
>
> http://ezlab.zrc-sazu.si/
>
> http://bijh.zrc-sazu.si/sl/sodelavci/simona-kralj-fi%C5%A1er-sl#v
>
>
>
> 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

	[[alternative HTML version deleted]]



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