17  인과 매개 분석 (Causal mediation analysis)

경고작업 진행 중 🚧

여러분은 현재 작성 중인 R을 이용한 인과 추론의 초판본을 읽고 계십니다. 이 장은 활발히 작업 중이며 구조가 변경되거나 수정될 수 있습니다. 또한 내용이 불완전할 수 있습니다.

지금까지 우리는 노출이 결과에 미치는 총 효과(total effect)를 추정하는 데 초점을 맞추어 왔습니다. 그러나 많은 연구에서는 “왜(why)?” 또는 “어떻게(how)?”라는 질문에 답을 구하고자 합니다. 예를 들어 엑스트라 매직 아워가 대기 시간에 영향을 준다면, 그 효과 중 얼마만큼이 방문객 수 변화를 거쳐 발생하는 것일까요?

이러한 질문에 답하는 분석 방법이 바로 인과 매개 분석(causal mediation analysis)입니다. 매개 분석은 처치 효과를 두 성분으로 나눕니다:

• 직접 효과(direct effect): 매개변수를 거치지 않고 노출이 결과에 직접 미치는 효과
• 간접 효과(indirect effect): 매개변수를 거쳐 노출이 결과에 미치는 효과

library(ggdag)
library(ggokabeito)

dagify(
  Y ~ X + M + C,
  M ~ X + C,
  coords = list(
    x = c(X = 0, M = 1, Y = 2, C = 1),
    y = c(X = 0, M = 0, Y = 0, C = 1)
  ),
  labels = c(
    X = "노출",
    M = "매개변수",
    Y = "결과",
    C = "교란 요인"
  ),
  exposure = "X",
  outcome = "Y"
) |>
  tidy_dagitty() |>
  node_status() |>
  ggplot(aes(x, y, xend = xend, yend = yend, color = status)) +
  geom_dag_edges() +
  geom_dag_point() +
  geom_dag_label_repel(seed = 42) +
  scale_color_okabe_ito(na.value = "grey90") +
  theme_dag() +
  theme(legend.position = "none")

그림 17.1: 매개 분석의 기본 구조를 보여주는 DAG. 노출(X)은 매개변수(M)를 거쳐 간접적으로, 그리고 직접적으로 결과(Y)에 영향을 미칩니다.

17.1 직접 효과 추정하기 (Estimating direct effects)

인과 매개 분석에서 직접 효과와 간접 효과를 정의하는 방식은 잠재적 결과(potential outcomes) 프레임워크에 기반합니다.

자연 직접 효과(Natural Direct Effect, NDE)는 매개변수가 노출을 받지 않았을 때의 값(자연적 수준)으로 고정된 상태에서 노출을 변화시켰을 때의 효과입니다:

\[\text{NDE} = E[Y(1, M(0)) - Y(0, M(0))]\]

여기서 \(Y(1, M(0))\)는 노출이 있고 매개변수는 노출이 없을 때의 자연적 수준인 경우의 잠재적 결과입니다.

17.1.1 순차적 회귀 접근법

직접 효과를 추정하는 가장 일반적인 방법 중 하나는 순차적 회귀법(sequential regression approach)입니다. 이 접근법은 두 단계로 이루어집니다:

1. 노출과 교란 요인을 사용하여 매개변수를 예측하는 모델 적합
2. 노출, 매개변수, 그리고 교란 요인을 사용하여 결과를 예측하는 모델 적합

library(broom)
library(touringplans)
library(dplyr)

seven_dwarfs_9 <- seven_dwarfs_train_2018 |>
  filter(wait_hour == 9)

# 결과 모델: 노출, 가상의 매개변수(게시 대기 시간), 교란 요인 포함
# 여기서는 설명 목적으로 오전 8시 게시 대기를 매개변수로 활용

eight <- seven_dwarfs_train_2018 |>
  filter(wait_hour == 8) |>
  select(park_date, wait_minutes_posted_avg) |>
  rename(wait_8am = wait_minutes_posted_avg)

mediation_data <- seven_dwarfs_9 |>
  left_join(eight, by = "park_date") |>
  drop_na(wait_8am)

# 매개변수 모델: 노출 -> 매개변수
mediator_model <- lm(
  wait_8am ~ park_extra_magic_morning +
    park_ticket_season + park_close + park_temperature_high,
  data = mediation_data
)

# 결과 모델: 노출 + 매개변수 -> 결과 (직접 효과 추정)
outcome_model_direct <- lm(
  wait_minutes_posted_avg ~
    park_extra_magic_morning + wait_8am +
    park_ticket_season + park_close + park_temperature_high,
  data = mediation_data
)

tidy(outcome_model_direct, conf.int = TRUE)

# A tibble: 7 × 7
  term    estimate std.error statistic p.value conf.low
  <chr>      <dbl>     <dbl>     <dbl>   <dbl>    <dbl>
1 (Inter…  4.24e+1 16.2           2.61 9.59e-3  1.04e+1
2 park_e…  4.68e+0  2.72          1.72 8.67e-2 -6.81e-1
3 wait_8…  9.77e-2  0.0578        1.69 9.24e-2 -1.62e-2
4 park_t… -1.12e+1  2.78         -4.02 8.08e-5 -1.67e+1
5 park_t… -1.22e+1  3.88         -3.14 1.96e-3 -1.98e+1
6 park_c…  6.75e-4  0.000155      4.35 2.13e-5  3.69e-4
7 park_t… -2.75e-1  0.119        -2.31 2.17e-2 -5.10e-1
# ℹ 1 more variable: conf.high <dbl>

결과 모델에서 park_extra_magic_morning의 계수가 직접 효과의 추정치입니다. 이는 매개변수(오전 8시 대기 시간)를 고정시킨 상태에서 엑스트라 매직 아워가 오전 9시 대기 시간에 미치는 효과를 나타냅니다.

경고매개 분석의 가정

인과 매개 분석을 유효하게 수행하려면 네 가지 형태의 교란이 없어야 합니다:

1. 노출-결과 교란이 없을 것
2. 매개변수-결과 교란이 없을 것
3. 노출-매개변수 교란이 없을 것
4. 노출에 의해 유발되는 매개변수-결과 교란이 없을 것 — 이 가정이 가장 까다롭습니다

특히 네 번째 가정은 종종 충족하기 어렵습니다. 이 경우 sec-g-comp에서 소개한 G-공식과 같은 고급 방법이 필요합니다.

17.2 간접 효과 추정하기 (Estimating indirect effects)

자연 간접 효과(Natural Indirect Effect, NIE)는 노출이 없을 때의 값으로 고정된 상태에서 매개변수가 노출로 인해 변화했을 때의 효과입니다:

\[\text{NIE} = E[Y(0, M(1)) - Y(0, M(0))]\]

간접 효과는 매개변수 모델의 계수와 결과 모델에서 매개변수의 계수를 곱하여 추정할 수 있습니다. 이 곱의 방법(product of coefficients method)은 간단하고 직관적입니다:

\[\hat{\text{NIE}} = \hat{\alpha}_1 \times \hat{\beta}_2\]

여기서 \(\hat{\alpha}_1\)은 노출이 매개변수에 미치는 효과의 추정치이고, \(\hat{\beta}_2\)는 매개변수가 결과에 미치는 효과의 추정치입니다.

# 매개변수 모델에서 노출 계수
alpha_1 <- tidy(mediator_model) |>
  filter(term == "park_extra_magic_morning") |>
  pull(estimate)

# 결과 모델에서 매개변수 계수
beta_2 <- tidy(outcome_model_direct) |>
  filter(term == "wait_8am") |>
  pull(estimate)

# 직접 효과 (결과 모델의 노출 계수)
direct_effect <- tidy(outcome_model_direct) |>
  filter(term == "park_extra_magic_morning") |>
  pull(estimate)

# 간접 효과 = 곱의 방법
indirect_effect <- alpha_1 * beta_2

# 총 효과 = 직접 효과 + 간접 효과
total_effect <- direct_effect + indirect_effect

cat(
  "직접 효과:", round(direct_effect, 3), "\n",
  "간접 효과:", round(indirect_effect, 3), "\n",
  "총 효과:", round(total_effect, 3), "\n",
  "매개 비율:", round(indirect_effect / total_effect * 100, 1), "%\n"
)

직접 효과: 4.679 
 간접 효과: 2.52 
 총 효과: 7.199 
 매개 비율: 35 %

17.2.1 부트스트랩을 이용한 불확실성 추정

매개 효과의 신뢰 구간은 부트스트랩을 사용하여 계산할 수 있습니다. 직접 효과와 간접 효과 모두 두 개의 모델에서 계산되기 때문에, 델타 방법(delta method)보다 부트스트랩이 더 일반적으로 사용됩니다.

library(rsample)

compute_mediation <- function(.split) {
  .df <- as.data.frame(.split)

  # 매개변수 모델
  m_model <- lm(
    wait_8am ~ park_extra_magic_morning +
      park_ticket_season + park_close + park_temperature_high,
    data = .df
  )

  # 결과 모델
  y_model <- lm(
    wait_minutes_posted_avg ~
      park_extra_magic_morning + wait_8am +
      park_ticket_season + park_close + park_temperature_high,
    data = .df
  )

  a1 <- tidy(m_model) |>
    filter(term == "park_extra_magic_morning") |>
    pull(estimate)
  b2 <- tidy(y_model) |>
    filter(term == "wait_8am") |>
    pull(estimate)
  de <- tidy(y_model) |>
    filter(term == "park_extra_magic_morning") |>
    pull(estimate)
  ie <- a1 * b2

  tibble(
    term = c("직접효과", "간접효과", "총효과"),
    estimate = c(de, ie, de + ie)
  )
}

boots <- bootstraps(mediation_data, times = 500, apparent = TRUE) |>
  mutate(results = map(splits, compute_mediation))

int_pctl(boots, results)

# A tibble: 3 × 6
  term     .lower .estimate .upper .alpha .method   
  <chr>     <dbl>     <dbl>  <dbl>  <dbl> <chr>     
1 간접효과 -0.285      2.48   5.42   0.05 percentile
2 직접효과  0.523      4.75   9.49   0.05 percentile
3 총효과    3.08       7.23  11.6    0.05 percentile

17.3 상호작용 고려하기 (Allowing for interactions)

sec-interaction에서 살펴본 것처럼, 노출과 매개변수 사이에 상호작용이 존재할 때 상황이 더 복잡해집니다. 노출-매개변수 상호작용이 있는 경우, 자연 직접 효과와 자연 간접 효과는 더 이상 간단한 회귀 계수로 추정할 수 없습니다.

이 경우 결과 모델에 상호작용 항을 포함시키고, 반사실적 예측값을 직접 계산해야 합니다.

\[E[Y \mid X, M, C] = \beta_0 + \beta_1 X + \beta_2 M + \beta_3 (X \times M) + \beta_4 C\]

# 노출-매개변수 상호작용을 포함한 결과 모델
outcome_with_interaction <- lm(
  wait_minutes_posted_avg ~
    park_extra_magic_morning * wait_8am +
    park_ticket_season + park_close + park_temperature_high,
  data = mediation_data
)

tidy(outcome_with_interaction, conf.int = TRUE)

# A tibble: 8 × 7
  term    estimate std.error statistic p.value conf.low
  <chr>      <dbl>     <dbl>     <dbl>   <dbl>    <dbl>
1 (Inter…  4.29e+1 15.5          2.76  6.35e-3  1.22e+1
2 park_e… -2.63e+1  7.47        -3.52  5.31e-4 -4.10e+1
3 wait_8…  3.72e-2  0.0571       0.651 5.16e-1 -7.54e-2
4 park_t… -8.67e+0  2.73        -3.18  1.71e-3 -1.40e+1
5 park_t… -1.06e+1  3.74        -2.83  5.04e-3 -1.80e+1
6 park_c…  6.71e-4  0.000149     4.51  1.09e-5  3.78e-4
7 park_t… -2.74e-1  0.114       -2.41  1.70e-2 -4.99e-1
8 park_e…  5.58e-1  0.126        4.42  1.56e-5  3.09e-1
# ℹ 1 more variable: conf.high <dbl>

상호작용이 존재할 때 자연 직접 효과와 간접 효과는 반사실적 시뮬레이션을 사용하여 계산할 수 있습니다. 이 경우 G-공식(장 13) 방식이 자연스럽게 확장됩니다.

힌트mediation 패키지

R의 mediation 패키지는 상호작용이 있는 경우를 포함하여 다양한 인과 매개 분석 방법을 지원합니다. 부트스트랩 기반의 신뢰 구간과 민감도 분석 기능도 포함하고 있어 실무에서 유용합니다.

# install.packages("mediation")
library(mediation)

med_out <- mediate(
  mediator_model,
  outcome_model_direct,
  treat = "park_extra_magic_morning",
  mediator = "wait_8am",
  boot = TRUE,
  sims = 500
)
summary(med_out)