Skip to content

Commit 612c335

Browse files
committed
add plot to vignette
1 parent db3dae1 commit 612c335

File tree

1 file changed

+21
-4
lines changed

1 file changed

+21
-4
lines changed

vignettes/practical_its.Rmd

Lines changed: 21 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -25,7 +25,7 @@ options(
2525
modelbased_select = "minimal"
2626
)
2727
28-
pkgs <- c("datawizard", "marginaleffects")
28+
pkgs <- c("datawizard", "marginaleffects", "ggplot2", "see")
2929
3030
if (!all(insight::check_if_installed(pkgs, quietly = TRUE))) {
3131
knitr::opts_chunk$set(eval = FALSE)
@@ -86,6 +86,23 @@ outcome <-
8686
rnorm(n, mean = 0, sd = 100) # Add some random noise
8787
8888
dat <- data.frame(outcome, time, event)
89+
# make event a factor for easier interpretation
90+
dat$event <- factor(dat$event, labels = c("Pre-Event", "Post-Event"))
91+
```
92+
93+
First, let's visualize the data to see the intervention's effect. We can clearly see both the immediate jump at day 200 and the change in trend afterward.
94+
95+
```{r}
96+
library(ggplot2)
97+
library(see)
98+
99+
# Visualize the simulated data
100+
ggplot(dat, aes(x = time, y = outcome, colour = event, group = event)) +
101+
geom_point(alpha = 0.5) +
102+
geom_smooth(method = "lm") +
103+
labs(title = "Interrupted Time Series Analysis", x = "Time (Days)", y = "Outcome") +
104+
theme_modern(show.ticks = TRUE) +
105+
scale_color_flat()
89106
```
90107

91108
# Modeling the Time Series
@@ -106,13 +123,13 @@ With our model fitted, we can now use `modelbased` to quantify and test the inte
106123

107124
## 1. The Interruption: Level Change
108125

109-
First, let's examine what the model estimates for the outcome right before and right after the intervention at `time = 200`. We can use `estimate_means()` to get the predicted values at `time = 199` and `time = 200` for both the factual (`event = 1`) and counterfactual (`event = 0`) scenarios.
126+
First, let's examine what the model estimates for the outcome right before and right after the intervention at `time = 200`. We can use `estimate_means()` to get the predicted values at `time = 199` and `time = 200` for both the factual (`event = Post-Event`) and counterfactual (`event = Pre-Event`) scenarios.
110127

111128
```{r}
112129
estimate_means(mod, by = c("time=c(199,200)", "event"))
113130
```
114131

115-
To directly test if the immediate "jump" or "level change" at the moment of the intervention is statistically significant, we can use `estimate_contrasts()`. We ask for the difference between `event = 1` and `event = 0` specifically at `time = 200`.
132+
To directly test if the immediate "jump" or "level change" at the moment of the intervention is statistically significant, we can use `estimate_contrasts()`. We ask for the difference between `event = Post-Event` and `event = Pre-Event` specifically at `time = 200`.
116133

117134
```{r}
118135
estimate_contrasts(mod, contrast = "event", by = "time=200")
@@ -122,7 +139,7 @@ The output shows a large, statistically significant difference of about 1020. Th
122139

123140
## 2. The Change in Trend: Slope Change
124141

125-
Next, we want to know if the intervention changed the long-term trend of the outcome. We can use `estimate_slopes()` to compute the slope of `time` for the pre-intervention period (`event = 0`) and the post-intervention period (`event = 1`).
142+
Next, we want to know if the intervention changed the long-term trend of the outcome. We can use `estimate_slopes()` to compute the slope of `time` for the pre-intervention period (`event = Pre-Event`) and the post-intervention period (`event = Post-Event`).
126143

127144
```{r}
128145
estimate_slopes(mod, trend = "time", by = "event")

0 commit comments

Comments
 (0)