# [R-sig-ME] heterogeneous residual variances using rcov in MCMCglmm?

David Atkins datkins at u.washington.edu
Tue Oct 25 19:16:38 CEST 2011

```Malcolm--

I haven't used MCMCglmm for this previously, though I suspect your
intuition is correct re. rcov.  However, you do need to provide a prior
for R that matches the dimensions of your (heterogeneous) residual variance.

For example, it looks to me like the following fits heterogeneous
residual variances by *sex* (below).

However, fitting a residual variance by *dam*, which has 106 unique
levels seems... unparsimonious?  (Caveat: I didn't look at the paper you
referenced, so perhaps that is precisely what is suggested.)

Hope that helps.

cheers, Dave

> prior2 <- list(R = list(V = diag(3), nu = 1.002), G = list(G1 =
list(V = 1, nu = 0.002)))

> m2 <- MCMCglmm(tarsus ~ sex, random = ~dam, rcov=~idh(sex):units,
prior = prior2, data = BTdata)

MCMC iteration = 0
[snip]
MCMC iteration = 13000
> summary(m2)

Iterations = 3001:12991
Thinning interval  = 10
Sample size  = 1000

DIC: 2015.288

G-structure:  ~dam

post.mean l-95% CI u-95% CI eff.samp
dam    0.2554    0.166   0.3539    886.7

R-structure:  ~idh(sex):units

post.mean l-95% CI u-95% CI eff.samp
Fem.units     0.5891   0.4970   0.6780    830.7
Male.units    0.6367   0.5430   0.7309   1000.0
UNK.units     0.5125   0.2885   0.7435   1000.0

Location effects: tarsus ~ sex

post.mean l-95% CI u-95% CI eff.samp  pMCMC
(Intercept)  -0.39849 -0.51948 -0.26687     1000 <0.001 ***
sexMale       0.76884  0.65555  0.87882     1000 <0.001 ***
sexUNK        0.15959 -0.06393  0.39012     1000  0.166
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

--
Dave Atkins, PhD
Research Associate Professor
Department of Psychiatry and Behavioral Science
University of Washington
datkins at u.washington.edu

Center for the Study of Health and Risk Behaviors (CSHRB)
1100 NE 45th Street, Suite 300
Seattle, WA  98105
206-616-3879
http://depts.washington.edu/cshrb/
(Mon-Wed)

Center for Healthcare Improvement, for Addictions, Mental Illness,
Medically Vulnerable Populations (CHAMMP)
325 9th Avenue, 2HH-15
Box 359911
Seattle, WA 98104
http://www.chammp.org
(Thurs)

Dear list,

Am I right in thinking that the "rcov" term in a call to MCMCglmm can
allow for heterogeneous residual (level-1) variances, where the
heterogeneity is across level-2 units? If so, how do I make use of rcov
in this way?

Consider the "BTdata" dataset included in the MCMCglmm package. The
following worked fine to estimate a two-level model, with observations
nested within dams, with a single residual variance across all dams:

prior <- list(R = list(V = 1, nu = 0.002), G = list(G1 = list(V = 1, nu
= 0.002)))
m1 <- MCMCglmm(tarsus ~ sex, random = ~dam, prior = prior, data = BTdata)

However, the variance of tarsus varies quite a bit by dam:

summary(by(BTdata, BTdata\$dam, function(x) var(x\$tarsus)))

(The variance is NA for one dam because there is only one observation
for that dam.) How therefore would I allow for heterogeneous residual
variances? I tried the following to estimate a unique residual variance
for each dam:

m2 <- MCMCglmm(tarsus ~ sex, random = ~dam, rcov=~idh(dam):units, prior
= prior, data = BTdata)

However, this fails and asks for a different set of priors, which I have
not been able to figure out. So two questions:
(1) Does this call to MCMCglmm make sense? (Or should I be trying some
other approach? I don't believe this can be done using lme4 and (RE)ML.)
(2) And how should I be specifying the priors for the MCMCglmm call?

If you're wondering why I want to do this, I am following the advice
given in: http://pan.oxfordjournals.org/content/15/2/165.

Any assistance would be much appreciated.

Cheers,
Malcolm

```