[R-sig-ME] Meta-analysis for heritability using MCMCglmm?

Sat Feb 14 09:34:41 CET 2015

Hi Jackie,

1/ For the residuals you can do:

my_data$y-predict(my_model)

In the code above the predictions are marginal with respect to the  
random effects. If you want the predictions to include the random  
effects then use:

my_data$y-predict(my_model, marginal=NULL)

2/ I don't have a good solution for the weighted analysis. In fact we  
didn't do a weighted analysis in the end because only about 50% of  
studies provided standard errors. We figured that we lost more power  
by discarding the 50% than we would gain by incorporating information  
about the relative precision of different estimates. REML standard  
errors are based on the sampling distribution of estimates tending to  
the normal as n goes to infinity. In reality the sampling distribution  
will be skewed (particularly when h2 and/or sample size is low) like  
the posterior distribution. Its far from perfect but you could use the  
standard error given by:

0.5*((u95-pm)/qnorm(0.975)+(l95-pm)/qnorm(0.025))

where l95 and u95 are the lower and upper 95% credible intervals and  
pm is the posterior mean/mode.  The motivation behind this is *if* the  
posterior was normal and the asymmetry in the 95% credible intervals  
was just due to Monte Carlo error then this would be the best estimate  
of the posterior standard deviation (sort of like a standard error)  
given the information at hand. Not great justification!

3/ DIC (as focussed in MCMCglmm) is not a reliable model selection  
criterion for scientific inference. For Gaussian resopnses you could  
refocus DIC at the correct level post-analysis and then it would be  
more reliable. I should probably write a function to do this ....

Cheers,

Jarrod

Quoting Jackie Wood <jackiewood7 at gmail.com> on Sun, 1 Feb 2015 12:10:51 -0500:

> Hi Jarrod,
>
> We were finally able to dig into the statistical analysis of our
> heritability data (as a reminder, we are conducting a meta-analysis
> investigating heritability in relation to population size) and of course a
> few questions have come up.
>
> You had mentioned examining the residuals for the model. This may seem like
> a "beginner" question, how does one extract the residuals from MCMCglmm?
> The residuals.mcmcglmm function does not work. Summary plots of the MCMC
> parameter estimates appear to be roughly normally distributed, and the
> traces seem fine. Is this what you were referring to? We ran an unweighted
> analysis (heritability estimates without SEs) in both MCMCglmm and lmer to
> see if they give similar results, and they are basically concordant as
> well. However, the residual distribution for the models run in lmer
> (treating heritability as gaussian) is slightly skewed. These are from an
> unweighted analysis, though, as I have also read that lme4 is unsuitable
> for conducting formal weighted meta-analyses.
>
> We also had a few questions that we thought would be worth discussing about
> some methodological issues relating to incorporating common estimates of
> heritability in the literature. Bayesian methodologies have become
> increasingly popular to use when estimating trait heritabilities, but
> bayesian estimates do not provide typical standard error or variance
> estimates, as parent-offspring/ANOVA/REML methods do. Published bayesian
> heritability estimates typically only include asymmetric confidence
> intervals, and we unsure whether these can be translated into variance
> estimates that can be used to weight our meta-analysis. For now, we plan on
> performing a weighted meta-analysis using heritability estimates that
> provide S.E.s, and an additional unweighted analysis that will include the
> bayesian point-estimates we have collected from the literature. We were
> wondering if you resolved this issue in your own heritability meta-analysis
> and knew of a way to incorporate bayesian estimates (which form a
> considerable proportion of the suitable heritability estimates available,
> at least in recent history) into a formal weighted meta-analysis.
>
> Additionally, we were wondering about the suitability of DIC to conduct
> model selection for our analysis of heritability. I recall reading on this
> SIG list that you had mentioned that there were potential issues using DIC
> for hierarchical models as well as non-gaussian data. Since we're treating
> heritability as gaussian, would it still be appropriate?
>
> Any advise would be much appreciated!
>
>
> On Thu, Jan 15, 2015 at 3:41 PM, Jackie Wood <jackiewood7 at gmail.com> wrote:
>
>> Hi Jarrod and Ken,
>>
>> Hope you had a great New Year! Thanks so much for your responses to my
>> inquiry. Given that we've been using MCMCglmm all along, we'll probably
>> stick with it unless there's a compelling reason to change programs. We'll
>> be running the h2 models in the coming days and will specify a Gaussian
>> distribution as Jarrod suggested; we have quite a bit of data so hopefully
>> the residuals will behave!
>>
>> The advice is much appreciated as always!
>> Jackie
>>
>> On Fri, Dec 26, 2014 at 1:58 AM, Jarrod Hadfield <j.hadfield at ed.ac.uk>
>> wrote:
>>
>>> Hi Jackie,
>>>
>>> The data are not binomial they are continuous: a beta distribution is
>>> probably most appropriate for continuos observations bounded by 0 and 1.
>>> However, although heritabilities are bounded by 0 and 1, heritability
>>> estimates are not necessarily so, depending on the method of inference (for
>>> example it would be possible to get a negative parent-offspring regression,
>>> either by chance or through certain types of maternal effect).
>>>
>>> We have just finished a meta-analysis of h2 estimates and just treated
>>> them as Gaussian. The distribution of the residuals wasn't far off and I
>>> think the conclusions are robust to the distributional assumptions. Have
>>> you checked your residuals - do they look badly non-normal?
>>>
>>>
>>> Cheers,
>>>
>>> Jarrod
>>>
>>>
>>>
>>>
>>>
>>> Quoting Ken Beath <ken.beath at mq.edu.au> on Wed, 24 Dec 2014 12:30:03
>>> +1100:
>>>
>>>  If you have the original data giving the numerator and denominator for
>>>> the
>>>> proportion then it is binomial data, and can be modelled in a
>>>> met-analysis.
>>>> I don't know if this can be done with MCMCglmm but should be possible
>>>> with
>>>> STAN, JAGS or BUGS. All will require a bit of effort in setting up the
>>>> model.
>>>>
>>>> On 24 December 2014 at 07:17, Jackie Wood <jackiewood7 at gmail.com> wrote:
>>>>
>>>>  Dear R-users,
>>>>>
>>>>> I am attempting to conduct a meta-analysis to investigate the
>>>>> relationship
>>>>> of narrow-sense heritability with population size. In previous work, I
>>>>> have
>>>>> used MCMCglmm to conduct a formal meta-analysis which allowed me to
>>>>> account
>>>>> for the effect of sampling error through the argument "mev". This was
>>>>> relatively easy to do for a continuous response variable, however,
>>>>> heritability is presented as a proportion and is therefore bounded by 0
>>>>> and
>>>>> 1 which clearly changes the situation.
>>>>>
>>>>> In fact, I am not actually certain if it possible to conduct a formal
>>>>> weighted meta-analysis on the heritability data using MCMCglmm. I have
>>>>> seen
>>>>> elsewhere where data presented as a proportion (survival,
>>>>> yolk-conversion
>>>>> efficiency for example) has been logit transformed and fitted using a
>>>>> Gaussian error distribution (though this was done using REML rather than
>>>>> Bayesian modelling) but I don't know if this is a legitimate strategy
>>>>> for a
>>>>> formal meta-analysis using heritability as a response variable since any
>>>>> transformation applied to the heritability data would also need to be
>>>>> applied to the standard errors?
>>>>>
>>>>> I would greatly appreciate any advice on this matter!
>>>>>
>>>>> Cheers,
>>>>> Jackie
>>>>>
>>>>> --
>>>>> Jacquelyn L.A. Wood, PhD.
>>>>> Biology Department
>>>>> Concordia University
>>>>> 7141 Sherbrooke St. West
>>>>> Montreal, QC
>>>>> H4B 1R6
>>>>> Phone: (514) 293-7255
>>>>>
>>>>>         [[alternative HTML version deleted]]
>>>>>
>>>>> _______________________________________________
>>>>> R-sig-mixed-models at r-project.org mailing list
>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>
>>>>>
>>>>
>>>>
>>>> --
>>>>
>>>> *Ken Beath*
>>>> Lecturer
>>>> Statistics Department
>>>> MACQUARIE UNIVERSITY NSW 2109, Australia
>>>>
>>>> Phone: +61 (0)2 9850 8516
>>>>
>>>> Building E4A, room 526
>>>> http://stat.mq.edu.au/our_staff/staff_-_alphabetical/staff/beath,_ken/
>>>>
>>>> CRICOS Provider No 00002J
>>>> This message is intended for the addressee named and may...{{dropped:9}}
>>>>
>>>> _______________________________________________
>>>> R-sig-mixed-models at r-project.org mailing list
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>
>>>>
>>>>
>>>
>>>
>>> --
>>> The University of Edinburgh is a charitable body, registered in
>>> Scotland, with registration number SC005336.
>>>
>>>
>>>
>>
>>
>> --
>> Jacquelyn L.A. Wood, PhD.
>> Biology Department
>> Concordia University
>> 7141 Sherbrooke St. West
>> Montreal, QC
>> H4B 1R6
>> Phone: (514) 293-7255
>>
>>
>
>
> --
> Jacquelyn L.A. Wood, PhD.
> Biology Department
> Concordia University
> 7141 Sherbrooke St. West
> Montreal, QC
> H4B 1R6
> Phone: (514) 293-7255
>

-- 
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.