I have data for an interrupted time series (ITS) design where I have multiple sites and multiple interventions per site. The problem that I am facing is that the interventions were all implemented at differing times at the different sites. The example plot that I show depicts this by showing the 4 different sites with the vertical lines showing when each intervention (A and B) was introduced at the site. With a single intervention, I would simply create a time variable for each site such that time=0 corresponded to the start of the intervention at each site. In this way, I could create a segmented regression model in brms that would partially pool across all sites, like this:
brm(outcome ~ time_zeroed + Intervention + time_zeroed * Intervention + (time_zeroed + Intervention + time_zeroed * Intervention | site_id)
Each site would have different numbers of observations before an after zero, because they had the intervention at different points in time, but as far as the model was concerned, time would mean the same thing to all sites.
However, with multiple interventions, I canât create a zeroed time variable except for one of the interventions. It doesnât makes sense to me to make multiple time variables to include for each intervention in the model. I would like to be able to include an interaction with time for each intervention, but I am not sure how to do this without the zeroed time for the second intervention.
Should I run a separate model for each site and then use meta analysis? It wouldnât be hard to run a model with time interactions for each intervention if I only did one site at a time.
So maybe my thinking isnât correct that time must be normalized to the start of the interventions for each site. I created a dataset using the above code and then added normalized time. I have attached it here. datatest.csv (2.3 KB)
I ran a few different models on it:
m1 â brm(outcome ~ time * InterventionA + (1|site_id),data=datatest)
m2 â brm(outcome ~ time_offset * InterventionA + (1|site_id),data=datatest)
m3 â brm(outcome ~ time * InterventionA + (time * InterventionA|site_id),data=datatest)
m4 â brm(outcome ~ time_offset * InterventionA + (time_offset * InterventionA|site_id),data=datatest)
m1 and m2 are the same except one uses time and the other the normalized time.
m3 and m4 are the same except one uses time and the other the normalized time.
It appears that the slope of time and the interaction are comparable for m1 vs m2, and for m3 vs m4. The parameter for âInterventionAâ is different for m1 vs m2 and for m3 vs m4.
Also, the sd for the group level effect for the InterventionA parameter is much higher for m3 vs m4, which makes sense, because this parameter value apparently varies depending on what point in time the intervention occurs (based on the result of m1 vs m2).
So, if I am interested only in the interactions, then maybe I donât need to normalize time???
And in that case, the answer to my original question would be more straightforward.
I have seen a variety of ways to specify ITS models, including this paper that specifies the model like this:
brm(outcome ~ time + intervention + time_post + site)
where âtime_postâ would be 0 before the intervention and 1:end of study after the intervention. This works for multi-site, but when adding multiple interventions to each site, I think you would need a different âtime_postâ for each intervention⌠it seems like you could do the same thing with âtime + intervention + time:interventionâ
Anyone have ideas? @martinmodrak any thoughts? Thanks
I am definitely not an expert on this, so just some very brief observations:
That looks roughly correct to me.
More generally, I think you are running into limitations of what brms and simple linear models can do for you. In the data you shouwed, it looks like the interventions donât change the value directly, but rather change the derivative (i.e. course of the trajectory / increments per time / âŚ). So a possible another approach would be to model this as a state space model via the ctsem package or similar but that adds some complexity you might not want to get into.
Also the data look pretty strong (flat before intervention, steady decline after), so it probably doesnât matter that much how you analyze it.
Thanks @martinmodrak for the thoughtful response!
I donât have experience with state space models. I would prefer to keep it close to something seen in the literature regarding ITS model and epidemiology.
I guess my main question is this - if the sites have different times at which the implementation occurred, is it ok to use the âtime*interventionâ formulation rather than some sort of time variable for each intervention that denotes the start of the intervention for that site?
I guess it is really no different than a longitudinal study among subjects who each get the treatment at a different time.
I think what is throwing me off was looking at the âconditional_effects()â plot in brms. I donât think that plot would be valid for the interaction unless you specified the site since the interventions happened at different timepoints for each site.
Assuming the baseline level of intervention is âno interventionâ, then having time*intervention vs. a specialized variable that is 0 before intervention and then increments by time, say time_int1 in a setup like time_int1 + intervention + time would be basically equivalent in a maximum-likelihood setting. (Note that time_int1 is very similar to the time:intervention term). One would expect the time:intervention and time_int1 coefficients to be the same as both describe change post intervention for one unit of time. The difference will be that in time*intervention the first time the intervention is included in the predictor, it would be with a possibly large time, so the coefficent for the intervention simple effect would change to accomodate this. So time*intervention will IMHO lead to coefficients that are harder to interpret, although model predictions would be the same.
In a Bayesian context, thist would further interact with your priors and I think it will be easier to elicit good prior for the time_int1 + intervention + time, so another reason to prefer this.
I am not sure I understand the issue here: conditional_effects always plots the change given all the predictors (if you donât give them explicitly, brms will plug-in some defaults for you), so you are looking how the predictions change as you keep all predictors except one fixed.
Thanks! That helps. That is also what I was thinking. If you use the âtime_int1 + intervention + timeâ setup, then if you have multiple sites, and multiple interventions and want to partially pool across sites for the effect of all interventions, then I think you would have to do something like this:
brm(outcome ~ time + time_int1 + intervention1 +time_int2 + intervention2 + (time + time_int1 + intervention1 +time_int2 + intervention2 | site))
I thought all of those âtime_intâ variables would get confusing (especially if more than two interventions), but I guess it is really no different than doing âtime * interventionâ for every different intervention. And, as you say, the coefficient for intervention is not so directly interpretable for the âtime*interventionâ setup.
Also, FYI I have decided to look into the state-space models as you suggested. It looks like the âwalkerâ package is helpful.
