[R-sig-ME] priors for univariate model with very low expected variance contributions (mcmcglmm)
Jarrod Hadfield
j.hadfield at ed.ac.uk
Sat Mar 28 08:55:59 CET 2015
Hi Rob,
Parameter expanded priors are probably a better option:
prior = list(G = list( G1 = list(V = 1, nu=1, alpha.mu=0, alpha.V=1000),
G2 = list(V = 1, nu=1, alpha.mu=0, alpha.V=1000),
G3 = list(V = 1, nu=1, alpha.mu=0, alpha.V=1000)),
R = list(V = 1, nu=0.002))
This prior is approximately flat for the standard deviation. Note that
this does not imply it is flat for the variance; it still puts quite a
bit of mass at zero, but less than the inverse-Wishart with low degree
of belief. If you want something flat on the variance scale (e.g. a
uniform distribution) you'll have to move to WinBugs or JAGS.
Parameter expansion is not implemented for the residuals so I've left
the inverse-Wishart prior on the residual variance. Usually the data
overwhelm the prior for the residual variance so you can probably be
pretty relaxed about that.
Cheers,
Jarrod
Quoting Robert Griffin <robgriffin247 at hotmail.com> on Fri, 27 Mar 2015
15:33:15 +0100:
> Dear list,
>
> I have sampled ~35 copies of a chromosome from a population because I want
> to estimate how much contribution that part of the genome makes to the
> variance in the traits. Therefore I want to estimate the additive genetic
> variance. I will do this by using making a univariate response model in
> MCMCglmm. The data for the trait was collected from 300-500 offspring per
> sampled chromosome, and measured in males. This was done across 4
> experimental *blocks *and within each *line *and *block *there were 4 *vials
> *of many individuals, each sourced from one of two *sets *of parents.
>
> Within the model there are *block*, parental set (*set*), *vial*, and *line
> *effects to model. I have done this in the following way:
>
> chain1 = MCMCglmm(trait ~1 + block,
> random = ~line + set + vial,
> rcov = ~units,
> nitt = nitt,
> thin = thin,
> burnin = burnin,
> prior = prior,
> family = "gaussian",
> start = list(QUASI = FALSE),
> data = df1)
>
>
> However, the phenotypic variance in this trait is large [var(trait) =
> ~150], and I am expecting an extremely large part of the variance to be
> environmental & measurement error (residual), and the variables of line,
> set, block, and vial to contribute very little (probably <5% of total
> variation each) - visual examination of the data suggests that there is
> almost no variance among lines, blocks, vials, or parental sets. Which
> leads me to my call for help.
>
> I am mainly concerned about how to choose priors for variances which are
> expected to be near zero (when the aim is to test if line variance is not
> 0) - can this affect the outcome of the model? How should I define my
> priors in such a case? Currently my best estimate from reading the
> literature is to use the following:
>
> prior = list(G = list( G1 = list(V = var(trait)/4, nu=0.002),
> G2 = list(V = var(trait)/4, nu=0.002),
> G3 = list(V = var(trait)/4, nu=0.002)),
> R = list(V = var(trait)/4, nu=0.002))
>
>
> Advice about the priors (and the model in general if you happen spot
> anything- e.g. should the family be Gaussian?) would be greatly appreciated,
>
> Rob
>
>
> -----------------------------
> Robert Griffin
> PhD candidate, Uppsala University
> griffinevo.wordpress.com
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> 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.
More information about the R-sig-mixed-models
mailing list