§ 9.3–9.4 — Stratified Cox Model & Left Truncation

1 도입 — 두 가지 확장의 동기

§9.2에서 시간의존 공변량을 활용하여 PH 가정의 위반을 검정했다. 검정 결과 PH가 깨진 변수가 발견되면, 다음 두 전략이 가능하다:

전략	적용 조건	핵심 아이디어
층화(Stratification)	PH 위반 변수가 범주형(또는 범주화 가능)이다	해당 변수의 각 수준에 별도의 기저 위험을 할당한다
좌절단(Left Truncation)	특정 중간 사건 이후의 하위 집단만 분석 대상이다	위험집합 정의를 “진입 시점 이후”로 수정한다

두 기법은 독립적으로도, 결합적으로도 사용된다. BMT 데이터에서는 MTX 사용 여부로 층화하고, 혈소판 회복 환자만 대상으로 좌절단 분석을 수행하는 것이 대표적이다 (Klein & Moeschberger, 2003, Ch.9).

2 층화 비례위험 모형 (§ 9.3)

2.1 모형 정의 — 식 (9.3.1)

$s$ 개의 층(stratum)이 있을 때, $j$ 번째 층에 속한 개인의 위험률은 다음과 같다:

\[ h_j[t \mid \mathbf{Z}(t)] = h_{0j}(t)\, \exp[\boldsymbol{\beta}^\top \mathbf{Z}(t)], \quad j = 1, \ldots, s \]

직관 — “기저 위험은 따로, 효과는 함께”

각 층은 고유한 기저 위험함수 $h_{0j}(t)$ 를 가진다 — 층 간 위험의 형태(shape)가 달라도 된다.

그러나 공변량의 회귀 계수 $\boldsymbol{\beta}$ 는 모든 층에서 동일하다. 즉, “FAB 등급이 위험을 2배 높인다”는 효과는 MTX 투여 여부와 무관하게 일정하다고 가정한다.

이는 분산분석의 블록(block) 설계와 유사하다: 블록 효과(기저 위험)는 자유롭게 변하게 두고, 처치 효과(β)만 추정한다.

2.2 층화 편로그우도 — 식 (9.3.2)

전체 로그 편우도는 각 층의 로그 편우도를 단순히 합산한다:

\[ LL(\boldsymbol{\beta}) = LL_1(\boldsymbol{\beta}) + LL_2(\boldsymbol{\beta}) + \cdots + LL_s(\boldsymbol{\beta}) \]

여기서 $LL_j(\boldsymbol{\beta})$ 는 $j$ 번째 층의 데이터만 사용하여 구성한 Ch.8의 로그 편우도(식 8.3.2)이다.

직관 — 왜 합산이 성립하는가

각 층에서 위험집합은 해당 층 내부의 개인들로만 구성된다. 층 $j$ 의 사건에 대해 층 $k$ 의 개인은 위험에 노출되지 않는다. 따라서 각 층의 편우도는 독립적으로 구성되며, 전체 정보량은 단순합이 된다.

이것이 “각 층에 별도의 기저 위험”을 허용하는 메커니즘이다 — $h_{0j}(t)$ 는 추정할 필요 없이 편우도에서 소거된다(Ch.8의 Cox 편우도 소거 원리와 동일).

2.3 추정 절차

각 층에서 사건 시점별 위험집합을 구성한다 (층 내부만)
각 층의 편우도 기여를 합산하여 $LL(\boldsymbol{\beta})$ 를 형성한다
Newton-Raphson으로 $\boldsymbol{\beta}$ 를 추정한다 — 도함수도 층별 합산으로 구한다
각 층에서 개별적으로 생존함수를 추정한다 (§8.8 Breslow 방법)

2.4 핵심 가정 — 공통 β의 타당성 검정

층화 모형의 핵심 가정은 공변량 효과가 층 간에 동일하다는 것이다. 이를 검정하는 두 가지 방법이 있다.

2.4.1 우도비 검정

단계	모형	로그 편우도
귀무	층화 모형 (공통 $\boldsymbol{\beta}$ )	$LL(\mathbf{b})$
대립	각 층 별도 모형 ( $\boldsymbol{\beta}_j$ 별도)	$\sum_{j=1}^s LL_j(\mathbf{b}_j)$

검정통계량:

\[ \chi^2 = -2\!\left[LL(\mathbf{b}) - \sum_{j=1}^{s} LL_j(\mathbf{b}_j)\right] \sim \chi^2_{(s-1)p} \]

직관

대립 모형에서는 각 층이 고유한 $\boldsymbol{\beta}_j$ 를 갖는다 — 자유도가 $s \times p$ 개.

귀무 모형은 $\boldsymbol{\beta}_1 = \cdots = \boldsymbol{\beta}_s$ 로 제한하므로 자유도가 $p$ 개.

차이 자유도 $(s-1)p$ 만큼의 카이제곱으로 “층 간 효과 차이가 우연인가”를 검정한다.

2.4.2 Wald 검정 (개별 계수 비교)

각 층의 추정치가 점근적으로 독립이므로(정보 행렬이 block diagonal), $i$ 번째 공변량에 대해:

\[ \chi^2_i = \frac{(b_{1i} - b_{2i})^2}{\text{SE}^2(b_{1i}) + \text{SE}^2(b_{2i})} \]

이는 각 공변량별로 “두 층에서의 효과 크기가 같은가”를 1 자유도로 검정한다.

2.5 예제: BMT 데이터 — MTX 층화 (Example 9.1 계속)

§9.2에서 MTX가 PH 가정을 위반함을 확인했다. 이를 층화 변수로 처리한다.

모형: MTX 사용(40명) / 미사용(97명) 두 층, 공변량 = Group + FAB + Age + $Z_P(t)$

결과 (Table 9.7):

공변량	$b$	SE$(b)$	Wald $\chi^2$	$p$
$Z_1$ : AML low-risk	$-0.990$	0.367	7.30	0.0069
$Z_2$ : AML high-risk	$-0.363$	0.371	0.96	0.3280
$Z_3$ : FAB Grade 4/5	$0.892$	0.284	9.90	0.0017
$Z_4$ : 환자 연령$-$28	$0.010$	0.020	0.23	0.6305
$Z_5$ : 공여자 연령$-$28	$-0.001$	0.018	0.01	0.9373
$Z_6$ : $Z_4 \times Z_5$	$0.003$	0.001	7.43	0.0064
$Z_P(t)$ : 혈소판 회복	$-1.003$	0.345	8.48	0.0036

질병 유형 전체 검정: $H_0: \beta_1 = \beta_2 = 0$ → Wald $\chi^2 = 8.92$ ( $p = 0.012$ ) — 유의.

해석

층화 전후의 추정치가 매우 유사하다 — MTX에 의한 기저 위험의 형태 차이는 존재하지만, 그것이 다른 공변량의 효과 추정을 크게 왜곡하지는 않았다. 그럼에도 층화 모형이 올바른 이유는, PH 위반 상태에서의 편우도는 이론적으로 일관 추정량을 보장하지 않기 때문이다.

2.6 공통 β 가정 검증

우도비 검정:

\[ -2\bigl[-303.189 - ((-219.677) + (-80.467))\bigr] = 6.09, \quad \text{df} = 7, \quad p = 0.53 \]

개별 Wald 검정(Table 9.8)에서도 7개 공변량 모두 두 층 간 계수 차이가 유의하지 않다 (최대 $p = 0.17$ ).

결론: 층화 모형의 공통 $\boldsymbol{\beta}$ 가정은 타당하다.

2.7 층화 모형의 또 다른 활용 — 대응표본(Matched Pairs)

층화 모형은 대응 설계(matched design)에도 적용된다. 각 대응 쌍을 하나의 층으로 정의하면:

대응 변수(병원, 관해 상태 등)의 효과는 $h_{0j}(t)$ 로 자유롭게 흡수된다
치료 효과만 $\beta$ 로 추정된다

2.7.1 예제: 6-MP vs 위약 (Example 9.3)

급성 백혈병 소아 42명을 관해 상태로 대응시키고, 쌍 내에서 6-MP 또는 위약에 무작위 배정했다.

21개 쌍 → 21개 층
단일 공변량: $Z = 1$ (6-MP), $Z = 0$ (위약)
추정: $\hat{\beta} = -1.792$ , SE $= 0.624$
우도비 검정: $\chi^2 = 11.89$ ( $p = 0.0006$ )
상대 위험: $e^{-1.792} = 0.167$ , 95% CI: [0.049, 0.566]

직관

위약군은 6-MP군보다 재발 위험이 2~20배 높다. Score 검정 통계량은 §7.7의 층화 로그순위 검정과 정확히 동일하다 — 층화 Cox 모형의 점수(score) 검정은 층화 로그순위 검정을 일반화한 것이다.

3 좌절단 회귀 (§ 9.4)

3.1 좌절단이란 무엇인가

좌절단(left truncation)은 개인이 시간 원점보다 늦게 관찰을 시작하는 상황이다. 다른 이름으로 지연 진입(delayed entry) 이라 한다.

구분	정의	예시
좌절단	진입 시점 $V$ 이전에 사건이 발생한 개인은 표본에 포함되지 않는다	은퇴 커뮤니티 입소 후 관찰 시작 — 입소 전 사망자는 연구에 나타나지 않는다
중도절단(right censoring)	관찰 종료 시점까지 사건이 발생하지 않는다	연구 종료 시 생존 중인 환자

직관 — 생존 편향

좌절단은 본질적으로 생존 편향(survival bias) 을 내포한다. “진입 시점까지 살아 있어야 관찰에 포함된다”는 조건이 붙기 때문이다.

Channing House 예시: 은퇴 커뮤니티에 입소할 만큼 오래 산 사람만 연구에 들어간다. 68세에 입소한 사람은 68세 이전의 사망 위험에 기여하지 않는다 — 왜냐하면 그 시점에는 아직 “관찰 대상”이 아니었기 때문이다.

3.2 두 가지 전형적 시나리오

시나리오 A: 연령 척도 분석

사건 시간 $X$ = 사망 연령
진입 시점 $V$ = 연구 등록 연령
예: Channing House — 입소 연령 $V_i$ , 사망 또는 중도절단 연령 $T_i$

시나리오 B: 중간 사건 후 하위 집단 분석

사건 시간 $X$ = 랜드마크(이식일 등)로부터 사망까지의 시간
진입 시점 $V$ = 중간 사건(혈소판 회복) 발생 시점
예: BMT — 혈소판 회복 환자만의 무병생존 분석

3.3 조건부 위험률과 독립 가정

좌절단 데이터에서의 조건부 위험률은 다음과 같다:

\[ h(t \mid \mathbf{Z}, X > V) = \frac{P(X = t \mid \mathbf{Z}, X > V)}{P(X \ge t \mid \mathbf{Z}, X > V)} \]

핵심 가정: 조건부 독립

사건 시간 $X$ 와 진입 시점 $V$ 가 공변량 $\mathbf{Z}$ 가 주어졌을 때 독립이면:

\[ h(t \mid \mathbf{Z}, X > V) = h(t \mid \mathbf{Z}) \]

즉, 조건부 위험률이 비조건부 위험률과 동일해진다 (Andersen et al., 1993).

이 가정이 깨지면 — 예를 들어 “오래 살아서 입소한 사람일수록 더 건강하다”는 선택 효과가 존재하면 — 좌절단 보정만으로는 편향을 제거할 수 없다.

3.4 위험집합의 재정의

좌절단에서의 핵심 수정은 위험집합 $R(t)$ 의 정의이다:

\[ R(t) = \{j \mid V_j < t \le T_j\} \]

직관 — 표준 위험집합과의 차이

표준 (중도절단만): $R(t) = \{j \mid t \le T_j\}$ — 시점 0부터 관찰 중인 모든 개인.

좌절단 보정: $R(t) = \{j \mid V_j < t \le T_j\}$ — 시점 $t$ 에서 이미 진입했으며 아직 관찰 종료되지 않은 개인만 포함.

이 단순한 변경이 편우도 전체에 파급된다 — 분모의 합산 범위가 달라지므로, “아직 나타나지 않은 사람은 위험에 노출된 것으로 세지 않는다”는 논리가 반영된다.

이 수정된 위험집합을 Ch.8과 §9.2–9.3의 편우도에 대입하면, 나머지 추정·검정 절차는 동일하게 적용된다.

3.5 예제 A: Channing House — 성별과 사망 연령 (Example 9.4)

캘리포니아 은퇴 커뮤니티 Channing House의 거주자 사망 데이터이다.

사건 시간: 사망 연령
좌절단 시점: 입소 연령
공변량: 성별 ( $Z = 1$ 남성, $Z = 0$ 여성)

위험집합은 시간(연령) $t$ 에서 “입소 연령이 $t$ 미만이며 아직 사망/퇴소하지 않은” 거주자로 구성된다. 집합 크기는 Figure 4.10처럼 연령에 따라 비단조적으로 변한다 — 70대까지 증가하다가 80대 이후 급감한다.

결과:

$\hat{\beta}_{\text{gender}} = 0.316$ , SE $= 0.173$
Wald $p = 0.068$ — 남녀 간 유의한 생존 차이 없음

직관

만약 좌절단을 무시하고 연령을 단순 공변량으로 넣었다면, “일찍 입소한 건강한 고령자”가 전체 연령 효과를 왜곡할 수 있다. 좌절단 처리는 각 연령 시점에서 실제로 관찰 가능했던 사람들만으로 위험집합을 구성함으로써 이 편향을 교정한다.

3.6 예제 B: BMT — 혈소판 회복 후 무병생존 (Example 9.5)

혈소판이 회복된 환자만을 대상으로 무병생존을 분석한다.

사건 시간: 이식일로부터 사망/재발까지의 시간
좌절단 시점: 혈소판 회복 시점 $V_i$
층화 변수: MTX 사용 여부
공변량: Group, FAB, Age

모형:

\[ h(t \mid \mathbf{Z}, \text{MTX}) = h_{0,\text{MTX}}(t)\, \exp(\boldsymbol{\beta}^\top \mathbf{Z}), \quad t > V_i \]

결과 (Table 9.9):

공변량	$b$	SE$(b)$	Wald $\chi^2$	$p$
$Z_1$ : AML low-risk	$-1.752$	0.438	16.03	<0.0001
$Z_2$ : AML high-risk	$-0.750$	0.408	3.39	0.0657
$Z_3$ : FAB Grade 4/5	$1.278$	0.325	15.46	<0.0001
$Z_4$ : 환자 연령$-$28	$0.042$	0.022	3.51	0.0611
$Z_5$ : 공여자 연령$-$28	$-0.035$	0.021	2.80	0.0943
$Z_6$ : $Z_4 \times Z_5$	$0.002$	0.001	3.49	0.0617

질병 유형 전체 검정: Wald $\chi^2 = 18.27$ (2 df, $p < 0.0001$ ) — 강하게 유의.

해석 — §9.2 시간의존 결과와의 비교

§9.2의 시간의존 공변량 분석에서는 혈소판 회복 전후를 모두 포함하여 $Z_P(t)$ 의 효과를 추정했다.

여기서는 혈소판 회복 후의 하위 집단만 분석한다 — 서로 다른 질문에 답하는 것이다:

§9.2: “혈소판 회복이 위험에 영향을 주는가?” → Yes, HR = 0.32
§9.4: “회복 후 환자들 사이에서, 질병 유형이 예후를 결정하는가?” → Yes, AML low-risk의 HR = 0.17

좌절단 분석은 이미 이정표를 통과한 집단 내부의 예후 인자를 식별하는 데 적합하다.

3.7 Practical Notes

항목	내용
연령과 좌절단	연령이 시간 척도일 때 좌절단으로 처리하면, 연령을 공변량으로 넣을 필요가 없다
소프트웨어	R `Surv(V, T, event)` — 왼쪽에 진입 시점, 오른쪽에 관찰 종료 시점
생존함수 추정	좌절단 + 고정 공변량 모형에서는 §8.8의 Breslow 방법을 그대로 적용할 수 있다
독립 가정 위반	$X \perp V \mid \mathbf{Z}$ 가 깨지면 추정량이 편향된다 (Keiding, 1992 참조)

4 층화 + 좌절단의 결합

BMT Example 9.5는 두 기법을 동시에 사용하는 전형적인 사례이다:

전체 BMT 환자 (n=137)
    │
    ├─ 층화: MTX 사용 (n=40) / 미사용 (n=97)
    │       → 각 층에 별도의 h₀(t)
    │
    └─ 좌절단: 혈소판 회복 환자만 진입
              → R(t) = {j | V_j < t ≤ T_j, 같은 층}

편우도 구성:

\[ LL(\boldsymbol{\beta}) = \sum_{j=1}^{s} \sum_{i: \text{사건 in 층 } j} \left[\boldsymbol{\beta}^\top \mathbf{Z}_{(i)}(t_i) - \log\!\sum_{k \in R_j(t_i)} e^{\boldsymbol{\beta}^\top \mathbf{Z}_k(t_i)}\right] \]

여기서 $R_j(t_i) = \{k \mid V_k < t_i \le T_k,\; k \in \text{층 } j\}$ — 위험집합이 층 내부 + 좌절단 보정을 동시에 반영한다.

5 코드 예시

5.1 Step 1: 순수 Python으로 층화 편우도 계산

import numpy as np
from collections import defaultdict

# 간소화 데이터: (층, 시작, 끝, 사건, Z1)
# 층 0 = No MTX, 층 1 = MTX
data = [
    # 층 0 (No MTX)
    (0, 0, 10, 1, 1),   # AML low-risk, 사건
    (0, 0, 15, 0, 0),   # ALL, 중도절단
    (0, 0, 20, 1, 0),   # ALL, 사건
    (0, 0, 25, 1, 1),   # AML low-risk, 사건
    # 층 1 (MTX)
    (1, 0, 12, 1, 0),   # ALL, 사건
    (1, 0, 18, 0, 1),   # AML low-risk, 중도절단
    (1, 0, 30, 1, 1),   # AML low-risk, 사건
]

def stratified_log_partial_likelihood(beta, data):
    """층화 로그 편우도를 계산한다."""
    # 층별로 데이터 그룹화
    strata = defaultdict(list)
    for stratum, start, stop, event, z in data:
        strata[stratum].append((start, stop, event, z))
    
    total_ll = 0.0
    for s_id, records in strata.items():
        # 이 층의 사건 시점
        event_times = sorted(
            stop for start, stop, event, z in records if event == 1
        )
        for t in event_times:
            # 사건 발생자 공변량
            z_event = next(
                z for start, stop, event, z in records
                if stop == t and event == 1
            )
            # 위험집합 (좌절단 보정 포함)
            risk_z = [
                z for start, stop, event, z in records
                if start < t <= stop
            ]
            num = np.exp(beta * z_event)
            den = sum(np.exp(beta * z) for z in risk_z)
            total_ll += np.log(num / den)
    
    return total_ll

# 격자 탐색
betas = np.linspace(-3, 1, 200)
lls = [stratified_log_partial_likelihood(b, data) for b in betas]
beta_hat = betas[np.argmax(lls)]
print(f"beta_hat = {beta_hat:.3f}")
print(f"hazard ratio = {np.exp(beta_hat):.3f}")

5.2 Step 2: R survival 패키지 — 층화 + 좌절단

library(survival)

# 층화 Cox 모형 (MTX로 층화)
# strata() 함수로 층화 변수 지정
fit_strat <- coxph(
  Surv(time, status) ~ z1 + z2 + z3 + z4 + z5 + z6 + zp_td + strata(mtx),
  data = bmt
)
summary(fit_strat)

# 공통 beta 가정 검정: 층별 개별 모형 vs 통합 모형
fit_no_mtx <- coxph(Surv(time, status) ~ z1+z2+z3+z4+z5+z6+zp_td,
                    data = bmt, subset = (mtx == 0))
fit_mtx    <- coxph(Surv(time, status) ~ z1+z2+z3+z4+z5+z6+zp_td,
                    data = bmt, subset = (mtx == 1))

# 우도비 검정
lr_stat <- -2 * (fit_strat$loglik[2] -
                 (fit_no_mtx$loglik[2] + fit_mtx$loglik[2]))
p_val <- pchisq(lr_stat, df = 7, lower.tail = FALSE)
cat("LR chi-sq =", lr_stat, ", p =", p_val, "\n")

# 좌절단 분석: 혈소판 회복 후 환자만
# Surv(진입시점, 종료시점, 사건) 형식
bmt_recovered <- subset(bmt, platelet_recovered == 1)
fit_trunc <- coxph(
  Surv(time_platelet_recovery, time_event, status) ~
    z1 + z2 + z3 + z4 + z5 + z6 + strata(mtx),
  data = bmt_recovered
)
summary(fit_trunc)

# Python: lifelines 좌절단
import pandas as pd
from lifelines import CoxPHFitter

# 좌절단 데이터: entry_time = 혈소판 회복 시점
df = pd.DataFrame({
    'entry': [15, 20, 10, 25, 18],   # 좌절단 시점 (혈소판 회복)
    'stop':  [45, 60, 35, 80, 50],   # 관찰 종료
    'event': [1, 0, 1, 0, 1],
    'group': [1, 1, 0, 0, 1],        # AML low-risk indicator
    'mtx':   [0, 0, 1, 1, 0],        # 층화 변수
})

cph = CoxPHFitter()
cph.fit(
    df,
    duration_col='stop',
    event_col='event',
    entry_col='entry',       # 좌절단
    strata=['mtx'],          # 층화
    formula='group',
)
cph.print_summary()

6 요약 — §9.3–9.4 핵심 흐름

PH 가정 위반 변수 발견 (§9.2에서 MTX 검출)
    │
    ▼
층화 Cox 모형 (§9.3)
    ├─ 위반 변수 → 층화 변수: h₀ⱼ(t) 별도, β 공통
    ├─ 편우도 = Σ(층별 편우도)
    ├─ 공통 β 가정 검정: LR 또는 Wald
    └─ 대응표본: 각 쌍을 층으로 → Score 검정 = 층화 로그순위

좌절단 / 지연 진입 (§9.4)
    ├─ 위험집합 재정의: R(t) = {j | Vⱼ < t ≤ Tⱼ}
    ├─ 핵심 가정: X ⊥ V | Z
    ├─ 시나리오 A: 연령 척도 → 좌절단 = 입소 연령
    └─ 시나리오 B: 중간 사건 → 좌절단 = 회복 시점

결합: 층화 + 좌절단 동시 적용 (Example 9.5)
    → R_j(t) = {k | V_k < t ≤ T_k, k ∈ 층 j}

7 관련 주제

선행 지식

후속 주제

관련 개념