はじめに
Cox比例ハザードモデル(Cox回帰とも呼ばれます)は、「何がイベントの起きやすさに影響しているか」を調べられる便利な手法です。医療をはじめ、解約予測や故障分析など、幅広い分野で使われています。ただ、きちんと理解しようとすると難しい手法です。
そこでこの記事では、細かい理論はいったん脇に置いて、まずCox比例ハザードモデルの「結果が読めるようになる」ことを目標にします。統計の知識はほとんどないという方でも読めるように、数式の導出には踏み込まず、「この手法で何が分かるのか」「出てきた結果をどう解釈するのか」に絞って説明していきます。
Cox比例ハザードモデルとは?
ひとことで言うと、「あるイベントの起きやすさに、どの要因がどれくらい影響しているか」を調べる手法です。生存時間解析と呼ばれる分野の手法の一つになります。
分析すると、各要因について「ハザード比」という数字が得られます。ハザード比とは、その要因が「イベントの起きやすさ」を何倍にするかを表す値です。たとえば「新薬を使うと死亡リスクが0.23倍になる」といった形で結果が出てきます。この記事のゴールは、このハザード比という数字を正しく読めるようになることです。
生存時間解析とは?
生存時間解析は、「あるイベントが発生するまでの時間」を分析することを目的とした分野です。名前に「生存」とありますが、生死に限らず、2択で表現できるイベントなら何にでも使えます。たとえば、生存/死亡、発症しない/発症、解約しない/解約、購入しない/購入といったものです。
「イベントが起きるまでどれくらいかかるか」「何が起きやすさに影響するか」を扱いますが、その大きな特徴は、次に説明する「打ち切り」を考慮できる点にあります。
打ち切りとは?
打ち切り(censoring)とは、「観察期間中にイベントが起きなかった(または追えなくなった)」状態のことです。
たとえば「新薬の効果を10年間追跡する」という研究を考えてみます。ある人は3年目に亡くなったとすれば、これはイベント発生(time=3)です。一方、ある人は10年経っても生存していたとすると、10年時点で観察が終わってしまうため、その後いつ亡くなるかは分かりません。これが打ち切りです。また別の人が5年目に引っ越して連絡が取れなくなった場合も、そこから先は追えないので打ち切りになります。
ここで大事なのは、打ち切りは「情報がない」わけではなく、「少なくともその時点まではイベントが起きていなかった」という情報を持っているということです。この打ち切りを捨てずに活かせることが、生存時間解析の強みです。もし単純に「イベントが起きた人だけ」で平均を取ってしまうと、打ち切りの人の情報が無駄になり、結果が歪んでしまいます。
実際に走らせてみる
ここからは、架空のダミーデータ(患者60人分)でCox比例ハザードモデルを実行した結果を見てみます。
使うデータは次のような形式で、イベントは死亡とします。
| 列名 | 意味 |
|---|---|
| time | 観察期間(年)。イベント発生 or 打ち切りまでの時間 |
| age | 年齢(共変量:影響を調べたい要因) |
| treat | 治療の種類(0=従来治療, 1=新薬) |
| event | 1=死亡(イベント発生), 0=打ち切り |
このデータをdummy.csvとして保存しておき、PythonのlifelinesというライブラリでCox比例ハザードモデルにかけます。なお、dummy.csvの中身はこの記事の最後に付録として載せてあるので、手元で試したい場合はコピーして使ってください。
コードは以下のようにとてもシンプルです。
import pandas as pd
from lifelines import CoxPHFitter
# データの読み込み
df = pd.read_csv("dummy.csv")
# Cox比例ハザードモデルの実行
cph = CoxPHFitter()
cph.fit(df, duration_col="time", event_col="event")
# 結果の表示(decimals=4 で小数点以下4桁まで表示)
cph.print_summary(decimals=4)
duration_colには時間の列(time)を、event_colにはイベント発生かどうかの列(event)を指定します。共変量(age, treat)は、指定した2列以外が自動的に使われます。また、decimals=4はp値などを小数点以下4桁まで表示するための指定で、これを付けないと<0.005のように丸められてしまいます。
実行すると、次のような結果が出力されます。
モデルの基本情報
| 項目 | 値 |
|---|---|
| model | lifelines.CoxPHFitter |
| duration col | ‘time’ |
| event col | ‘event’ |
| baseline estimation | breslow |
| number of observations | 60 |
| number of events observed | 28 |
| partial log-likelihood | -88.4716 |
各共変量の推定結果
| covariate | coef | exp(coef) | se(coef) | coef lower 95% | coef upper 95% | exp(coef) lower 95% | exp(coef) upper 95% | cmp to | z | p | -log2(p) |
|---|---|---|---|---|---|---|---|---|---|---|---|
| age | 0.0808 | 1.0841 | 0.0223 | 0.0371 | 0.1244 | 1.0378 | 1.1325 | 0.0000 | 3.6288 | 0.0003 | 11.7778 |
| treat | -1.4542 | 0.2336 | 0.4453 | -2.3270 | -0.5814 | 0.0976 | 0.5591 | 0.0000 | -3.2657 | 0.0011 | 9.8388 |
モデル全体の評価
| 項目 | 値 |
|---|---|
| Concordance | 0.7602 |
| Partial AIC | 180.9432 |
| log-likelihood ratio test | 28.5998 on 2 df |
| -log2(p) of ll-ratio test | 20.6304 |
数字がたくさん並んでいて戸惑うかもしれませんが、次の章でこの出力の「どこを見ればいいか」を順番に説明していきます。
結果の読み方
たくさんの数字が並んでいますが、まず見るべきは「ハザード比(exp(coef))」の列です。ここだけ読めれば、結果の大部分は理解できます。
① ハザード比(exp(coef)):影響の大きさと向き
ハザード比は「その要因がイベント(死亡)の起きやすさを何倍にするか」を表します。読み方のルールはシンプルで、ハザード比が1より大きければ、その要因が大きいほどイベントが起きやすい(リスクが上がる)ことを意味します。逆に1より小さければ、その要因が大きいほどイベントが起きにくい(リスクが下がる=保護的)ということです。ちょうど1のときは、その要因はイベントの起きやすさに影響しないと読めます。
今回の結果を見てみます。age(年齢)は1.0841なので、年齢が1歳上がるごとに、死亡のリスクが約1.08倍になります。1歳では小さく見えますが、10歳差なら1.0841の10乗で約2.2倍になり、積み重なると大きな差です。一方、treat(新薬)は0.2336なので、新薬を使うと死亡リスクが従来治療の約0.23倍にまで下がります。つまり新薬は死亡リスクを大きく下げており、保護的に働いていると解釈できます。
② p値:その影響は「偶然ではなさそう」か
p値は、その影響が「偶然ではなさそうか」を判断するための数字です。p値が小さければ(一般に0.05未満)、「偶然とは考えにくい=統計的に意味のある影響」とみなします。逆にp値が大きければ、「偶然の範囲かもしれず、効果があるとは言い切れない」と判断します。ざっくり「小さいほど、その影響は信頼できる」と覚えておけば十分です。
今回はage(p=0.0003)もtreat(p=0.0011)もどちらも0.05未満なので、両方とも統計的に意味のある影響だと解釈できます。
③ 95%信頼区間:ハザード比の「ぶれ幅」
ハザード比は1つの推定値ですが、データが変われば多少ぶれます。その「もっともらしい範囲」を示すのが95%信頼区間です。
見るポイントは、信頼区間が1をまたいでいるかどうかです。1をまたいでいなければ(たとえばtreatの0.0976〜0.5591は全て1未満)「効果あり」、下限0.8・上限1.3のように1をまたいでいれば「リスクを上げるとも下げるとも言い切れない」と判断します。今回はageもtreatも1をまたいでいないので「効果あり」と読めます。
勘のいい方は「p値と同じことを言っているのでは?」と気づくかもしれません。そのとおりで、「信頼区間が1をまたがない」ことと「p値が0.05未満」は必ず一致します(どちらも同じ計算から出てくるためです)。ではなぜ両方あるのかというと、信頼区間は有意かどうかに加えて「効果がどのくらいの大きさになりそうか」まで分かるからです。たとえばtreatの0.0976〜0.5591からは、「新薬は死亡リスクを1/10くらいまで下げるかもしれないし、半分程度かもしれない」という幅が読み取れます。
まとめ:最低限これだけ押さえれば読める
Cox比例ハザードモデルの結果は、まずハザード比(exp(coef))を見て、1より大きいか小さいかで影響の向きを判断します。そのうえで、p値(0.05未満か)と信頼区間(1をまたがないか)を見て、その影響が信頼できるかどうかを判断します。この3つを押さえておけば、Cox比例ハザードモデルの結果はひとまず読みこなせます。
おわりに
今回は「結果の読み方」に絞って、Cox比例ハザードモデルを使えるところまでを見てきました。一方で、「そもそもハザードとは何か」「なぜ比を取るのか」「打ち切りをどう計算に入れているのか」といった裏側の仕組みには触れていません。このあたりは、機会があれば別の記事で扱いたいと思います。
付録:使用したダミーデータ
この記事で使ったdummy.csvの中身は以下のとおりです(患者60人分)。手元で試す場合は、これをコピーしてdummy.csvという名前で保存してください。
time,age,treat,event
0.7,78,0,0
10.7,66,0,1
11.6,68,0,1
1.1,76,0,1
20.0,65,1,0
6.7,72,0,0
7.6,74,0,1
0.9,52,0,1
8.5,46,1,0
20.0,55,1,0
20.0,54,1,0
4.6,75,0,1
20.0,76,1,0
20.0,45,0,0
11.5,62,0,1
3.2,73,0,1
20.0,49,1,0
17.8,72,1,1
20.0,49,1,0
2.1,61,0,1
6.1,73,1,1
10.0,55,0,0
20.0,56,1,0
1.9,54,1,1
14.6,70,0,1
20.0,53,1,0
4.4,79,1,0
20.0,60,1,0
13.8,61,1,0
20.0,62,1,0
0.3,65,0,1
20.0,64,1,0
8.3,62,0,1
13.2,79,0,1
1.3,73,0,1
20.0,72,1,0
0.5,69,1,1
20.0,66,1,0
6.3,56,1,1
3.2,79,1,1
4.6,61,1,0
20.0,52,0,0
15.9,74,1,1
12.1,50,1,0
8.8,75,0,1
6.7,66,0,1
20.0,49,1,0
20.0,46,0,0
9.1,60,1,0
20.0,46,0,0
4.0,49,0,0
16.4,63,0,0
0.6,78,0,1
11.0,61,1,0
0.1,73,0,1
10.1,77,0,1
20.0,73,1,0
10.2,67,1,1
13.6,60,0,1
8.4,62,1,0

