Hi All,

I'm relatively new to R and I've been playing around with lmer recently regarding which I have run into some basic doubts. My objective is to write an R function for a given lmer model such that when provided with a set of values for the predictors it should return the expected value of the dependent variable.

To this end, I created a dummy dataset with the following attributes:

y: binary (0,1)
x1: continuous
x2: continuous, range depends on value of group
group: factor with levels A,B,C and D

'data.frame':     200 obs. of  4 variables:
$ y    : int  0 1 1 1 1 1 0 0 1 0 ...
$ x1   : int  12 17 12 18 13 12 16 16 12 14 ...
$ x2   : int  774 1254 2 1121 1565 4 1093 597 1019 167 ...
$ Group: Factor w/ 4 levels "A","B","C","D": 1 1 1 1 1 1 1 1 1 1 ...

I then used lmer to fit a model to this using the call:

a <- lmer(y ~ x1 + x2 + ((1 + x1 + x2)|Group),data = modelData , family = binomial)

Now to the best of my ability to use the internet, I found that the following commands all seemed relevant:

coef(a) - returns a 4x3 matrix of coefficients: intercept, x1 and x2 coefficients by Group
ranef(a) - returns a 1x3 vector of fixed effect coefficients: intercept, x1 and x2
fixef(a) - returns a 4x3 matrix of random effect coefficients: intercept, x1 and x2 coefficients by Group

As I suspected, the sum of results of ranef and fixef equals the result of coef -

ie.

cf <- as.matrix(coef(a)[[1]])
rf <- ranef(a)[[1]]
ff <- fixef(a)

> rf + matrix(ff,nrow=4,ncol=3,byrow=TRUE) == cf
  (Intercept)   x1   x2
A        TRUE TRUE TRUE
B        TRUE TRUE TRUE
C        TRUE TRUE TRUE
D        TRUE TRUE TRUE

So, based on the values obtained using coef I created a model function which implements the following logic:

applyModel <- function(row){
      inVec <- as.numeric(row[2:3])
      group <- row[4]
      logx <- log1p(inVec + 1)
      if(group == "A"){
            ret <- sum(cf[1,1], cf[1,2:3]*logx)
      }else if(group == "B"){
            ret <- sum(cf[2,1], cf[2,2:3]*logx)
      }else if(group == "C"){
            ret <- sum(cf[3,1], cf[3,2:3]*logx)
      }else{
            ret <- sum(cf[4,1], cf[4,2:3]*logx)
      }
      return(ret)
}

af <- as.matrix(fitted(a))
b <- as.matrix(apply(as.matrix(binomData),1,applyModel))
c <- exp(b)/(1+exp(b))
binomOut <- cbind(binomData,b,af,c)
names(binomOut) <- c(names(binomData),"computed","fitted","computed_y")

Having done this, I was expecting that values in columns "fitted" and "computed_y" would match exactly. This however has turned out not to be the case and  computed_y seems to be consistently higher than the corresponding fitted value.

Pasted below are the first few observations from binomOut:


y

x1

x2

Group

computed

fitted

computed_y

1

0

12

774

A

-0.30243

0.411052

0.424962

2

1

17

1254

A

-0.14744

0.452844

0.463206

3

1

12

2

A

0.589756

0.641315

0.643309

4

1

18

1121

A

-0.08875

0.46796

0.477826

5

1

13

1565

A

-0.36802

0.396194

0.40902

6

1

12

4

A

0.521084

0.62118

0.627401

7

0

16

1093

A

-0.16609

0.44764

0.458573

8

0

16

597

A

-0.06392

0.473048

0.484025

9

1

12

1019

A

-0.34891

0.399837

0.413648

10

0

14

167

A

0.059154

0.502542

0.514784


I have a few guesses about why this might be happening:


1)      I'm coverting the logit incorrectly

2)      Numerical-errors are throwing values off

3)      My understanding of the outputs of ranef/fixef/coef is wrong

4)      I shouldn't be considering/shouldn't be considering only the output of coef while reconstructing the model equation

Could this mismatch be due to one/any of the guesses above, or is it due to something entirely different? Please help!

Thanks,

Krish.

________________________________
This email message may contain proprietary, private and confidential information. The information transmitted is intended only for the person(s) or entities to which it is addressed. Any review, retransmission, dissemination or other use of, or taking of any action in reliance upon, this information by persons or entities other than the intended recipient is prohibited and may be illegal. If you received this in error, please contact the sender and delete the message from your system.

Mu Sigma takes all reasonable steps to ensure that its electronic communications are free from viruses. However, given Internet accessibility, the Company cannot accept liability for any virus introduced by this e-mail or any attachment and you are advised to use up-to-date virus checking software.

	[[alternative HTML version deleted]]


