1 Linear Mixed Model (2): 모델 구조

1.1 3가지 기본 구조

LMM은 어떤 Random Effect를 포함하느냐에 따라 크게 세 가지 구조로 나뉜다.

구조	수식 표기 (R lme4)	의미
Random Intercept Only	`(1 \\| group)`	그룹마다 절편이 다름
Random Slope Only	`(0 + X \\| group)`	그룹마다 기울기가 다름
Random Intercept + Slope	`(1 + X \\| group)`	절편과 기울기 모두 다름

1.2 Structure 1: Random Intercept Model

1.2.1 수식

\[Y_{ij} = (\beta_0 + u_i) + \beta_1 X_{ij} + \epsilon_{ij}\]

\[u_i \sim N(0, \sigma^2_u), \quad \epsilon_{ij} \sim N(0, \sigma^2_e)\]

기울기 \(\beta_1\)은 모든 그룹이 동일 (고정)
절편은 그룹마다 다름: 그룹 \(i\)의 절편 = \(\beta_0 + u_i\)

1.2.2 직관적 이해

만족도
5 │              /  ← 사용자 3 (u₃ = +0.8, 원래 만족도 높음)
  │             /
4 │            /   ← 사용자 1 (u₁ = 0, 평균)
  │           /
3 │          /     ← 사용자 2 (u₂ = -0.8, 원래 만족도 낮음)
  │         /
2 │
  └─────────────── 시간
     기울기는 세 사람 모두 동일 (같은 slope)
     시작점(절편)만 다름

1.2.3 언제 사용하는가

그룹마다 기본 수준이 다르지만 (사람마다 기본 만족도 다름)
처치/시간에 대한 반응은 모두 동일하다고 가정할 때
가장 단순하고 일반적인 Mixed Model

1.2.4 예시: AI Agent 개인화 실험

import statsmodels.formula.api as smf

# Random Intercept: 사용자마다 기본 만족도가 다름
model_ri = smf.mixedlm(
    "satisfaction ~ personalized + week",
    data=df,
    groups=df["user_id"]   # (1 | user_id)
).fit()

print(model_ri.summary())

library(lme4)

# (1 | user_id): 사용자별 랜덤 절편
model_ri <- lmer(
  satisfaction ~ personalized + week + (1 | user_id),
  data = df
)
summary(model_ri)

결과 해석:

Random Effects:
 Groups   Name        Variance  Std.Dev.
 user_id  (Intercept)  0.42     0.65    ← 사용자 간 기본 만족도 차이 (σ²_u)
 Residual              0.58     0.76    ← 세션 간 랜덤 변동 (σ²_e)

Fixed Effects:
               Estimate Std.Error t-value
(Intercept)      3.50     0.08    43.8   ← 전체 평균 기본 만족도
personalized     0.48     0.04    12.0   ← 개인화 효과 (관심 대상)
week             0.02     0.01     2.0   ← 주당 자연 증가

ICC = 0.42 / (0.42 + 0.58) = 0.42
→ 만족도 분산의 42%가 사람 차이에서 비롯됨

1.3 Structure 2: Random Intercept + Random Slope Model

1.3.1 수식

\[Y_{ij} = (\beta_0 + u_{0i}) + (\beta_1 + u_{1i}) X_{ij} + \epsilon_{ij}\]

\[\begin{pmatrix} u_{0i} \\ u_{1i} \end{pmatrix} \sim N\left(\begin{pmatrix} 0 \\ 0 \end{pmatrix}, \begin{pmatrix} \sigma^2_{u0} & \sigma_{u01} \\ \sigma_{u01} & \sigma^2_{u1} \end{pmatrix}\right)\]

\(u_{0i}\): 그룹 \(i\)의 절편 편차 (기본 수준이 다름)
\(u_{1i}\): 그룹 \(i\)의 기울기 편차 (반응 속도가 다름)
\(\sigma_{u01}\): 절편과 기울기 편차 간의 공분산 (아래에서 설명)

1.3.2 직관적 이해

만족도
5 │    ╱ ← 사용자 3 (u₀₃=+0.5, u₁₃=+0.3: 시작도 높고 더 빠르게 증가)
  │  ╱
4 │ ╱── ← 사용자 1 (u₀₁=0, u₁₁=0: 평균적 사용자)
  │╱
3 │─── ← 사용자 2 (u₀₂=-0.3, u₁₂=-0.1: 시작도 낮고 증가도 느림)
  │
  └─────────────── 시간
     절편(시작점)도 다르고 기울기(증가 속도)도 다름

1.3.3 언제 사용하는가

그룹마다 기본 수준이 다를 뿐 아니라
처치/시간에 대한 반응 속도도 다를 것으로 예상될 때
예: 어떤 사용자는 개인화에 빠르게 반응, 어떤 사용자는 느리게 반응

1.3.4 예시: 개인화 반응 속도 차이

# Random Intercept + Slope: 사용자마다 기본 만족도도 다르고
# 시간에 따른 반응 속도도 다름
model_rs = smf.mixedlm(
    "satisfaction ~ personalized + week",
    data=df,
    groups=df["user_id"],
    re_formula="~week"   # week에 대한 랜덤 기울기 추가
).fit()

# (1 + week | user_id): 절편과 week 기울기 모두 랜덤
model_rs <- lmer(
  satisfaction ~ personalized + week + (1 + week | user_id),
  data = df
)

결과 해석:

Random Effects:
 Groups   Name        Variance  Std.Dev.  Corr
 user_id  (Intercept)  0.42     0.65
          week         0.01     0.10     -0.30  ← 절편-기울기 상관
 Residual              0.35     0.59

해석:
- σ²_u0 = 0.42: 사용자 간 기본 만족도 차이
- σ²_u1 = 0.01: 사용자 간 주당 변화량 차이 (기울기 분산)
- Corr = -0.30: 기본 만족도가 높은 사람일수록 주당 증가폭이 작은 경향
  (ceiling effect: 이미 높으니까 더 올라갈 여지가 적음)

1.3.5 절편-기울기 공분산의 해석

\(\sigma_{u01}\)의 부호에 따라 다른 패턴을 의미한다:

\(\sigma_{u01}\)	의미	예시
양수 (+)	기본 수준이 높을수록 더 빠르게 성장	원래 성적이 좋은 학생이 교육 효과도 더 큼
음수 (-)	기본 수준이 높을수록 더 느리게 성장	이미 만족도가 높은 사용자는 개선폭이 작음 (ceiling effect)
0	기본 수준과 성장 속도 무관	개인화 효과가 기본 만족도와 무관

1.4 Structure 3: 복잡한 계층 구조

실무에서는 여러 수준의 계층이 중첩될 수 있다.

1.4.1 3수준 모델 (학생 → 학급 → 학교)

\[Y_{ijk} = \beta_0 + u_{0k} + v_{0jk} + \beta_1 X_{ijk} + \epsilon_{ijk}\]

기호	의미
\(Y_{ijk}\)	학교 \(k\), 학급 \(j\), 학생 \(i\)의 성적
\(u_{0k}\)	학교 \(k\)의 랜덤 절편
\(v_{0jk}\)	학교 \(k\) 내 학급 \(j\)의 랜덤 절편
\(\beta_1\)	교수법의 고정 효과

# 3수준 모델
model_3level <- lmer(
  score ~ teaching_method + (1 | school/class),
  # (1 | school/class) = (1 | school) + (1 | class:school)
  data = df
)

분산 분해:

전체 분산 = 학교 간 분산 + 학급 간 분산(학교 내) + 학생 간 분산(학급 내)
         = σ²_school   + σ²_class            + σ²_student

예:
σ²_school  = 50  (15%)  ← 학교 간 차이
σ²_class   = 80  (24%)  ← 학교 내 학급 간 차이
σ²_student = 203 (61%)  ← 학급 내 학생 간 차이
합계        = 333 (100%)

1.4.2 Cross-classified 구조

중첩(nested)이 아니라 교차(crossed)되는 경우도 있다.

예: 학생이 여러 과목을 수강하고, 과목도 여러 학생에게 들음 (학생 ↔︎ 과목 교차)

# Cross-classified: 학생과 과목이 교차
model_crossed <- lmer(
  score ~ study_hours + (1 | student_id) + (1 | subject_id),
  data = df
)

1.5 공분산 구조 (Covariance Structure)

Random Effect의 공분산 행렬 \(\mathbf{G}\)를 어떻게 설정하느냐에 따라 모델이 달라진다.

1.5.1 Random Intercept Only

\[\mathbf{G} = \begin{pmatrix} \sigma^2_u \end{pmatrix}\]

추정할 파라미터: 1개 (\(\sigma^2_u\))

1.5.2 Random Intercept + Slope

\[\mathbf{G} = \begin{pmatrix} \sigma^2_{u0} & \sigma_{u01} \\ \sigma_{u01} & \sigma^2_{u1} \end{pmatrix}\]

추정할 파라미터: 3개 (\(\sigma^2_{u0}\), \(\sigma^2_{u1}\), \(\sigma_{u01}\))

1.5.3 Diagonal (절편-기울기 독립 가정)

공분산 없이 분산만 추정 (파라미터 절약):

\[\mathbf{G} = \begin{pmatrix} \sigma^2_{u0} & 0 \\ 0 & \sigma^2_{u1} \end{pmatrix}\]

# R에서 ||로 공분산 제거
model_diag <- lmer(
  satisfaction ~ week + (1 || user_id) + (0 + week || user_id),
  data = df
)

1.5.4 잔차의 공분산 구조

반복 측정에서는 잔차끼리도 상관될 수 있다 (시간적으로 가까운 측정치끼리 더 비슷).

구조	설명	적합 상황
Compound Symmetry (CS)	모든 시점 간 상관 동일	시간 순서 무관한 반복 측정
AR(1)	인접 시점일수록 상관 높음	시계열적 반복 측정
Unstructured	모든 시점 쌍의 상관 자유 추정	가장 유연, 파라미터 많음

library(nlme)

# AR(1) 잔차 구조
model_ar1 <- lme(
  satisfaction ~ personalized + week,
  random = ~1 | user_id,
  correlation = corAR1(form = ~week | user_id),
  data = df
)

1.6 모델 구조 선택 가이드

1단계: Random Intercept Only 로 시작
       → ICC 확인, 기본 구조 파악

2단계: Random Slope 추가 필요한가?
       → 이론적 근거: "사람마다 반응 속도가 다를 것인가?"
       → 데이터: 개인별 time plot에서 기울기 차이 확인
       → 통계: LRT로 Intercept Only vs Intercept+Slope 비교

3단계: 절편-기울기 공분산 필요한가?
       → 기본은 포함 (1 + X | group)
       → 수렴 문제 시 제거 (1 || group) + (0 + X || group)

4단계: 잔차 구조 조정
       → 시간 데이터면 AR(1) 고려
       → 기본은 독립 가정으로 시작

1.6.1 실무 권장 사항

# Step 1: Null model (ICC 계산)
m0 = smf.mixedlm("Y ~ 1", data=df, groups=df["id"]).fit()

# Step 2: Fixed effects 추가
m1 = smf.mixedlm("Y ~ X + time", data=df, groups=df["id"]).fit()

# Step 3: Random slope 추가
m2 = smf.mixedlm("Y ~ X + time", data=df, groups=df["id"],
                  re_formula="~time").fit()

# Step 4: LRT로 비교 (더 복잡한 모델이 유의하게 나은지)
from scipy.stats import chi2
lr_stat = 2 * (m2.llf - m1.llf)
p_value = 1 - chi2.cdf(lr_stat, df=2)  # 파라미터 차이 = 2
print(f"LRT: χ² = {lr_stat:.2f}, p = {p_value:.4f}")

1.7 개인별 fitted line 시각화

Mixed Model의 강점 중 하나는 개인별 예측선을 뽑을 수 있다는 것이다.

import matplotlib.pyplot as plt
import numpy as np

# 개인별 랜덤 효과 추출
random_effects = model_rs.random_effects
# {user_id: {'Intercept': u0i, 'week': u1i}, ...}

# 개인별 예측선 그리기
weeks = np.linspace(0, 4, 50)
fig, ax = plt.subplots(figsize=(10, 6))

for user_id, re in list(random_effects.items())[:10]:  # 10명만
    intercept = model_rs.fe_params['Intercept'] + re['Intercept']
    slope = model_rs.fe_params['week'] + re['week']
    y_pred = intercept + slope * weeks
    ax.plot(weeks, y_pred, alpha=0.4, color='steelblue')

# 전체 평균 (Fixed Effect만)
y_mean = model_rs.fe_params['Intercept'] + model_rs.fe_params['week'] * weeks
ax.plot(weeks, y_mean, 'r-', linewidth=3, label='Population mean')
ax.set_xlabel('Week')
ax.set_ylabel('Satisfaction')
ax.set_title('Individual trajectories + Population mean')
ax.legend()

시각화 결과:
만족도
5 │ ╱╱╱╱╱╱╱╱ ← 일부 사용자는 빠르게 증가
  │╱──────── ← 전체 평균 (빨간 선)
4 │──────── ← 일부 사용자는 완만하게 증가
  │─────
3 │────
  └─────────── 주

1.8 수식 정리

Random Intercept Model: \[Y_{ij} = \beta_0 + \beta_1 X_{ij} + u_i + \epsilon_{ij}, \quad u_i \sim N(0, \sigma^2_u)\]

Random Intercept + Slope Model: \[Y_{ij} = \beta_0 + \beta_1 X_{ij} + u_{0i} + u_{1i} X_{ij} + \epsilon_{ij}\] \[\begin{pmatrix} u_{0i} \\ u_{1i} \end{pmatrix} \sim N\left(\mathbf{0}, \mathbf{G}\right), \quad \mathbf{G} = \begin{pmatrix} \sigma^2_{u0} & \sigma_{u01} \\ \sigma_{u01} & \sigma^2_{u1} \end{pmatrix}\]

3수준 모델: \[Y_{ijk} = \beta_0 + u_k + v_{jk} + \beta_1 X_{ijk} + \epsilon_{ijk}\]

다음: [03-mixed-model-estimation.qmd] — 추정 방법 (ML vs REML), 모델 선택 (LRT, AIC, BIC)