[R-sig-ME] prior specification in MCMCglmm
Ned Dochtermann
ned.dochtermann at gmail.com
Wed Apr 27 00:35:56 CEST 2011
There's a mistake in my email below.
The generalization of an inverse-gamma prior for multi-response models
should, I think, just be:
list(V=diag(x), nu=1.002)
Sorry about that,
Ned
-----Original Message-----
From: Ned Dochtermann [mailto:ned.dochtermann at gmail.com]
Sent: Tuesday, April 26, 2011 1:12 PM
To: 'r-sig-mixed-models at r-project.org'; 'bonamy at horus.ens.fr'
Subject: Re: [R-sig-ME] prior specification in MCMCglmm
Celine and Pierre,
I too am still very unclear on priors however a month ago Jarrod replied to
some questions of mine regarding multi-response models with this:
"From experience I find
list(V=diag(2), nu=2, alpha.mu=c(0,0), alpha.V=diag(2)*a))
where a is something large (e.g. 1000, depending on the scale of the
data) works well for the two standard deviations and the
correlation, in terms of informativeness.
You can't use parameter expanded priors for the residual term yet,
so you will have to stick with the standard inverse-Wishart (or use
another program). Generally, data are highly informative for the residual
part so often the posterior is not very sensitive to the prior
specification. Nevertheless, you should check alternatives:
V=diag(2), nu=1.002 gives the inverse-gamma prior for the variances
with shape=scale=0.001 V=diag(2)*1e-6, nu=3 is flat for the correlation
from -1 to 1"
I think that this would generalize to G1=(list(V=diag(x), nu=1.00x)) for x
response variables for an inverse-gamma prior for G but I'm not entirely
positive that such is the case. I'm also not positive how to generalize the
flat prior (mainly would nu stay at three?), the first one I assume
generalizes to (V=diag(x), nu=x, alpha.mu=c(0,0...0x),alpha.V=diag(x)*a).
The whole thread for this discussion starts here:
https://stat.ethz.ch/pipermail/r-sig-mixed-models/2011q1/005694.html
The whole prior issue is still a mystery to me so I continue to be
uncomfortable with those approaches where they're necessary; nonetheless for
some data I'm analyzing they're basically necessary (e.g. multi-response
generalized mixed models). By and large the univariate analyses I've
conducted produce highly concordant results regardless of whether fitting
via REML, ML, Laplace or MCMC and regardless of priors for the latter. This
gives me a bit of confidence in the quantitative results and quite a bit of
confidence in the inferences. Unfortunately such comparisons aren't as
feasible (for me) with multi-response models.
Ned
--
Ned Dochtermann
Department of Biology
University of Nevada, Reno
ned.dochtermann at gmail.com
http://wolfweb.unr.edu/homepage/mpeacock/Ned.Dochtermann/
http://www.researcherid.com/rid/A-7146-2010
--
Message: 5
Date: Tue, 26 Apr 2011 09:41:58 +0200
From: "Pierre B. de Villemereuil" <bonamy at horus.ens.fr>
To: Celine Teplitsky <teplitsky at mnhn.fr>
Cc: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] prior specification in MCMCglmm
Message-ID: <4DB67746.4000607 at horus.ens.fr>
Content-Type: text/plain
Dear Ciline,
I'm not very comfortable with covariance priors, but my guess is that,
is this case, you've got to really specify a inverse-Wishart as a prior.
You should check into Hadfield's article introducing MCMCglmm, they use
something like V=diag(dimV)/4 and n=dimV. Why ? I have no idea. If your
data are not standardized, I don't think you should divide by 4 (but
then, you specify your prior as if your guess for variance components is
that they all equal 1), but for the rest...
Sorry I can not help further. Maybe somebody else would be able to help
on this subject ?
Pierre.
Le 26/04/2011 09:17, Celine Teplitsky a icrit :
> Dear Pierre,
>
> thanks a lot! It does help, but I will need time to fully understand
> the paper.
>
> Just one further question if I may, what prior would you use for a
> covariance then?
>
> Many thanks again,
>
> All the best
>
> Celine
>
>> Dear Ciline !
>>
>> One usual "non informative" prior on variance component is V=1, and
>> nu=0.002, which correspond to a inverse-Gamma(0.001,0.001). This is
>> usual, but that is not to say that it is really non informative.
>> Indeed, inverse-Gamma(e,e) is weakly informative, since the posterior
>> can depend on the choice of e.
>>
>> Concerning the WAMwiki suggestion to use the phenotypic variance to
>> set the prior, this quite not orthodox since no information coming
>> from your current dataset should be used to define the prior (but you
>> could use previous data to parametrize your prior).
>>
>> I would suggest you to refit your model considering V=1 and nu=0.002
>> as a (so-called) non informative prior. Other solutions exist like
>> using the parameter expansion and a chi2 distribution by setting
>> V=1,nu=1000,alpha.mu=0 and alpha.V=1, which is also weakly
>> informative (it has more weight in variance values less than 10).
>>
>> For more information about priors on variance component, and
>> parameter expansion, it would suggest you to read :
>>
>> 1. A. Gelman, + Prior distributions for variance parameters in
>> hierarchical models ;, /Bayesian analysis/ 1, n^o . 3 (2006):
>> 515--533.
>>
>> In the hope I'm helping.
>>
>> Pierre de Villemereuil.
>>
>>
>> Le 25/04/2011 13:26, Celine Teplitsky a icrit :
>>> Dear all,
>>>
>>> I realise that Jarrod is doing field work, but I'm really hoping
>>> someone can answer my question while he's not around.
>>>
>>> I am running animal models estimating covariances between life
>>> history traits, and I'm having trouble knowing which prior to use.
>>>
>>> Thing is, if I use a prior as described on the Wam wiki site with
>>> V=PhenotypicVar/4 (as I have 3 random effects + residual), I have
>>> very nice results, with some significant genetic correlations
>>> between some life history traits.
>>>
>>> However, one reviewer asked about prior sensitivity because CI were
>>> pretty large, so I went back to MCMCglmm course notes and saw that
>>> non informative prior were supposed to be V=diag(nbDimV)*0 and
>>> n=nbDimV-3. This led to an error message about G being ill
>>> conditioned, so I tried with diag(nbDimV)*0.001 and
>>> diag(nbDimV)*0.01 instead of diag(nbDimV)*0, and diag(nbDimV)*0.01
>>> worked... But then I have the posterior of additive genetic variance
>>> collapsing on 0 for some trait. So my guess would be that I should
>>> use those latest priors, and believe my nice results did not exist.
>>> But as Hadfield et al paper and the Wam wiki website do not
>>> recommend those priors, I am a bit confused. Could someone help me
>>> figure out what would be the right thing to do?
>>>
>>> All my apologies if this is a silly question, but I'm feeling a bit
>>> lost here
>>>
>>> Thanks a lot in advance
>>>
>>> Celine
>>>
>>
>
>
> --
>
> Celine Teplitsky
> Dipartement Ecologie et Gestion de la Biodiversiti UMR 7204
> Uniti Conservation des Esphces, Restauration et Suivi des Populations
> Case Postale 51
> 55 rue Buffon 75005 Paris
>
> Webpage :http://www2.mnhn.fr/cersp/spip.php?rubrique96
> Fax : (33-1)-4079-3835
> Phone: (33-1)-4079-3443
[[alternative HTML version deleted]]
------------------------------
Message: 6
Date: Tue, 26 Apr 2011 11:53:02 +0200
From: peter dalgaard <PDalgd at gmail.com>
To: Junqian Gordon Xu <xjqian at gmail.com>
Cc: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] interaction term in null hypothesis
Message-ID: <E1153437-757B-42E9-9A72-89EF9A368101 at gmail.com>
Content-Type: text/plain; charset="us-ascii"
On Apr 26, 2011, at 00:05 , Junqian Gordon Xu wrote:
> I have a quick question for a simple model as below:
>
>> Fix + (1 | Rand) + (1 | Rand : Fix)
>
> Which one is the null hypothesis:
>
>> 1 + (1 | Rand)
>
> or
>
>> 1 + (1 | Rand) + (1 | Rand : Fix)
>
> To me the interaction term (1 | Rand : Fix) does not make much sense if
> no fixed effect term is present in the model, but I'm not sure.
>
It does make sense, at least sometimes. For one thing, such interactions are
often aliased to aspects of the experimental design; e.g., if you have
randomized different treatment to left side and the right side of test
subjects, then the random interaction is equivalent to the (random)
difference between the two sides within the same subject. Also, you could
conceivably have a kill-or-cure drug with positive effects for some and
negative effects for others, with the question being whether the total
effect is positive or negative.
> Regards
> Gordon
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
--
Peter Dalgaard
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: (+45)38153501
Email: pd.mes at cbs.dk Priv: PDalgd at gmail.com
------------------------------
_______________________________________________
R-sig-mixed-models mailing list
R-sig-mixed-models at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
End of R-sig-mixed-models Digest, Vol 52, Issue 47
More information about the R-sig-mixed-models
mailing list