[prev in list] [next in list] [prev in thread] [next in thread] 

List:       r-sig-mixed-models
Subject:    [R-sig-ME] =?gb2312?b?tPC4tDogIEhvdyBjYW4gYSBtaXhlZCBtb2RlbCBk?= =?gb2312?b?aWZmZXJlbnRpYXRlIHRoZSBm
From:       Chen Chun <talischen () hotmail ! com>
Date:       2016-07-28 10:18:27
Message-ID: VI1PR01MB0989FB5C3D0615012C4D253AAC000 () VI1PR01MB0989 ! eurprd01 ! prod ! exchangelabs ! com
[Download RAW message or body]

[Attachment #2 (text/plain)]

Thank you Aaron for the help. It's clear to me now.


Regards,

Chun


________________________________
·¢¼þÈË: Aaron Mackey <ajmackey@gmail.com>
·¢ËÍʱ¼ä: 2016Äê7Ô 27ÈÕ 21:07
ÊÕ¼þÈË: Chen Chun
³­ËÍ: r-sig-mixed-models@r-project.org
Ö÷Ìâ: Re: [R-sig-ME] How can a mixed model differentiate the fixed effect and random \
effect when there is only one group of control in the data?

because the random effect in your model is independent of the fixed effect, the mixed \
model has six groups from which to estimate the variance of the random effect. If you \
had asked for the random effect group variances to be different between treatments  ~ \
(1|treatment:group) then you'd be in trouble.

hope that helps,
-Aaron

On Wed, Jul 27, 2016 at 4:44 PM, Chen Chun \
<talischen@hotmail.com<mailto:talischen@hotmail.com>> wrote: Dear all,


I have a data set where experiments were conducted in groups of animals and each \
group was assigned to one treatment (treatment vs. control). In my data, \
unfortunately I have only one control group, so the model would be:

lmer(out ~ treatment + (1| group),  REML = TRUE, data=dat)


I am wondering whether in this case the mixed model would still be able to provide \
unbiased estimate of the fixed effect and random group effect for the control group. \
I have written a simulation code for 5 treatment groups vs. 1 control group (see \
below). Seems that the linear mixed model is still able to provide unbiased estimate. \
Can someone give me more insight about why lmer could identify the fixed and random \
effect when there is only one group from one treatment arm?


Thanks


Chun Chen


library(lme4)
set.seed(320)
N_group <- 6
N_per_group <- 20
NO_con_group <-1
beta_t <- 3
beta_c <- 0
intercept <- 2

res <- matrix(NA, nrow=1000, ncol=2)
for(k in 1:1000) {
    random_error1 <- rnorm(N_per_group*(N_group-NO_con_group), mean = 0, sd = 1)
    random_error2 <- rnorm(N_per_group*NO_con_group, mean = 0, sd = 1)
    group_error <- rnorm(N_group, mean = 0, sd = 5)

    yt <- intercept + beta_t + rep(group_error[1:(N_group-NO_con_group)], \
each=N_per_group) + random_error1  y_c <- intercept + beta_c + \
rep(group_error[(N_group-NO_con_group+1):N_group], each=N_per_group) + random_error2  \
y <- c(yt, y_c)

    dat <- data.frame(out=y, treatment=c(rep("treat", \
N_per_group*(N_group-NO_con_group)), rep("con",N_per_group*NO_con_group)), \
group=rep(letters[1:N_group], each=N_per_group))

    fit <- lmer(out ~ treatment + (1| group),  REML = TRUE, data=dat)
    summary(fit)
    res[k,1] <- fixef(fit)[1]
    res[k,2] <- fixef(fit)[2]
}

colMeans(res)  ##---similar to intercept  and beta_t


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-mixed-models@r-project.org<mailto:R-sig-mixed-models@r-project.org> mailing \
list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models


	[[alternative HTML version deleted]]



_______________________________________________
R-sig-mixed-models@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

[prev in list] [next in list] [prev in thread] [next in thread]