# C. elegans lifespan: a mixed-effects Cox model with a random effect for each plate and the trial as a stratum,
# the alternative to analyses C and D. R 4.3.3 with survival and coxme. R only.
suppressPackageStartupMessages({library(survival); library(coxme)})
d <- read.csv("lifespan.csv")
fit <- coxme(Surv(day, died) ~ condition + strata(trial) + (1 | plate), data = d)
b <- fixef(fit)[["conditiontreated"]]; se <- sqrt(vcov(fit)[1, 1]); p <- 2 * pnorm(-abs(b / se))
cat(sprintf("Mixed-effects Cox model, plate random effect, stratified by trial: hazard ratio treated / control %.2f (95%% CI %.2f to %.2f), p = %s; plate SD %.2f\n",
            exp(b), exp(b - 1.96 * se), exp(b + 1.96 * se), formatC(signif(p, 2), format = "g"), sqrt(unlist(VarCorr(fit))[[1]])))