[R-sig-ME] fit correlation between two random effects

[Nicolas Rode]

Dear all,

I am trying to fit the correlation between two random effects (intercepts
sampled from a bivariate normal distribution).
I was wondering if there would be a way to do this in lmer?

Thanks,

Nicolas

Here is a small example with the correlation I am tryign to estimate:

library(lme4)
library(MASS)


numStrain <- 10
data <- expand.grid(Strain=1:numStrain,NeighStrain=1:numStrain,rep=1:3)
data$Strain <- as.factor(data$Strain)
data$NeighStrain <- as.factor(data$NeighStrain)

sig1 <- 5
sig2 <- 2
rho <- 0.7
Sigma <- matrix(c(sig1^2, rho*sig1*sig2, rho*sig1* sig2, sig2^2),2,2)

## Group level random intercept for the focal strain and its effect as a
neighbour
BV <- mvrnorm(n=numStrain,mu=rep(0,2),Sigma)

Fi <- BV[,1]
names(Fi) <- levels(data$Strain)

Nj <- BV[,2]
names(Nj) <- levels(data$Strain)

## Creating a matrix with observations as rows and strain identity as
columns (300 observations x 41 strains)
L1 <- model.matrix(~-1+Strain,data)
dim(L1)
data$U_j1 <- L1 %*% Fi # strain genotypic value


## Creating a matrix with observations as rows and neighbour strain
identity as columns (240 observations x 41 strains)
L2 <- model.matrix(~-1+NeighStrain,data)
dim(L2)
data$U_j2 <- L2 %*% Nj # neighbour strain genotypic value

## Response variable with the two random effects and a random residual error
data$y <-  data$U_j1+data$U_j2+rnorm(nrow(data))

m1 <- lmer(y ~ 1 + (1 | Strain)+ (1 | NeighStrain), data)
summary(m1)

## Correlation I would like to fit directly in lmer
cor.test(unlist(ranef(m1)$Strain),unlist(ranef(m1)$NeighStrain))

	[[alternative HTML version deleted]]