[R-sig-ME] estimating variance components for arbitrarily defined var/covar matrices

[Viechtbauer Wolfgang (STAT)]

Good point. Using sparse matrix methods is one way of handling that. Set 'sparse=TRUE' for rma.mv() to do so.

Best,
Wolfgang

> -----Original Message-----
> From: Jarrod Hadfield [mailto:j.hadfield at ed.ac.uk]
> Sent: Thursday, February 26, 2015 17:42
> To: Viechtbauer Wolfgang (STAT)
> Cc: r-sig-mixed-models at r-project.org; Matthew Keller
> Subject: Re: [R-sig-ME] estimating variance components for arbitrarily
> defined var/covar matrices
> 
> Hi,
> 
> If the `genetic' matrix is pedigree derived (?) the matrix is a
> correlation matrix (in the absence of inbreeding) or proportional to
> the covariance matrix (with inbreeding). However, the matrix has a lot
> of structure that can be exploited by *not* treating it as just
> another correlation matrix. Specifically, the inverse of the matrix is
> all zeros, except for the diagonal elements and those that correspond
> to mates and/or parent-offspring. ASReml, WOMBAT, MCMCglmm, pedigreem
> (?) + others all exploit this structure leading to faster more robust
> results. For small problems this might not be an issue.
> 
> Cheers,
> 
> Jarrod
> 
> Quoting "Viechtbauer Wolfgang (STAT)"
> <wolfgang.viechtbauer at maastrichtuniversity.nl> on Thu, 26 Feb 2015
> 17:24:39 +0100:
> 
> > And just to throw another package in the mix, you can also use the
> > metafor package for this. While this doesn't quite sound like a
> > meta-analysis, in the end, meta-analysis models are just
> > mixed-effects models. And in a phylogenetic meta-analysis, we add
> > random effects with known correlations matrices to the model.
> >
> > The syntax would be:
> >
> > id.e <- 1:nrow(dat)
> > id.r <- 1:nrow(dat)
> > rma.mv(y ~ <fixed effects>, V = 0, random = list(~ 1 | id.r, ~ 1 |
> > id.e), R = list(id.r = GRM), data=dat)
> >
> > where 'GRM' is the n x n matrix of similarities. The V = 0 part
> > seems a bit strange, but in meta-analytic models, we usually don't
> > estimate the error variance and instead have known sampling
> > variances (or even a known variance-covariance matrix of the
> > sampling errors). Here, we don't, so we just set that part to 0
> > (you'll get a warning that 'V appears to be not positive definite.'
> > but you can safely ignore this). This will fit the model that you
> > specified below.
> >
> > Besides the estimates of the fixed effects, the results will include
> > two variance components, one for id.r (this is VG) and one for id.e
> > (this is sigma^2_error). The default is REML estimation (method="ML"
> > is you want MLEs).
> >
> > By default, rma.mv() tries to rescale the matrix supplied via the R
> > argument into a correlation matrix. Not sure if your matrix of
> > similarities is really a correlation matrix or not. In case you
> > don't want the function to mess with the matrix, set 'Rscale=FALSE'
> > and it won't touch it.
> >
> > Best,
> > Wolfgang
> >
> > --
> > Wolfgang Viechtbauer, Ph.D., Statistician
> > Department of Psychiatry and Psychology
> > School for Mental Health and Neuroscience
> > Faculty of Health, Medicine, and Life Sciences
> > Maastricht University, P.O. Box 616 (VIJV1)
> > 6200 MD Maastricht, The Netherlands
> > +31 (43) 388-4170 | http://www.wvbauer.com
> >
> >> -----Original Message-----
> >> From: R-sig-mixed-models [mailto:r-sig-mixed-models-bounces at r-
> >> project.org] On Behalf Of Anthony R Ives
> >> Sent: Thursday, February 26, 2015 13:04
> >> To: Jarrod Hadfield; Matthew Keller
> >> Cc: r-sig-mixed-models at r-project.org
> >> Subject: Re: [R-sig-ME] estimating variance components for arbitrarily
> >> defined var/covar matrices
> >>
> >> Matthew,
> >>
> >> You should be able to do this in communityPGLMM in {pez}. Also, Steve
> >> Walker is currently working on a way to do this in lmer/glmer.
> >>
> >> Cheers, Tony
> >>
> >>
> >> On 02/26/15, Jarrod Hadfield  wrote:
> >> > Hi Matthew,
> >> >
> >> > Both MCMCglmm and asreml-r fit these models in R.
> >> >
> >> > Cheers,
> >> >
> >> > Jarrod
> >> >
> >> > Quoting Matthew Keller <mckellercran at gmail.com> on Wed, 25 Feb 2015
> >> 16:42:32 -0700:
> >> >
> >> > >Hi all,
> >> > >
> >> > >This is a typical problem in genetics and I'm trying to figure out
> >> whether
> >> > >there's any way to solve it using lmer or similar, and if not, why
> it
> >> isn't
> >> > >possible.
> >> > >
> >> > >Often in genetics, we have an n-by-n matrix (n=sample size) of
> genetic
> >> > >relationships, where the diagonal is how related you are to
> yourself
> >> (~1,
> >> > >depending on inbreeding) and off-diagonals each pairwise
> relationship.
> >> I'd
> >> > >like to be able to use lmer or some other function in R to estimate
> >> the
> >> > >variance attributable to this genetic relationship matrix. Thus:
> >> > >y = b0 + b*X + g*Z + error
> >> > >where y is a vector of observations, b is a vector of fixed
> covariate
> >> > >effects and g is a vector of random genetic effects. X and Z are
> >> incidence
> >> > >matrices for b & g respectively, and we assume g ~ N(0, VG). The
> >> variance
> >> > >of y is therefore
> >> > >var(y) = Z*Z' * VG + I*var(e)
> >> > >
> >> > >Z*Z' is the observed n-by-n genetic relationship matrix. Given an
> >> observed
> >> > >Z*Z' genetic relationship matrix, is there a way to estimate VG?
> >> > >
> >> > >I guess this boils down to, if we have an observed n-by-n matrix of
> >> > >similarities, can we use mixed models in R to get the variance in y
> >> that is
> >> > >explained by that similarity?
> >> > >
> >> > >Thanks in advance!