[R-sig-ME] binary trait correlation across environments (experimental trials) using MCMCglmm?
HADFIELD Jarrod
j@h@dfield @ending from ed@@c@uk
Tue Jan 8 12:15:55 CET 2019
Hi,
Is it possible to share the data - it is hard to troubleshoot otherwise?
Cheers,
Jarrod
On 06/01/2019 19:40, Plough, Louis wrote:
> Hi Jarrod et al.
> Thanks for these suggestions. The model coding the interaction as
> ~us(Trait):animal is working (running), but convergence appears to be
> a way off even after 2 or 3 million iterations, especially for the
> TrialB variance ( survival in the second environment). Would another
> family (ordinal?) be worth considering over the threshold family? Or
> changing the prior somehow, e.g. bumping up the alpha.V=diag(2)*1000?
> We do expect some positive correlation between the two trials...i.e.
> should not be zero.
>
> Coding the model the other way you suggested as animal+animal:Trait
> gives me the following error:
>
>> bi_model_trial_pos <- MCMCglmm(phen ~ 1, random = ~animal+animal:Trial, family = "threshold",prior = prior_b, pedigree = pedigree.t, data = trials_up, nitt = 5e+05, burnin = 25000, thin = 20)
> Error in buildZ(rmodel.terms[r], data = data, nginverse = names(ginverse)) :
>
> interactions not permitted between standard random effects and those
> associated with ginverse
>
> I tried used a similar prior as for the above model, but changed the
> number of structures... (not sure if that is the right prior here):
>
>> prior_b <- list(R = list(V = 1, fix=1), G = list(G1 = list(V = 1, nu = 2),G2=list(V=1,alpha.mu=0, nu=2)))
> For either model, I have been studying your tutorials for how to tweak
> the priors in a sensible way (for a binary trait), but I still
> struggle with how the various supplementary parameters alter the prior
> distribution or CDF... you had mentioned that my previous priors were
> 'very informative', but in which way? Towards zero or too high?
>
> LVP
>
>
> On Thu, Jan 3, 2019 at 4:24 AM HADFIELD Jarrod <j.hadfield using ed.ac.uk> wrote:
>> Hi,
>>
>> The prior is not fixed at one. Also, the priors you are using are quite informative. I would use
>>
>> prior <- list(R = list(V = 1, fix=1), G = list(G1 = list(V = diag(2), nu = 2, alpha.mu=c(0,0), alpha.V=diag(2)*100)))
>>
>> Your equation for the genetic correlation is correct.
>>
>> You shouldn't expect the correlation in family means to equal the model based estimate of the genetic correlation for a number of reasons; most importantly they are on different scales (data versus latent) and the variances of the family means contains residual variation but the covariance doesn't so the correlation in family means is a (downwardly) biased estimator.
>>
>> Cheers,
>>
>> Jarrod
>>
>>
>>
>>
>>
>>
>> On 02/01/2019 22:18, Plough, Louis wrote:
>>
>> Thanks Jarrod for the prompt response!
>>
>> I think the residual variance was set at 1. My prior looks like this:
>>
>> prior <- list(R = list(V = 1, nu = 3), G = list(G1 = list(V = diag(2), nu = 2)))
>>
>> Not sure why the difference for nu in R vs G (3 vs 2)...might have been a typo. Is there guidance on setting the gamma parameter for this kind of binary trait, cross environment correlation model with the 'threshold' family? In the low iteration (toy) run I did for the R-Sig-ME post, the correlation was a bit lower than I would expect based on simply phenotypic correlation of family means, so I think I might need to tweak nu for both?
>>
>> Can I also confirm that I am estimating the rG correctly with the following code:
>>
>> corr.gen<-model_trial4$VCV[,'TrialA:TrialB.animal']/sqrt(model_trial4$VCV[,'TrialA:TrialA.animal']*model_trial4$VCV[,'TrialB:TrialB.animal'])
>>
>> Any reason to use 'TrialA:TrialB.animal' vs 'TrialB:TrialA.animal' in the numerator? Or are they basically equivalent? I wouldn't want add them, would I?
>>
>> Thanks for your help!!!
>>
>> LVP
>>
>>
>>
>> On Wed, Jan 2, 2019 at 3:17 PM HADFIELD Jarrod <j.hadfield using ed.ac.uk> wrote:
>>> Hi,
>>>
>>> Your model bi_model_trial is the correct one. However you must fix the residual variance at one in the prior. You have estimated the genetic (co)variance matrix for the two trials, from which you can obtain the genetic correlation. Alternatively, you could fit animal+animal:Trial which assumes the genetic variances in the two trials are the same and the correlation is positive. The correlation in this latter model is obtained as VAR(animal)/(VAR(animal)+VAR(animal:Trial)). Also, there seems to be a problem with your Surv data as it has three levels rather than 2 and so a cutpoint is being estimated.
>>>
>>> Cheers,
>>>
>>> Jarrod
>>>
>>>
>>>
>>>
>>>
>>>
>>> On 2 Jan 2019, at 16:33, Plough, Louis <lplough using umces.edu> wrote:
>>>
>>> bi_model_trial
>>>
>>>
>>> The University of Edinburgh is a charitable body, registered in Scotland, with registration number SC005336.
--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.
More information about the R-sig-mixed-models
mailing list