From: Dominik Grathwohl <DominikGrathwohl_at_gmx.de>

Date: Tue 08 Mar 2005 - 12:17:46 EST

Date: Tue 08 Mar 2005 - 12:17:46 EST

Hello R-help group,

I need some suggestions in stating a mixed effects model. I would like to
fit one mixed effect model to two or more variables and than investigating
treatment effects by applying the multcomp library. I will explain my
problem along an example of weight and length measurements on infants.
Weight and length are measured at two occasions in time and for two
treatment groups, one group got some experimental formula and the other some
control formula. It is also reasonable to assume that weight and length are
correlated. Aim is to estimate the treatment effect on weight and on length
taking multiplicity into account.

Below I generated some data of the proposed properties and I applied also
two mixed effects models, one model for weight and one for length. Up to now
I’m not able to state a model for both variables together, any
suggestions?

library(nlme)

n <- 200 # n = number of subjectsid <- rep(1:n, each=2) # id = subject identification number t <- rep(0:1, times=n) # two occations in time trt <- rep(sample(rep(0:1, times=n/2)), each=2) # two treatment groups b1 <- rnorm(n, mean=0, sd=0.5) # b1 = random effect of weight

y1 <- 3.5 + rep(b1, each=2) + 0.7 * t + 0.01 * trt + 0.1 * t * trt +
rnorm(2*n,mean=0, sd=0.09)

*# y1 = weight measurements, 3.5 kg at baseline,
**# 0.7 kg more at the second time occation
**# 0.01 = the treatment effect at baseline,
**# the treatment effect is 0.1 kg for the
*

# experimental group at at the second time occation.

# The whithin subject standard deviation is 0.09.

b2 <- 0.9 * 3/sd(b1) * b1 + rnorm(length(b1), mean=0, sd=sqrt(1-0.9**2)*3)

# b2 = random effect for length with standard deviation = 3

# and correlated with the random effect of weight (b1) by r=0.9

y2 <- 49 + rep(b2, each=2) + 2 * t - 0.05 * trt + 0.5 * t * trt +
rnorm(2*n,mean=0, sd=1)

*# y2 = length measurements, 49 cm at baseline,
**# 2 cm more at the second time occation
**# 0.05 = the treatment effect at baseline,
**# the treatment effect is 0.5 cm for the
*

# experimental group at at the second time occation.

# The whithin subject standard deviation is 1.

# data frame of the data:

df <- data.frame(var=as.factor(rep(0:1, each=2*n)),
id=c(id, id), t=as.factor(c(t, t)), trt=as.factor(c(trt, trt)), y=c(y1,
y2))

# grouped data object:

gd <- groupedData(y ~ t | id, data=df)

# mixed effects model on weight:

fm1weight <- lme(y ~ t * trt, random = ~ 1 | id, data = gd, subset=var==0)
summary(fm1weight)

# mixed effect model on length:

fm1length <- lme(y ~ t * trt, random = ~ 1 | id, data = gd, subset=var==1)
summary(fm1length)

-- DSL Komplett von GMX +++ Supergünstig und stressfrei einsteigen! AKTION "Kein Einrichtungspreis" nutzen: http://www.gmx.net/de/go/dsl ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.htmlReceived on Tue Mar 08 12:32:57 2005

*
This archive was generated by hypermail 2.1.8
: Fri 03 Mar 2006 - 03:30:41 EST
*