Cox比例ハザードモデルの解説｜生物統計家への道

記事の目次

はじめに

生存時間解析は、製薬業界の臨床試験や観察研究において欠かせない統計手法です。特に、治療効果や予後因子を評価する際に、打ち切りデータ（censored data）を適切に扱える点が大きな強みです。
その中でもCox比例ハザードモデル（Cox Proportional Hazards Model）は、共変量を考慮しつつイベント発生リスクを推定できる柔軟な回帰モデルとして広く用いられています。
本稿では、数式による導入から、比例ハザード性の仮定、製薬業界での典型的な利用場面、そして図を用いた直感的理解までを網羅します。

Cox比例ハザードモデルとは？

Cox比例ハザードモデルは、時間 t におけるハザード関数 h(t|X) を次のように表します。

\[h(t|X) = h_0(t) \exp(\beta_1 X_1 + \beta_2 X_2 + \dots + \beta_p X_p)\]

$h(t|X)$：共変量 X を持つ被験者の時点 t におけるハザード（単位時間あたりのイベント発生率）
$h_0(t)$：ベースラインハザード（共変量がゼロのときのハザード関数）
$X_j$：共変量（例：治療群、年齢、性別、バイオマーカー値など）
$\beta_j$：回帰係数$（\exp(\beta_j)$ がハザード比）

セミパラメトリックモデル：$h_0(t)$ の形は仮定せず、$\beta$ のみを推定するため、柔軟性が高い。
ハザード比（Hazard Ratio, HR）：$\exp(\beta_j)$ は、他の変数を固定したときの $X_j$ の単位変化に伴う相対的リスク。

ここで、2群のハザード関数(共変量はそれぞれ$(X_{1},…,X_{p}) , (Z_{1},…,Z_{p})$)とすると、ハザード比は下記のようになる。

$$
\begin{eqnarray}
\frac{h(t|X)}{h(t|Z)} &=& \frac{exp(\beta_{1}X_{1} + \dots + \beta_{p}X_{p})}{exp(\beta_{1}Z_{1} + \dots + \beta_{p}Z_{p})}\\
&=& exp(\beta_{1}(X_{1} – Z_{1}) + \dots + \beta_{p}(X_{p} – Z_{p}))
\end{eqnarray}
$$

上記式からも分かるように、時点tに関係なく一定の値をとるため、比例ハザードと呼ばれています。

Cox比例ハザードモデルの仮説検定、信頼区間

ここで、$X_{1}$を治療変数として、A治療群：$X_{1}=1$、B治療群：$X_{1}=0$とすると、$\beta_{1}$が治療Aの治療Bに対する治療効果を表します。なぜなら、他の共変量$(X_{2},…,X_{p})$が同じであるA治療群とB治療群のハザード比は時間に関係なく一定であり、下記のようになる。

\[\frac{h(A治療群,X_{2}, … ,X_{p})}{h(B治療群,X_{2}, … ,X_{p})} = exp(\hat{\beta})\]

A治療群とB治療群を比較するために、以下のような帰無仮説、対立仮説を立てる。

\[H_{0}: \beta_{1} = 0\]

\[H_{1}: \beta_{1} \neq 0\]

この時の検定は、\(H_{0})の下でZ値が近似的に正規分布に従うことを利用する。

\[Z = \frac{\hat{\beta_{1}}}{SE(\hat{\beta_{1}})} \mathcal {N}(0,1)\]

これより、95%信頼区間は以下のようになる。

\[exp(\hat{\beta_{1}} \pm 1.96SE(\hat{\beta_{1}}))\]

比例ハザード性の仮定

Cox比例ハザードモデルで重要なのは、2群間のハザードが時間によらず一定である、比例ハザード性が成り立つことです。

\[\frac{h(t|X_a)}{h(t|X_b)} = \text{一定} \quad (\forall t)\]

例えば、２群間でカプランマイヤー曲線が交差したり、交絡因子が調整されていない場合は比例ハザード性が成り立たない場合があります。そのため、下記のように比例ハザード性が成り立っているか検証することが大切です。

検証方法

Schoenfeld残差の検定・プロット
log(-log S(t)) プロットの平行性確認

Rによる実践例

#パッケージ読み込み
install.packages(“survival”)
install.packages(“survminer”) # 図示用
install.packages(“gtsummary”) # 表形式出力用

library(survival)
library(survminer)
library(gtsummary)

＃データ確認
head(lung)

# Survオブジェクト作成
# status: 1=打ち切り, 2=死亡 → 死亡をイベント(=1)に変換
lung$status <- ifelse(lung$status == 2, 1, 0)

#単変量モデル（治療群はないので性別で例示）
cox_uni <- coxph(Surv(time, status) ~ sex, data = lung)
summary(cox_uni)

# 多変量モデル
cox_multi <- coxph(Surv(time, status) ~ age + sex + ph.ecog, data = lung)
summary(cox_multi)

# ハザード比を表形式で出力
tbl_regression(cox_multi, exponentiate = TRUE)

# カプラン–マイヤー曲線（性別）
fit <- survfit(Surv(time, status) ~ sex, data = lung)
ggsurvplot(fit, data = lung, pval = TRUE, risk.table = TRUE)

# 比例ハザード性の検証
cox.zph(cox_multi)
plot(cox.zph(cox_multi))

変数	coef	exp(coef)	se(coef)	z	p値
age	0.019	1.019	0.009	2.11	0.035
sex（男性=1）	-0.531	0.588	0.167	-3.18	0.001
ph.ecog	0.454	1.575	0.114	3.98	<0.001

年齢：10歳増加でHRは約1.02（リスク増加）
性別：男性のHRは女性の約0.59（リスク低下）
ECOG PS：1ポイント悪化でHRは約1.58（リスク増加）

実務上の注意点

1. 比例ハザード性の検証は必須

Coxモデルの根幹仮定である比例ハザード性（PH仮定）が破れると、推定されたハザード比は時間全体の平均的な効果を表すに過ぎず、臨床的な意味づけが難しくなります。

検証方法
- Schoenfeld残差：時間との相関を検定（cox.zph()関数など）
- log(-log S(t))プロット：群間の曲線が平行であればPH仮定はおおむね成立
PH仮定が破れた場合の対応
- 時間依存共変量モデル：例）coxph(Surv(time, status) ~ trt + tt(trt), tt=function(x,t,...) x*(t>6))
- 層別Coxモデル：PH仮定が怪しい変数を層別化因子として扱い、ベースラインハザードを群ごとに別推定
製薬業界での注意

承認申請資料（CTD）では、主要解析のPH仮定検証結果を添付することが望ましい
PH仮定が破れている場合は、代替解析（RMSTなど）を併記することが多い

2. 共変量の選択は臨床的妥当性を優先

モデルに含める変数は、統計的有意性だけでなく、臨床的・薬理学的な妥当性を重視します。

選択基準
- 事前にプロトコルやSAP（Statistical Analysis Plan）で定義
- 臨床的に重要とされる既知の予後因子（例：年齢、性別、ECOG PS）
- 試験デザイン上の層別化因子
過剰な変数投入のリスク
- イベント数に対して変数が多すぎると推定が不安定（経験則：10イベント/変数以上）
- 多重共線性による係数の不安定化
製薬業界での注意