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

steve.walker at utoronto.ca steve.walker at utoronto.ca
Thu Feb 26 23:20:53 CET 2015

Good point Tony.  I'll add that to the README.


Quoting Anthony Ives <arives at wisc.edu>:

> Steve,
> Just a quick follow-up. Even though your work (and communityPGLMM)   
> does not allow transformations in the covariance matrices themselves  
>  (like the OU), it does in effect perform Pagel's lambda in the form  
>  sigma2*V + I ~= lambda*V + (1-lambda)*I.
> Cheers, Tony
> Anthony Ives
> Department of Zoology
> 459 Birge Hall (4th floor, E end of bldg)
> UW-Madison
> Madison, WI 53706
> 608-262-1519
>> On Feb 26, 2015, at 1:09 PM, Steve Walker <steve.walker at utoronto.ca> wrote:
>> Thanks Tony.  For those interested, you can check out the   
>> development of my work on phylogenetic models in lme4 here:
>> https://github.com/stevencarlislewalker/lme4ord
>> The section in the README file about phylogenetic models gives an   
>> example and briefly describes the kinds of models that can be fitted.
>> Regarding speed, the glmerc (with a 'c' for known covariance over   
>> grouping factor levels) function can fit models similar to the one   
>> in the README but with a 500 tip phylogeny in about 20 seconds on   
>> my pretty standard macbook pro.  lme4 uses sparse matrices for all   
>> random effect-related matrix computations, so most of the speed   
>> comes from this infrastructure.  However, it may be possible to   
>> speed things up more by exploiting Kronecker-product-type   
>> structure, but this would require some substantial additions to   
>> lme4, I think.
>> I should also note that the models fitted by glmerc do _not_ allow  
>> for fitting parameters that scale branch lengths (e.g.  
>> Ornstein-Uhlenbeck).  The fundamental challenge with these models  
>> in lme4 is that lme4 is based on a Cholesky factor parameterization  
>> of the random effects covariance matrix.  To parameterize the  
>> covariance matrix on the covariance or branch length scale, one  
>> would need to compute a Cholesky factor every time the deviance  
>> function is evaluated.  This could get costly.  The recent methods  
>> of Ho and Ané (2014) in the phylolm package might help here.   
>> However, in their current form these methods only apply to  
>> quadratic terms involving the inverse covariance matrix and to the  
>> log determinant of the covariance matrix.  These computations are  
>> certainly important in mixed modelling in general, but  
>> unfortunately not _so_ much for lme4, because of the Cholesky  
>> parameterization.  However, I bet that their ideas could be  
>> modified to be applicable to lme4.  In particular I
>>  bet one could compute a series of small Cholesky decompositions at  
>> the tips and then iteratively rank-one update them by traversing  
>> the tree to the root.  My guess is that this procedure would scale  
>> linearly, as do the other methods of Ho and Ané (2014).
>> Cheers,
>> Steve
>> On 2015-02-26 7:04 AM, Anthony R Ives wrote:
>>> 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!
>>>>> [[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.
>>>> _______________________________________________
>>>> R-sig-mixed-models at r-project.org mailing list
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>> _______________________________________________
>>> R-sig-mixed-models at r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

More information about the R-sig-mixed-models mailing list