兩組存活曲線比較 | 存活資料分析

兩組 Kaplan–Meier 曲線可以先用圖形描述，但曲線之間看起來有距離，不代表母體存活情形一定不同。Log-rank test 會在每一個事件時間比較兩組實際事件數與虛無假設下的期望事件數，再把整段追蹤期間的差異合併成一個檢定統計量。

先從虛無假設出發

令 S₁(t)、S₂(t) 分別代表兩組母體的存活函數。虛無假設是兩組在所有時間點的存活曲線相同；也就是在已經存活到任一時間 t 的條件下，兩組接下來發生事件的風險沒有系統性差異。

虛無假設

H_0:S_1(t)=S_2(t)\quad\text{for every }t

對立假設

H_1:S_1(t)\ne S_2(t)\quad\text{for at least one }t

若 H₀ 成立，在同一個事件時間內，每一位仍在風險集中的人都應面對相同的事件機會。因此，該時間的總事件數應依兩組風險集人數所占比例分配，而不是依研究最初的樣本數分配。這正是 log-rank test 的出發點。

例子：兩種骨髓移植的存活資料

以下比較成人接受 autologous bone marrow transplant 與 allogeneic bone marrow transplant 後的存活曲線。Autologous 組有 33 人，allogeneic 組有 21 人；月份後面的「+」代表該時間發生右截尾，沒有「+」則代表觀察到事件。

補充：完整事件與截尾資料

Autologous month	Deaths or lost to follow-up	Allogeneic month	Deaths or lost to follow-up
1	3	1	1
2	2	2	1
3	1	3	1
4	1	4	1
5	1	6	1
6	1	7	1
7	1	12	1
8	2	15+	1
10	1	20+	1
12	2	21+	1
14	1	24	1
17	1	30+	1
20+	1	60+	1
27	2	85+	2
28	1	86+	1
30	2	87+	1
36	1	90+	1
38+	1	100+	1
40+	1	119+	1
45+	1	132+	1
50	3	—	—
63+	1	—	—
132+	2	—	—
Total	33	Total	21

先分別估計兩組 Kaplan–Meier 曲線

每組都依上一單元的方法，在事件時間乘上 (nᵢ−dᵢ)/nᵢ；截尾時間只會讓後續風險集變小，不會讓曲線下降。完整計算保留如下，並將教材中的英文表格重建為可搜尋、可縮放的網站表格。

補充：兩組完整 Kaplan–Meier 計算表

Autologous bone marrow transplant（n=33）

Month, tᵢ	Deaths or lost, dᵢ	Number at risk, nᵢ	Conditional survival	Cumulative survival, Ŝ(t)
1	3	33	30/33 = 0.909	0.909
2	2	30	28/30 = 0.933	0.848
3	1	28	27/28 = 0.964	0.817
4	1	27	26/27 = 0.963	0.787
5	1	26	25/26 = 0.962	0.757
6	1	25	24/25 = 0.960	0.727
7	1	24	23/24 = 0.958	0.697
8	2	23	21/23 = 0.913	0.636
10	1	21	20/21 = 0.952	0.605
12	2	20	18/20 = 0.900	0.545
14	1	18	17/18 = 0.944	0.514
17	1	17	16/17 = 0.941	0.484
20+	1	16	Censored	Unchanged: 0.484
27	2	15	13/15 = 0.867	0.420
28	1	13	12/13 = 0.923	0.388
30	2	12	10/12 = 0.833	0.323
36	1	10	9/10 = 0.900	0.291
38+	1	9	Censored	Unchanged: 0.291
40+	1	8	Censored	Unchanged: 0.291
45+	1	7	Censored	Unchanged: 0.291
50	3	6	3/6 = 0.500	0.145
63+	1	3	Censored	Unchanged: 0.145
132+	2	2	Censored	Unchanged: 0.145

Allogeneic bone marrow transplant（n=21）

Month, tᵢ	Deaths or lost, dᵢ	Number at risk, nᵢ	Conditional survival	Cumulative survival, Ŝ(t)
1	1	21	20/21 = 0.952	0.952
2	1	20	19/20 = 0.950	0.904
3	1	19	18/19 = 0.947	0.857
4	1	18	17/18 = 0.944	0.809
6	1	17	16/17 = 0.941	0.762
7	1	16	15/16 = 0.938	0.714
12	1	15	14/15 = 0.933	0.666
15+	1	14	Censored	Unchanged: 0.666
20+	1	13	Censored	Unchanged: 0.666
21+	1	12	Censored	Unchanged: 0.666
24	1	11	10/11 = 0.909	0.605
30+	1	10	Censored	Unchanged: 0.605
60+	1	9	Censored	Unchanged: 0.605
85+	2	8	Censored	Unchanged: 0.605
86+	1	6	Censored	Unchanged: 0.605
87+	1	5	Censored	Unchanged: 0.605
90+	1	4	Censored	Unchanged: 0.605
100+	1	3	Censored	Unchanged: 0.605
119+	1	2	Censored	Unchanged: 0.605
132+	1	1	Censored	Unchanged: 0.605

自體與異體骨髓移植的 Kaplan–Meier 存活曲線；階梯只在事件時間下降，短直線表示右截尾 — 兩組 Kaplan–Meier 估計曲線。Allogeneic 組的樣本存活曲線大多高於 autologous 組；短直線是右截尾標記。來源：本站依原始筆記的骨髓移植範例數據重製

圖上 allogeneic 組的估計曲線在追蹤期間大多高於 autologous 組，但圖形只能描述樣本。要判斷這個差距是否足以反對 H₀，仍須把每個事件時間的觀察事件數與期望事件數合併比較。

H₀ 下，事件數應如何分配？

在第 j 個事件時間，令 nAj、nBj 為兩組事件發生前的風險集人數，dAj、dBj 為該時間的事件數。先將兩組合併，得到總風險集 nⱼ 與總事件數 dⱼ。

合併風險集

n_j=n_{Aj}+n_{Bj}

合併事件數

d_j=d_{Aj}+d_{Bj}

A 組期望事件數

e_{Aj}=E(D_{Aj}\mid H_0)=d_j\frac{n_{Aj}}{n_j}

該時間的觀察值減期望值

u_j=d_{Aj}-e_{Aj}

例如第 1 個月共有 54 人仍在風險集中，其中 autologous 組有 33 人；該月兩組合計發生 4 個事件。若 H₀ 成立，4 個事件中分配給 autologous 組的期望數就是 4×33/54=2.444。實際觀察到 3 個，因此該月 u₁=3−2.444=0.556。

第 1 個月的期望事件數

e_{A1}=4\times\frac{33}{54}=2.444

第 1 個月的 O−E

u_1=3-2.444=0.556

為什麼變異數是這個公式？

在 H₀ 下，先固定該時間共有 dⱼ 個事件，再問其中有幾個落在 nAj 位 autologous 受試者。這相當於從 nⱼ 位風險集成員中不放回抽出 dⱼ 位事件者，D_Aj 因而服從超幾何分配。超幾何分配的平均數正是上面的比例分配，而變異數包含不放回抽樣的有限母體修正。

事件配置的條件分配

D_{Aj}\mid n_{Aj},n_{Bj},d_j,H_0\sim\operatorname{Hypergeometric}(n_j,n_{Aj},d_j)

單一事件時間的變異數

v_j=\operatorname{Var}(D_{Aj}\mid H_0)=d_j\frac{n_{Aj}}{n_j}\frac{n_{Bj}}{n_j}\frac{n_j-d_j}{n_j-1}

整理後的形式

v_j=\frac{n_{Aj}n_{Bj}d_j(n_j-d_j)}{n_j^2(n_j-1)}

補充：從超幾何分配完整推到 log-rank 變異數

把第 j 個事件時間看成一個 2×2 配置問題：風險集內共有 nAj 位 A 組與 nBj 位 B 組受試者，其中固定有 dⱼ 人發生事件。H₀ 下，事件者在風險集中的所有等大小配置具有相同機會。

超幾何機率

P(D_{Aj}=x\mid H_0)=\frac{\binom{n_{Aj}}{x}\binom{n_{Bj}}{d_j-x}}{\binom{n_j}{d_j}}

超幾何平均數

E(D_{Aj}\mid H_0)=d_j\frac{n_{Aj}}{n_j}=e_{Aj}

超幾何變異數

\operatorname{Var}(D_{Aj}\mid H_0)=d_j\frac{n_{Aj}}{n_j}\left(1-\frac{n_{Aj}}{n_j}\right)\frac{n_j-d_j}{n_j-1}

用 nBj/nⱼ 取代括號

1-\frac{n_{Aj}}{n_j}=\frac{n_j-n_{Aj}}{n_j}=\frac{n_{Bj}}{n_j}

得到表格使用的變異數貢獻

v_j=\frac{n_{Aj}n_{Bj}d_j(n_j-d_j)}{n_j^2(n_j-1)}

各時間的 uⱼ 是依當時風險集計算的條件觀察值減期望值。在 H₀ 下，這些逐時條件增量的期望值都是 0；利用條件變異數逐步相加，可得到 U_L 的變異數。事件時間夠多且沒有單一時點完全支配總和時，標準化後的 U_L 會近似標準常態。

累加觀察值減期望值

U_L=\sum_j u_j=\sum_j(d_{Aj}-e_{Aj})

虛無假設下的平均數

E(U_L\mid H_0)=\sum_jE(u_j\mid H_0)=0

累加條件變異數

\operatorname{Var}(U_L\mid H_0)=\sum_jv_j

標準化後的近似分配

Z=\frac{U_L}{\sqrt{\sum_jv_j}}\overset{H_0}{\approx}N(0,1)

逐一合併所有事件時間

Log-rank test 只在至少一組發生事件的月份建立比較列；只有截尾而沒有事件的月份不產生 uⱼ，但先前的截尾者必須從後續風險集移除。下表依列出的風險集重新計算全部數值。最後一欄是對 Var(U_L) 的貢獻 vⱼ，不是標準誤本身；先將它們相加後再開根號，才得到 sd(U_L)。

Month	dA	nA	dB	nB	d total	n total	d/n	Expected A, eA	A: O−E	Variance contribution, vⱼ
1	3	33	1	21	4	54	0.074	2.444	0.556	0.897
2	2	30	1	20	3	50	0.060	1.800	0.200	0.691
3	1	28	1	19	2	47	0.043	1.191	−0.191	0.471
4	1	27	1	18	2	45	0.044	1.200	−0.200	0.469
5	1	26	0	17	1	43	0.023	0.605	0.395	0.240
6	1	25	1	17	2	42	0.048	1.190	−0.190	0.471
7	1	24	1	16	2	40	0.050	1.200	−0.200	0.469
8	2	23	0	15	2	38	0.053	1.211	0.789	0.465
10	1	21	0	15	1	36	0.028	0.583	0.417	0.243
12	2	20	1	15	3	35	0.086	1.714	0.286	0.691
14	1	18	0	14	1	32	0.031	0.563	0.438	0.246
17	1	17	0	13	1	30	0.033	0.567	0.433	0.246
24	0	15	1	11	1	26	0.038	0.577	−0.577	0.244
27	2	15	0	10	2	25	0.080	1.200	0.800	0.460
28	1	13	0	10	1	23	0.043	0.565	0.435	0.246
30	2	12	0	10	2	22	0.091	1.091	0.909	0.472
36	1	10	0	9	1	19	0.053	0.526	0.474	0.249
50	3	6	0	9	3	15	0.200	1.200	1.800	0.617
Total									U_L=6.572	Var(U_L)=7.884

由 U_L 得到 z、χ² 與 p 值

把每個事件時間的 O−E 相加得到 U_L，再用累加變異數的平方根標準化。兩組比較只有一個獨立的組間方向，所以 Z² 服從自由度 1 的卡方近似；用雙尾標準常態或右尾 χ²₁ 會得到相同 p 值。

累加結果

U_L=6.572,\qquad \operatorname{Var}(U_L)=7.884

標準差

sd(U_L)=\sqrt{7.884}=2.808

標準常態統計量

z=\frac{6.572}{2.808}=2.341

等價的卡方統計量

\chi^2=z^2=5.479,\qquad df=1

雙尾 p 值

p=2P(Z\ge2.341)\approx0.019

因為 p≈0.019<0.05，所以拒絕兩組存活曲線完全相同的 H₀。本例 U_L 為正，表示 autologous 組觀察到的事件總數高於 H₀ 下的期望值；配合 Kaplan–Meier 圖，可看出其估計存活曲線較低。檢定結果說明兩組曲線有差異，但差異大小仍應搭配特定時間的存活率、存活時間摘要或 Cox model 的 hazard ratio 與信賴區間報告。

連續性校正

教材另外將 0.5 連續性校正套用在離散的 O−E 總和上。雙尾檢定應從 |U_L| 朝 0 的方向修正，因此分子寫成 |U_L|−0.5。校正後的統計量較小，結論較保守。

連續性校正

z_{cc}=\frac{|U_L|-0.5}{sd(U_L)}

本例

z_{cc}=\frac{6.572-0.5}{2.808}=2.163

Log-rank test 的判讀界線

兩組受試者應彼此獨立；同一人的重複事件、配對或群聚資料不能直接當成獨立兩組處理。
截尾應近似非資訊性，而且每位受試者必須有明確的追蹤起點、時間與事件指標。
Log-rank test 比較整段追蹤期間的曲線，不只是某一個月份的存活率；顯著結果也不代表每一時間點都不同。
當兩組 hazard 大致維持固定比例時，普通 log-rank test 通常最有檢定力；若曲線交叉，早期與晚期的 O−E 可能互相抵消，必須同時查看曲線並考慮其他預先指定的方法。
曲線尾端只剩少數受試者時，估計很不穩定；比較圖應同時提供 number at risk，而不能只看兩條線的距離。