diff --git a/_freeze/chapters/chapter-09/execute-results/html.json b/_freeze/chapters/chapter-09/execute-results/html.json index 35cb85d..b17908b 100644 --- a/_freeze/chapters/chapter-09/execute-results/html.json +++ b/_freeze/chapters/chapter-09/execute-results/html.json @@ -1,8 +1,10 @@ { - "hash": "55b373691b598915f6ff3286b9b408d6", + "hash": "f583ec2f984196312b026132c88b28c6", "result": { - "markdown": "# Evaluating your propensity score model {#sec-eval-ps-model}\n\n\n\n\n::: {.cell}\n\n```{.r .cell-code}\nrnorm(5)\n```\n\n::: {.cell-output .cell-output-stdout}\n\n```\n[1] -1.2297 -0.4071 -1.8826 1.5237 -0.4425\n```\n\n\n:::\n:::\n\n\n## Calculating the standardized mean difference\n\n## Visualizing balance via Love Plots, boxplots, and eCDF plots\n", - "supporting": [], + "markdown": "# Evaluating your propensity score model {#sec-eval-ps-model}\n\n\n\n\n\nPropensity scores are inherently *balancing* scores. The goal is to *balance* the exposure groups across confounders. \n\n## Calculating the standardized mean difference\n\nOne way to assess balance is the *standardized mean difference*. This measure helps you assess whether the average value for the confounder is balanced between exposure groups. For example, if you have some continuous confounder, $Z$, and $\\bar{z}_{exposed} = \\frac{\\sum Z_i(X_i)}{\\sum X_i}$ is the mean value of $Z$ among the exposed, $\\bar{z}_{unexposed} = \\frac{\\sum Z_i(1-X_i)}{\\sum 1-X_i}$ is the mean value of $Z$ among the unexposed, $s_{exposed}$ is the sample standard deviation of $Z$ among the exposed and $s_{unexposed}$ is the sample standard deviation of $Z$ among the unexposed, then the standardized mean difference can be expressed as follows:\n\n$$\nd =\\frac{\\bar{z}_{exposed}-\\bar{z}_{unexposued}}{\\frac{\\sqrt{s^2_{exposed}+s^2_{unexposed}}}{2}}\n$$\nIn the case of a binary $Z$ (a confounder with just two levels), $\\bar{z}$ is replaced with the sample proportion in each group (e.g., $\\hat{p}_{exposed}$ or $\\hat{p}_{unexposed}$ ) and $s^2=\\hat{p}(1-\\hat{p})$. In the case where $Z$ is categorical with more than two categories, $\\bar{z}$ is the vector of proportions of each category level within a group and the denominator is the multinomial covariance matrix ($S$ below), as the above can be written more generally as:\n\n$$\nd = \\sqrt{(\\bar{z}_{exposed} - \\bar{z}_{unexposed})^TS^{-1}(\\bar{z}_{exposed} - \\bar{z}_{unexposed})}\n$$\n\nOften, we calculate the standardized mean difference for each confounder in the full, unadjusted, data set and then compare this to an *adjusted* standardized mean difference. If the propensity score is incorporated using *matching*, this adjusted standardized mean difference uses the exact equation as above, but restricts the sample considered to only those that were matched. If the propensity score is incorporated using *weighting*, this adjusted standardized mean difference *weights* each of the above components using the constructed propensity score weight. \n\nIn R, the `{halfmoon}` package has a function `tidy_smd` that will calculate this for a data set.\n\n\n::: {.cell}\n\n```{.r .cell-code}\nlibrary(halfmoon)\n\nsmds <- tidy_smd(\n df,\n .vars = c(confounder_1, confounder_2, ...),\n .group = exposure,\n .wts = wts # weight is optional\n)\n```\n:::\n\n\nLet's look at an example using the same data as @sec-using-ps. \n\n\n::: {.cell}\n\n```{.r .cell-code}\nlibrary(broom)\nlibrary(touringplans)\nlibrary(propensity)\n\nseven_dwarfs_9 <- seven_dwarfs_train_2018 |> filter(wait_hour == 9)\n\nseven_dwarfs_9_with_ps <-\n glm(\n park_extra_magic_morning ~ park_ticket_season + park_close + park_temperature_high,\n data = seven_dwarfs_9,\n family = binomial()\n ) |>\n augment(type.predict = \"response\", data = seven_dwarfs_9)\nseven_dwarfs_9_with_wt <- seven_dwarfs_9_with_ps |>\n mutate(w_ate = wt_ate(.fitted, park_extra_magic_morning))\n```\n:::\n\n\nNow, using the `tidy_smd` function, we can examine the standardized mean difference before and after weighting.\n\n\n::: {.cell}\n\n```{.r .cell-code}\nlibrary(halfmoon)\nsmds <- \n seven_dwarfs_9_with_wt |>\n mutate(park_close = as.numeric(park_close)) |>\n tidy_smd(\n .vars = c(park_ticket_season, park_close, park_temperature_high),\n .group = park_extra_magic_morning,\n .wts = w_ate\n )\nsmds\n```\n\n::: {.cell-output .cell-output-stdout}\n\n```\n# A tibble: 6 × 4\n variable method park_extra_magic_mor…¹ smd\n \n1 park_ticket_se… obser… 1 0.391 \n2 park_close obser… 1 0.126 \n3 park_temperatu… obser… 1 0.157 \n4 park_ticket_se… w_ate 1 0.0413\n5 park_close w_ate 1 -0.0602\n6 park_temperatu… w_ate 1 0.0613\n# ℹ abbreviated name: ¹​park_extra_magic_morning\n```\n\n\n:::\n:::\n\n\nFor example, we see above that the *observed* standardized mean difference (prior to incorporating the propensity score) for ticket season is 0.39, however after incorporating the propensity score weight this is attenuated, now 0.04.\n\nOne downside of this metric is it only quantifying balance *on the mean*, which may not be sufficient for continuous confounders, as it is possible to be balanced on the mean but severely imbalanced in the tails. At the end of this chapter we will show you a few tools for examining balance across the full distribution of the confounder.\n\n## Visualizing balance\n\n### Love Plots\n\nLet's start by visualizing these standardized mean differences. To do so, we like to use a *Love Plot* (named for Thomas Love, as he was one of the first to popularize them). The `{halfmoon}` package has a function `geom_love` that simplifies this implementation.\n\n\n::: {.cell}\n\n```{.r .cell-code}\nggplot(\n data = smds,\n aes(\n x = abs(smd), \n y = variable, \n group = method, \n color = method\n )\n) + \n geom_love()\n```\n\n::: {.cell-output-display}\n![](chapter-09_files/figure-html/unnamed-chunk-5-1.png){width=672}\n:::\n:::\n\n\n\n### Boxplots and eCDF plots\n\nAs mentioned above, one issue with the standardized mean differences is they only quantify balance on a single point for continuous confounders (the mean). It can be helpful to visualize the whole distribution to ensure that there is not residual imbalance in the tails. Let's first use a boxplot. As an example, let's use the `park_temperature_high` variable. When we make boxplots, we prefer to always jitter the points on top to make sure we aren't masking and data anomolies -- we use `geom_jitter` to accomplish this. First, we will make the unweighted boxplot.\n\n\n::: {.cell}\n\n```{.r .cell-code}\nggplot(\n seven_dwarfs_9_with_wt, \n aes(\n x = factor(park_extra_magic_morning), \n y = park_temperature_high,\n group = park_extra_magic_morning\n )\n) + \n geom_boxplot(outlier.color = NA) + \n geom_jitter() + \n labs(x = \"Extra magic morning\",\n y = \"Temperature high\")\n```\n\n::: {.cell-output-display}\n![Unweighted boxplot showing the difference in historical high temperature between days that had extra magic hours and those that did not.](chapter-09_files/figure-html/fig-boxplot-1.png){#fig-boxplot width=672}\n:::\n:::\n\n::: {.cell}\n\n```{.r .cell-code}\nggplot(\n seven_dwarfs_9_with_wt, \n aes(\n x = factor(park_extra_magic_morning), \n y = park_temperature_high,\n group = park_extra_magic_morning,\n weight = w_ate\n )\n) + \n geom_boxplot(outlier.color = NA) + \n geom_jitter() + \n labs(x = \"Extra magic morning\",\n y = \"Historic temperature high\")\n```\n\n::: {.cell-output-display}\n![Weighted boxplot showing the difference in historical high temperature between days that had extra magic hours and those that did not after incorporating the propensity score weight (ATE weight).](chapter-09_files/figure-html/fig-weighted-boxplot-1.png){#fig-weighted-boxplot width=672}\n:::\n:::\n\n\nSimilarly, we can also examine the empirical cumulative distribution function (eCDF) for the confounder stratified by each exposure group. The unweighted eCDF can be visualized using `geom_ecdf`\n\n\n::: {.cell}\n\n```{.r .cell-code}\nggplot(seven_dwarfs_9_with_wt, \n aes(x = park_temperature_high,\n color = factor(park_extra_magic_morning))) + \n geom_ecdf() + \n scale_color_manual(\n \"Extra Magic Morning\",\n values = c(\"#5154B8\", \"#5DB854\"),\n labels = c(\"Yes\", \"No\")\n ) + \n labs(x = \"Historic temperature high\",\n y = \"Proportion <= x\")\n```\n\n::: {.cell-output-display}\n![Unweighted eCDF examining the difference in distribution for historic high temperature among days that had extra magic morning hours (purple) compared to those that did not (green).](chapter-09_files/figure-html/fig-ecdf-1.png){#fig-ecdf width=672}\n:::\n:::\n\n\nThe `{halfmoon}` package allows for the additional `weight` argument to be passed to `geom_ecdf` to display a weighted eCDF plot.\n\n\n::: {.cell}\n\n```{.r .cell-code}\nggplot(seven_dwarfs_9_with_wt, \n aes(x = park_temperature_high,\n color = factor(park_extra_magic_morning))) + \n geom_ecdf(aes(weights = w_ate)) + \n scale_color_manual(\n \"Extra Magic Morning\",\n values = c(\"#5154B8\", \"#5DB854\"),\n labels = c(\"Yes\", \"No\")\n ) + \n labs(x = \"Historic temperature high\",\n y = \"Proportion <= x\")\n```\n\n::: {.cell-output-display}\n![Weighted eCDF examining the difference in distribution for historic high temperature among days that had extra magic morning hours (purple) compared to those that did not (green) after incorporating the propensity score weight (ATE).](chapter-09_files/figure-html/fig-weighted-ecdf-1.png){#fig-weighted-ecdf width=672}\n:::\n:::\n\n\nExamining @fig-weighted-ecdf, we can notice a few things. First, compared to @fig-ecdf there is improvement in the overlap between the two distributions. In @fig-ecdf, the green line is almost always noticeably above the purple, whereas in @fig-weighted-ecdf the two lines appear to mostly overlap until we reach slightly above 80 degrees. After 80 degrees, the lines appear to diverge in the weighted plot. This is why it can be useful to examine the full distribution rather than a single summary measure. If we had just used the standardized mean difference, for example, we would have likely said these two groups are balanced and moved on. Looking at @fig-weighted-ecdf suggests that perhaps there is a non-linear relationship between the probability of having an extra magic morning and the historic high temperature. Let's try refitting our propensity score model using a natural spline. We can use the function `splines::ns` for this.\n\n\n\n\n::: {.cell}\n\n```{.r .cell-code}\nseven_dwarfs_9_with_ps <-\n glm(\n park_extra_magic_morning ~ park_ticket_season + park_close +\n splines::ns(park_temperature_high, df = 5), # refit model with a spline\n data = seven_dwarfs_9,\n family = binomial()\n ) |>\n augment(type.predict = \"response\", data = seven_dwarfs_9)\nseven_dwarfs_9_with_wt <- seven_dwarfs_9_with_ps |>\n mutate(w_ate = wt_ate(.fitted, park_extra_magic_morning))\n```\n:::\n\n\nNow let's see how that impacts the weighted eCDF plot\n\n\n::: {.cell}\n\n```{.r .cell-code}\nggplot(seven_dwarfs_9_with_wt, \n aes(x = park_temperature_high,\n color = factor(park_extra_magic_morning))) + \n geom_ecdf(aes(weights = w_ate)) + \n scale_color_manual(\n \"Extra Magic Morning\",\n values = c(\"#5154B8\", \"#5DB854\"),\n labels = c(\"Yes\", \"No\")\n ) + \n labs(x = \"Historic temperature high\",\n y = \"Proportion <= x\")\n```\n\n::: {.cell-output-display}\n![Weighted eCDF examining the difference in distribution for historic high temperature among days that had extra magic morning hours (purple) compared to those that did not (green) after incorporating the propensity score weight where historic high temperature was modeled flexibly with a spline.](chapter-09_files/figure-html/fig-weighted-ecdf-2-1.png){#fig-weighted-ecdf-2 width=672}\n:::\n:::\n\n\nNow in @fig-weighted-ecdf-2 the lines appear to overlap across the whole space.\n", + "supporting": [ + "chapter-09_files" + ], "filters": [ "rmarkdown/pagebreak.lua" ], diff --git a/_freeze/chapters/chapter-09/figure-html/boxplot-1.png b/_freeze/chapters/chapter-09/figure-html/boxplot-1.png new file mode 100644 index 0000000..bcc06fa Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/boxplot-1.png differ diff --git a/_freeze/chapters/chapter-09/figure-html/ecdf-1.png b/_freeze/chapters/chapter-09/figure-html/ecdf-1.png new file mode 100644 index 0000000..31a5b2b Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/ecdf-1.png differ diff --git a/_freeze/chapters/chapter-09/figure-html/fig-boxplot-1.png b/_freeze/chapters/chapter-09/figure-html/fig-boxplot-1.png new file mode 100644 index 0000000..34dd4ef Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/fig-boxplot-1.png differ diff --git a/_freeze/chapters/chapter-09/figure-html/fig-ecdf-1.png b/_freeze/chapters/chapter-09/figure-html/fig-ecdf-1.png new file mode 100644 index 0000000..31a5b2b Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/fig-ecdf-1.png differ diff --git a/_freeze/chapters/chapter-09/figure-html/fig-weighted-boxplot-1.png b/_freeze/chapters/chapter-09/figure-html/fig-weighted-boxplot-1.png new file mode 100644 index 0000000..519abef Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/fig-weighted-boxplot-1.png differ diff --git a/_freeze/chapters/chapter-09/figure-html/fig-weighted-ecdf-1.png b/_freeze/chapters/chapter-09/figure-html/fig-weighted-ecdf-1.png new file mode 100644 index 0000000..8d9f859 Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/fig-weighted-ecdf-1.png differ diff --git a/_freeze/chapters/chapter-09/figure-html/fig-weighted-ecdf-2-1.png b/_freeze/chapters/chapter-09/figure-html/fig-weighted-ecdf-2-1.png new file mode 100644 index 0000000..3a62ca3 Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/fig-weighted-ecdf-2-1.png differ diff --git a/_freeze/chapters/chapter-09/figure-html/unnamed-chunk-5-1.png b/_freeze/chapters/chapter-09/figure-html/unnamed-chunk-5-1.png new file mode 100644 index 0000000..b085e7e Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/unnamed-chunk-5-1.png differ diff --git a/_freeze/chapters/chapter-09/figure-html/weighted-boxplot-1.png b/_freeze/chapters/chapter-09/figure-html/weighted-boxplot-1.png new file mode 100644 index 0000000..76a2f5e Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/weighted-boxplot-1.png differ diff --git a/_freeze/chapters/chapter-09/figure-html/weighted-ecdf-1.png b/_freeze/chapters/chapter-09/figure-html/weighted-ecdf-1.png new file mode 100644 index 0000000..8d9f859 Binary files /dev/null and b/_freeze/chapters/chapter-09/figure-html/weighted-ecdf-1.png differ diff --git a/chapters/chapter-03.qmd b/chapters/chapter-03.qmd index 352a6fd..026dc78 100644 --- a/chapters/chapter-03.qmd +++ b/chapters/chapter-03.qmd @@ -339,7 +339,7 @@ We know the *true* average causal effect of (unspoiled) chocolate in the sample ##### Interference -Interference would mean that an individual's exposure does not impact another's potential outcome. For example, let's say each individual has a partner, and their potential outcome depends on both what flavor of ice cream they ate *and* what flavor their partner ate. For example, in the simulation below, having a partner that received a different flavor of ice cream increases the happiness by two units. +Interference would mean that an individual's exposure impacts another's potential outcome. For example, let's say each individual has a partner, and their potential outcome depends on both what flavor of ice cream they ate *and* what flavor their partner ate. For example, in the simulation below, having a partner that received a different flavor of ice cream increases the happiness by two units. ```{r} data <- data.frame( diff --git a/chapters/chapter-09.qmd b/chapters/chapter-09.qmd index 49fe813..4e9c0bc 100644 --- a/chapters/chapter-09.qmd +++ b/chapters/chapter-09.qmd @@ -2,10 +2,207 @@ {{< include 00-setup.qmd >}} +Propensity scores are inherently *balancing* scores. The goal is to *balance* the exposure groups across confounders. + +## Calculating the standardized mean difference + +One way to assess balance is the *standardized mean difference*. This measure helps you assess whether the average value for the confounder is balanced between exposure groups. For example, if you have some continuous confounder, $Z$, and $\bar{z}_{exposed} = \frac{\sum Z_i(X_i)}{\sum X_i}$ is the mean value of $Z$ among the exposed, $\bar{z}_{unexposed} = \frac{\sum Z_i(1-X_i)}{\sum 1-X_i}$ is the mean value of $Z$ among the unexposed, $s_{exposed}$ is the sample standard deviation of $Z$ among the exposed and $s_{unexposed}$ is the sample standard deviation of $Z$ among the unexposed, then the standardized mean difference can be expressed as follows: + +$$ +d =\frac{\bar{z}_{exposed}-\bar{z}_{unexposued}}{\frac{\sqrt{s^2_{exposed}+s^2_{unexposed}}}{2}} +$$ +In the case of a binary $Z$ (a confounder with just two levels), $\bar{z}$ is replaced with the sample proportion in each group (e.g., $\hat{p}_{exposed}$ or $\hat{p}_{unexposed}$ ) and $s^2=\hat{p}(1-\hat{p})$. In the case where $Z$ is categorical with more than two categories, $\bar{z}$ is the vector of proportions of each category level within a group and the denominator is the multinomial covariance matrix ($S$ below), as the above can be written more generally as: + +$$ +d = \sqrt{(\bar{z}_{exposed} - \bar{z}_{unexposed})^TS^{-1}(\bar{z}_{exposed} - \bar{z}_{unexposed})} +$$ + +Often, we calculate the standardized mean difference for each confounder in the full, unadjusted, data set and then compare this to an *adjusted* standardized mean difference. If the propensity score is incorporated using *matching*, this adjusted standardized mean difference uses the exact equation as above, but restricts the sample considered to only those that were matched. If the propensity score is incorporated using *weighting*, this adjusted standardized mean difference *weights* each of the above components using the constructed propensity score weight. + +In R, the `{halfmoon}` package has a function `tidy_smd` that will calculate this for a data set. + ```{r} -rnorm(5) +#| eval: false +library(halfmoon) + +smds <- tidy_smd( + df, + .vars = c(confounder_1, confounder_2, ...), + .group = exposure, + .wts = wts # weight is optional +) ``` -## Calculating the standardized mean difference +Let's look at an example using the same data as @sec-using-ps. + +```{r} +library(broom) +library(touringplans) +library(propensity) + +seven_dwarfs_9 <- seven_dwarfs_train_2018 |> filter(wait_hour == 9) + +seven_dwarfs_9_with_ps <- + glm( + park_extra_magic_morning ~ park_ticket_season + park_close + park_temperature_high, + data = seven_dwarfs_9, + family = binomial() + ) |> + augment(type.predict = "response", data = seven_dwarfs_9) +seven_dwarfs_9_with_wt <- seven_dwarfs_9_with_ps |> + mutate(w_ate = wt_ate(.fitted, park_extra_magic_morning)) +``` + +Now, using the `tidy_smd` function, we can examine the standardized mean difference before and after weighting. + +```{r} +library(halfmoon) +smds <- + seven_dwarfs_9_with_wt |> + mutate(park_close = as.numeric(park_close)) |> + tidy_smd( + .vars = c(park_ticket_season, park_close, park_temperature_high), + .group = park_extra_magic_morning, + .wts = w_ate + ) +smds +``` + +For example, we see above that the *observed* standardized mean difference (prior to incorporating the propensity score) for ticket season is `r smds |> filter(variable == "park_ticket_season" & method == "observed") |> pull(smd) |> round(2)`, however after incorporating the propensity score weight this is attenuated, now `r smds |> filter(variable == "park_ticket_season" & method == "w_ate") |> pull(smd) |> round(2)`. + +One downside of this metric is it only quantifying balance *on the mean*, which may not be sufficient for continuous confounders, as it is possible to be balanced on the mean but severely imbalanced in the tails. At the end of this chapter we will show you a few tools for examining balance across the full distribution of the confounder. + +## Visualizing balance + +### Love Plots + +Let's start by visualizing these standardized mean differences. To do so, we like to use a *Love Plot* (named for Thomas Love, as he was one of the first to popularize them). The `{halfmoon}` package has a function `geom_love` that simplifies this implementation. + +```{r} +ggplot( + data = smds, + aes( + x = abs(smd), + y = variable, + group = method, + color = method + ) +) + + geom_love() +``` + + +### Boxplots and eCDF plots + +As mentioned above, one issue with the standardized mean differences is they only quantify balance on a single point for continuous confounders (the mean). It can be helpful to visualize the whole distribution to ensure that there is not residual imbalance in the tails. Let's first use a boxplot. As an example, let's use the `park_temperature_high` variable. When we make boxplots, we prefer to always jitter the points on top to make sure we aren't masking and data anomolies -- we use `geom_jitter` to accomplish this. First, we will make the unweighted boxplot. + +```{r} +#| label: fig-boxplot +#| fig.cap: "Unweighted boxplot showing the difference in historical high temperature between days that had extra magic hours and those that did not." +ggplot( + seven_dwarfs_9_with_wt, + aes( + x = factor(park_extra_magic_morning), + y = park_temperature_high, + group = park_extra_magic_morning + ) +) + + geom_boxplot(outlier.color = NA) + + geom_jitter() + + labs(x = "Extra magic morning", + y = "Temperature high") +``` + +```{r} +#| label: fig-weighted-boxplot +#| fig.cap: "Weighted boxplot showing the difference in historical high temperature between days that had extra magic hours and those that did not after incorporating the propensity score weight (ATE weight)." +ggplot( + seven_dwarfs_9_with_wt, + aes( + x = factor(park_extra_magic_morning), + y = park_temperature_high, + group = park_extra_magic_morning, + weight = w_ate + ) +) + + geom_boxplot(outlier.color = NA) + + geom_jitter() + + labs(x = "Extra magic morning", + y = "Historic temperature high") +``` + +Similarly, we can also examine the empirical cumulative distribution function (eCDF) for the confounder stratified by each exposure group. The unweighted eCDF can be visualized using `geom_ecdf` + +```{r} +#| label: fig-ecdf +#| fig.cap: "Unweighted eCDF examining the difference in distribution for historic high temperature among days that had extra magic morning hours (purple) compared to those that did not (green)." + +ggplot(seven_dwarfs_9_with_wt, + aes(x = park_temperature_high, + color = factor(park_extra_magic_morning))) + + geom_ecdf() + + scale_color_manual( + "Extra Magic Morning", + values = c("#5154B8", "#5DB854"), + labels = c("Yes", "No") + ) + + labs(x = "Historic temperature high", + y = "Proportion <= x") +``` + +The `{halfmoon}` package allows for the additional `weight` argument to be passed to `geom_ecdf` to display a weighted eCDF plot. + +```{r} +#| label: fig-weighted-ecdf +#| fig.cap: "Weighted eCDF examining the difference in distribution for historic high temperature among days that had extra magic morning hours (purple) compared to those that did not (green) after incorporating the propensity score weight (ATE)." + +ggplot(seven_dwarfs_9_with_wt, + aes(x = park_temperature_high, + color = factor(park_extra_magic_morning))) + + geom_ecdf(aes(weights = w_ate)) + + scale_color_manual( + "Extra Magic Morning", + values = c("#5154B8", "#5DB854"), + labels = c("Yes", "No") + ) + + labs(x = "Historic temperature high", + y = "Proportion <= x") +``` + +Examining @fig-weighted-ecdf, we can notice a few things. First, compared to @fig-ecdf there is improvement in the overlap between the two distributions. In @fig-ecdf, the green line is almost always noticeably above the purple, whereas in @fig-weighted-ecdf the two lines appear to mostly overlap until we reach slightly above 80 degrees. After 80 degrees, the lines appear to diverge in the weighted plot. This is why it can be useful to examine the full distribution rather than a single summary measure. If we had just used the standardized mean difference, for example, we would have likely said these two groups are balanced and moved on. Looking at @fig-weighted-ecdf suggests that perhaps there is a non-linear relationship between the probability of having an extra magic morning and the historic high temperature. Let's try refitting our propensity score model using a natural spline. We can use the function `splines::ns` for this. + + + +```{r} +seven_dwarfs_9_with_ps <- + glm( + park_extra_magic_morning ~ park_ticket_season + park_close + + splines::ns(park_temperature_high, df = 5), # refit model with a spline + data = seven_dwarfs_9, + family = binomial() + ) |> + augment(type.predict = "response", data = seven_dwarfs_9) +seven_dwarfs_9_with_wt <- seven_dwarfs_9_with_ps |> + mutate(w_ate = wt_ate(.fitted, park_extra_magic_morning)) +``` + +Now let's see how that impacts the weighted eCDF plot + +```{r} +#| label: fig-weighted-ecdf-2 +#| fig.cap: "Weighted eCDF examining the difference in distribution for historic high temperature among days that had extra magic morning hours (purple) compared to those that did not (green) after incorporating the propensity score weight where historic high temperature was modeled flexibly with a spline." + +ggplot(seven_dwarfs_9_with_wt, + aes(x = park_temperature_high, + color = factor(park_extra_magic_morning))) + + geom_ecdf(aes(weights = w_ate)) + + scale_color_manual( + "Extra Magic Morning", + values = c("#5154B8", "#5DB854"), + labels = c("Yes", "No") + ) + + labs(x = "Historic temperature high", + y = "Proportion <= x") +``` -## Visualizing balance via Love Plots, boxplots, and eCDF plots +Now in @fig-weighted-ecdf-2 the lines appear to overlap across the whole space. diff --git a/chapters/chapter-09_files/figure-html/fig-finite-sample-bias-1.png b/chapters/chapter-09_files/figure-html/fig-finite-sample-bias-1.png deleted file mode 100644 index 1d94ce7..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-finite-sample-bias-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-finite-sample-bias-2-1.png b/chapters/chapter-09_files/figure-html/fig-finite-sample-bias-2-1.png deleted file mode 100644 index 6abc919..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-finite-sample-bias-2-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-ate-hist-1.png b/chapters/chapter-09_files/figure-html/fig-sd-ate-hist-1.png deleted file mode 100644 index e6c963e..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-ate-hist-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-atm-hist-1.png b/chapters/chapter-09_files/figure-html/fig-sd-atm-hist-1.png deleted file mode 100644 index 2410084..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-atm-hist-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-ato-hist-1.png b/chapters/chapter-09_files/figure-html/fig-sd-ato-hist-1.png deleted file mode 100644 index 2ee4d03..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-ato-hist-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-att-hist-1.png b/chapters/chapter-09_files/figure-html/fig-sd-att-hist-1.png deleted file mode 100644 index bab808b..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-att-hist-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-atu-hist-1.png b/chapters/chapter-09_files/figure-html/fig-sd-atu-hist-1.png deleted file mode 100644 index d97b011..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-atu-hist-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-ate-1.png b/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-ate-1.png deleted file mode 100644 index dbbec6e..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-ate-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-atm-1.png b/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-atm-1.png deleted file mode 100644 index b574146..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-atm-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-ato-1.png b/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-ato-1.png deleted file mode 100644 index e704b96..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-ato-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-att-1.png b/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-att-1.png deleted file mode 100644 index b574146..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-att-1.png and /dev/null differ diff --git a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-atu-1.png b/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-atu-1.png deleted file mode 100644 index fad99b5..0000000 Binary files a/chapters/chapter-09_files/figure-html/fig-sd-mirror-hist-atu-1.png and /dev/null differ