Nothing
## ------------------------------------------------------------------------
library(slim)
## ------------------------------------------------------------------------
data(dialysis)
str(dialysis)
## ------------------------------------------------------------------------
library(ggplot2)
print(p <- ggplot(dialysis, aes(x = month, y = renalfn, group = id)) +
geom_line(alpha = 0.1) + facet_grid(~ group))
## ------------------------------------------------------------------------
naive_means <- dialysis[, list(id = group, renalfn = mean(renalfn)),
by = list(group, month)]
print(p <- p + geom_point(data = naive_means))
## ------------------------------------------------------------------------
slim_basic_fit <- slim(renalfn ~ group + month, dialysis)
summary(slim_basic_fit)
## ---- results = "hide"---------------------------------------------------
# identify one fully observed individual in each group
dialysis[, fully_observed := (length(month) == 5), by = id]
group_reps <- dialysis[(fully_observed), list(id = min(id)), by = group]
setkey(group_reps, id)
dialysis[, slim_basic := fitted(slim_basic_fit)]
## ------------------------------------------------------------------------
print(p <- p + geom_line(aes(y = slim_basic), data = dialysis[group_reps]))
## ---- results = "hide"---------------------------------------------------
slim_intercept_slope_fit <- slim(renalfn ~ group + month, dialysis,
limit = ~ 1 + month)
dialysis[, slim_intercept_slope := fitted(slim_intercept_slope_fit)]
## ------------------------------------------------------------------------
print(p <- p + geom_line(aes(y = slim_intercept_slope),
data = dialysis[group_reps], linetype = "dashed"))
## ------------------------------------------------------------------------
small_example <- data.table(month = c(0, 3, 18))
z <- model.matrix(~ 1, small_example)
z %*% t(z)
z <- model.matrix(~ 1 + month, small_example)
z %*% t(z)
## ------------------------------------------------------------------------
vcov(slim_basic_fit)
vcov(slim_basic_fit, empirical = FALSE) # a bad idea for unmodelled covariances!
summary(slim_basic_fit)
summary(slim_basic_fit, empirical = FALSE) # still a bad idea
## ------------------------------------------------------------------------
library(lme4)
lmer_fit <- lmer(renalfn ~ group + month + (1 + month | id), dialysis)
slim_lmer_fit <- slim(renalfn ~ group + month, dialysis, covariance = lmer_fit)
summary(slim_lmer_fit)
summary(slim_lmer_fit, empirical = FALSE) # this now makes sense
## ------------------------------------------------------------------------
library(jmcm)
jmcm_fit <- jmcm(renalfn | id | month ~ group | 1, data = dialysis,
triple = rep(2L, 3), cov.method = "mcd")
slim_jmcm_fit <- slim(renalfn ~ group + month, dialysis, covariance = jmcm_fit)
summary(slim_jmcm_fit)
summary(slim_jmcm_fit, empirical = FALSE)
## ---- results = "hide"---------------------------------------------------
dialysis[, slim_lmer := fitted(slim_lmer_fit)]
dialysis[, slim_jmcm := fitted(slim_jmcm_fit)]
## ------------------------------------------------------------------------
print(p <- p + geom_line(aes(y = slim_lmer), data = dialysis[group_reps],
colour = "darkblue") +
geom_line(aes(y = slim_jmcm), data = dialysis[group_reps],
colour = "darkred"))
## ---- results = "hide"---------------------------------------------------
slim_pascal_fit <- slim(renalfn ~ group + month, dialysis,
covariance = "pascal")
dialysis[, slim_pascal := fitted(slim_pascal_fit)]
## ------------------------------------------------------------------------
print(p <- p + geom_line(aes(y = slim_pascal), data = dialysis[group_reps],
colour = "cyan"))
## ------------------------------------------------------------------------
extract_se <- function(fit) sqrt(diag(vcov(fit)))
sapply(list(se_basic = slim_basic_fit, se_pascal = slim_pascal_fit), extract_se)
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.