Posterior predictive check for Bayesian mixed effects model for repeated measures

Interfaces rstanarm
1

I was looking into fitting a Bayesian mixed effects model for repeated measures (MMRM). The standard situation for this model would be that we have a clinical trial with multiple independent patients that are assigned into two treatment groups and the same variable is measured for each patient at multiple visits that occur at pre-specified times for the same patient (each time under the same treatment). As a result the observations within a patient are correlated, which in a standard frequentist analayis one would capture with an unstructured covariance matrix (either the same one for both treatment groups, or separate ones for each treatment group). I am attaching an plot of simulated data (see code below) illustrating what the data may look like.

mmrm_data_plot.jpeg762×765 71.7 KB

Any comments on the model implementation in my rstanarm code below are of course welcome.

What I am trying to figure out is how I can get the posterior predictive plot to be grouped by both visit and treatment group.

mmrm_pp_check.jpeg762×765 87.9 KB

patient ← rep(1:20, each=5)
visit ← rep(1:5, 20)

baseline (pre-treatment) value for each patient

base ← rep( rnorm(20, 0, 1), each=5)
epsilon ← rnorm(100, 0, 1)
treatment ← rep(0:1, each=50)

value at each visit

aval ← base + visit*treatment + treatment + epsilon

change grom baseline at each visit

chg ← aval - base
testdata ← data.frame(patient=as.factor(patient), treatment=as.factor(treatment), visit=as.factor(visit), chg=chg, aval=aval, base=base)

Illustrate what the data look like

library(ggplot2)
ggplot(testdata, aes(visit, chg, col=treatment)) + theme_grey(base_size = 24) + scale_colour_manual(values=setNames(c(rgb(0,0,1,0.5),rgb(1,0,0,0.5)), c(0,1))) + geom_jitter(width=0.1, size=3)

library(rstanarm)
stanreg1 ← stan_glmer(chg ~ 0 + visit + (visit|patient) + base + treatment + visitbase + visittreatment,
data=testdata, family=“gaussian”,
prior_intercept = student_t(df=4, 0, 1000),
prior=student_t(df=4, 0, 10),
prior_aux = student_t(4, 0, 10),
prior_covariance = decov(regularization=2, concentration=1, shape=1, scale=1),
adapt_delta=0.99)
#summary(stanreg1)
#pp_check(stanreg1,lwd=3)
pp_check(stanreg1, plotfun = “violin_grouped”, group = visit , y_draw=“points”)

Oh, and here is my brms implementation.

library(brms)
bsfit2 ← brm(chg ~ 0 + visit + (visit|patient) + base + treatment + visitbase + visittreatment,
data=testdata, family=“gaussian”,
prior=c( prior(“lkj(2)”, class = “cor”),
prior(“student_t(4,0,10)”, class = “b”),
prior(“student_t(4,0,10)”, class = “sd”)))

2

This sounds like a question for @jonah.

3

Currently pp_check is only going to let you group by one variable. But allowing two grouping variables for plots is on the to-do list. For now I think you’ll have to manually create the plot you want using ggplot2.

1 Like
4

A quick hack: Change the grouping variable to include information from both grouping variables.

pp_check(stanreg1, plotfun = "violin_grouped", group = interaction(visit, treatment), y_draw="points")

1 Like
5

Nice answer. One can also change around the order:
pp_check(stanreg1, plotfun = “violin_grouped”, group = factor( interaction(visit, treatment), levels = c(“1.0”,“1.1”,“2.0”, “2.1”, “3.0”, “3.1”, “4.0”, “4.1”, “5.0”, “5.1”)), y_draw=“points”)
But I did not manage to change to colors, but this is pretty nice already.