Reputation: 77
I am fairly new to R and mixed models analysis.
I wanted to produce a single estimate for the change in the variable ln_ahr
based on time
for every individual. I believe this can be thought of as a slope of change in time. Here is how my data is structured (long format):
v001 ln_ahr time
13404 28337 0.28438718 0
13405 28337 NA 3
13406 28337 NA 6
13407 28337 1.05015991 9
13408 28337 NA 12
13409 28337 1.33345188 15
13410 28337 NA 19
13413 28355 1.14904314 0
13414 28355 NA 3
13415 28355 1.06546008 6
13416 28355 NA 9
13417 28355 1.17865500 12
13418 28355 2.84949593 15
13423 29983 0.07015499 0
13424 29983 0.21056477 3
13426 29983 0.36125306 9
13427 29983 0.66139848 12
13428 29983 0.16962391 16
where v001
is the subject identifier.
I tried to calculate the slope using the nlme
package in R as:
slope <- lme(ln_ahr~time,random=~1+time|v001,
data=restructured,na.action="na.omit")
and I tried obtained the ranef(slope)
and coef(slope)
values. I read that the coef(slope)
values "computes the sum of the fixed and random effects coefficients for each explanatory variable of each grouping factor" thus I believed that printing out the coefficients for time (leaving out the intercept values) would give me an estimate of each individual's change in ln_ahr
over time and I can use that as my "slope" or estimate for change in ln_ahr
.
Time is calculated as years where time
0 indicates the first year of ln_ahr
measurement; everyone is measured every three years.
I am wondering if this a correct approach at all or if I did it correctly; if not what are your suggestions?
Upvotes: 0
Views: 594
Reputation: 226087
The basic answer is "yes"; the numbers returned by lme4::coef()
are the estimated subject-specific parameters.
Variants of this example are around on the internet, but:
Fit one of the basic random-slopes examples from lme4
:
library(lme4)
fm1 <- lmer(Reaction~Days+(Days|Subject),sleepstudy)
Extract estimated intercepts and slopes for each group (Subject):
d1 <- coef(fm1)$Subject
d1$Subject <- rownames(d1)
For comparison, fit a separate model for each group and extract subject-specific slopes and intercepts:
fm2 <- lmList(Reaction~Days|Subject,sleepstudy)
d2 <- coef(fm2)
d2$Subject <- rownames(d2)
Plot (random-effects estimates as solid lines; fixed-effects dashed)
library(ggplot2); theme_set(theme_bw())
gg0 <- ggplot(sleepstudy,aes(Days,Reaction,
colour=Subject))+
geom_point()+
geom_abline(data=d1,
mapping=aes(intercept=`(Intercept)`,
slope=Days,colour=Subject))+
geom_abline(data=d2,linetype=2,
mapping=aes(intercept=`(Intercept)`,
slope=Days,colour=Subject))
Doug Bates doesn't like these "spaghetti plots" with all groups in a single panel, he prefers facets:
gg0+facet_wrap(~Subject)+
theme(panel.spacing=grid::unit(0,"lines"))
(Ideally we would also order the subjects in some non-arbitrary way, e.g. by slope)
Upvotes: 2