ワイブル分布（Weibull distribution）を使った生存時間解析

前回、パラメトリック生存時間解析を勉強してみて、ワイブル分布について学び直す必要性を痛感したのでまとめました。

necostat.hatenablog.jp

指数分布（exponential distribution）
ワイブル分布（Weibull distribution）
- 形状パラメータ
- 尺度パラメータ
おわりに
参考資料

指数分布（exponential distribution）

まずは生存時間（あるいはイベントの発生間隔）のパラメトリックモデルの基本になる指数分布から。確率密度関数（probability density function, PDF）は次の式で表されます。

$\begin{aligned} f(t) &= \lambda e^{- \lambda t} \\ \end{aligned}$

ここで、 $\lambda$ は率パラメータ（rate parameter）と呼ばれ、大きいほどイベントの発生が速くなり、イベント発生までの平均時間 $\frac{1}{\lambda}$ は短くなります。パラメータとして $\lambda$ の逆数を使った書き方をする場合もありますが、Rで使う関数（rexp( )など）ではrateパラメータを指定するので、それに合わせました。

パラメータλの大きさによって、発生するまでの平均時間が変わりますが、PDFはいずれも単調減少です。

生存関数（時点tよりも長く生存する確率）は、PDFをtから∞まで積分して、

$\begin{aligned} S(t) &= \int_t^{\infty} f(s) ds \\ &= \left[ - e^{- \lambda s} \right]_t^{\infty} \\ &= e^{- \lambda t} \end{aligned}$

となります。また、ハザード関数は、

$\begin{aligned} h(t) &= \frac{f(t)}{S(t)} \\ &= \frac{\lambda e^{- \lambda t}}{e^{- \lambda t}} \\ &= \lambda \end{aligned}$

となるため、生存時間によらずハザードが一定であることがわかります。

ワイブル分布（Weibull distribution）

形状パラメータ

ワイブル分布の基本形のPDFは次のようになります。

$\begin{aligned} f(t) &= \alpha t^{\alpha - 1} e^{-t^{\alpha}} \end{aligned}$

α（>0）はワイブル分布の形状を決めるパラメータで、形状パラメータ（shape parameter）と呼ばれます。

このPDFから生存関数は、

$\begin{aligned} S(t) &= \int_t^{\infty} f(s) ds \\ &= \left[ - e^{- s^\alpha} \right]_t^{\infty} \\ &= e^{- t^\alpha} \end{aligned}$

となるので、ハザード関数は、

$\begin{aligned} h(t) &= \frac{f(t)}{S(t)} \\ &= \frac{ \alpha t^{\alpha - 1} e^{-t^{\alpha}}}{e^{- t^\alpha}} \\ &= \alpha t^{\alpha - 1} \end{aligned}$

となります。

ハザード関数を見ると、αの値によって分布が以下のような特徴を持つことが分かります。

0<α<1のとき：
時間が経過するにつれてハザード率が低下する。術後合併症など、早期に発生するイベントをモデル化するときに適している。
α=1のとき：
指数分布と同じく、ハザード率が時間に依存しない。ランダムな要因によるイベントをモデル化するときに適している。
α>1のとき：
時間が経過するにつれてハザード率が増加する。老化が原因のイベントなど、時間が経つほど発生しやすいイベントをモデル化するときに適している。

尺度パラメータ

ワイブル分布の基本形で t を $\frac{t}{\sigma}$ に置き換えると、分布の形状はそのままで、t軸方向にσ倍に拡大することができます。最初の係数 $\frac{\alpha}{\sigma}$ は、曲線下面積が1になるようにするために付け加えられています。

$\begin{aligned} f(t) &= \left( \frac{\alpha}{\sigma} \right) \left( \frac{t}{\sigma} \right)^{\alpha - 1} e^{-\left( \frac{t}{\sigma} \right)^\alpha} \end{aligned}$

このパラメータσ（>0）は尺度パラメータ（scale parameter）と呼ばれ、分布の平均や分散に影響します。例えば、発生時間の期待値は、

$\begin{aligned} E[T] &= \sigma \Gamma \left( 1 + \frac{1}{\alpha} \right) \end{aligned}$

となり、αとσの両方で決定されます（ここで $\Gamma$ はガンマ関数であり、階乗を実数に拡張したもの）。

尺度パラメータを追加したワイブル分布の生存関数は、

$\begin{aligned} S(t) &= e^{- \left( \frac{t}{\sigma} \right)^\alpha} \end{aligned}$

ハザード関数は、

$\begin{aligned} h(t) &= \frac{f(t)}{S(t)} \\ &= \frac{\left( \frac{\alpha}{\sigma} \right) \left( \frac{t}{\sigma} \right)^{\alpha - 1} e^{-\left( \frac{t}{\sigma} \right)^\alpha}}{e^{- \left( \frac{t}{\sigma} \right)^\alpha}} \\ &= \left( \frac{\alpha}{\sigma} \right) \left( \frac{t}{\sigma} \right)^{\alpha - 1} \end{aligned}$

となります。

rweibull( )などRで使う関数では、αを引数shapeで、σを引数scaleで指定します。

おわりに

ワイブル分布の特殊形としての指数分布は、最初に紹介した指数分布の書き方とパラメータが逆数（ $\lambda = \frac{1}{\sigma}$ ）になっているので、少し混乱しやすいかもしれません。ポアソン分布との関係を考えると、指数分布のパラメータは逆数でない方がいいですね。
ワイブル分布は一般的に指数分布が持つような「無記憶性（memorylessness）」を持ちません。
最強寒波のため猫トイレの掃除は1週間延期です。