アルゴリズム (phm cox)

i = 1, 2, ..., nで、\[t_i\,\!\]をi番目の観測値は、p 共変量\[Z_j(j=1,2,\ldots ,p)\] のベクターデータを持ち、これに対する故障時間または打ち切り時間にします。故障と打ち切りの過程は独立しているものとします。ハザード関数\[\lambda (z,t)\,\!\]は、共変量 z を持つ個々のデータが時間tで故障する確率で、時間tは、個々が生存する時間として与えられます。Con比例ハザードモデルは、次のような形式になっています。

\[\lambda (z,t)=\lambda _0(t)\exp (z^{T}\beta +\omega )\,\!\]

ここで\[\lambda _0\,\!\]は、ハザード関数のベースラインで時間関数ではなく、 \[\beta \,\!\]は、不明なパラメータのベクターデータで、\[\omega\,\!\]は分かっているオフセット値です。

時間\[t_{(i)}\,\!\]における個別の故障\[d_i\,\!\]のように、\[n_d < n\,\!\]で与えられる故障時間は、 明確な故障時間 \[t_{(1)} < t_{(2)} < ?< t_{(nd)}\]と結びついており、 \[\beta\]に対する周辺尤度は、次式で近似されます。

\[L=\prod_{i=1}^{n_d}\frac{\exp (s_i^{T}\beta +\omega _i)}{[\sum_{l\in R(t_{(1)})}\exp (z_i^{T}\beta +\omega _i)]^{d_{i}}}\]

(1)

ここで、\[s_i\,\!\]は、時間\[t_{(i)}\,\!\] における観測した個々の故障の共変量の合計であり、 は、\[R(t_{(i)})\,\!\] より前のリスクにおける個々の故障です。これは、時間\[t_{(i)}\,\!\]以上に生存した個々のデータに加えて、時間\[t_{(i)}\] での故障または打ち切りのデータすべてとなります。 \[\beta\,\!\]のMLE(最大尤度見積り)は、\[\hat \beta\,\!\]で与えられ、Newton-Raphson反復法を使って(1)を最大化することで取得されます。この反復法は、段階的に行われ、下記の(2)と(3)で与えられる(1)の一階および二階微分を利用します。

\[U_j(\beta )=\frac{\partial Ln(L)}{\partial \beta _j}=\sum_{i=1}^{n_d}[s_{ji}-d_i\alpha _{ji}(\beta )]=0\]

(2)

j = 1, 2,..., p, ここで\[s_{ji}\,\!\]は、ベクターデータ\[s_i\,\!\]の j番目の要素です。

\[\alpha _{ji}(\beta )=\frac{\sum_{l\in R(t_{(1)})}z_{jl}\exp (z_l^{T}\beta +\omega _l)}{\sum_{l\in R(t_{(1)})}\exp (z_l^{T}\beta +\omega _l)}\]

同様に、

\[I_{hj}(\beta )=-\frac{\partial ^2Ln(L)}{\partial \beta _h\partial \beta _j}=\sum_{i=1}^{n_d}d_i\gamma _{hji}\]

(3)

ここで \[\gamma _{hji}=\frac{\sum_{l\in R(t_{(1)})}z_{hl}z_{jl}\exp (z_l^{T}\beta +\omega _l)}{\sum_{l\in R(t_{(1)})}\exp (z_l^{T}\beta +\omega _l)}-\alpha _{hi}(\beta )\alpha _{ji}(\beta )\] h, j = 1, ..., p.p.

\[U_j(\beta )\,\!\] は、スコアベクターの j 番目の成分で、\[I_{hi}(\beta )\,\!\] は、観測情報行列 \[I(\beta )\,\!\]の(h, j)要素です。この行列の逆行列\[I(\beta )^{-1}=I_{hi}(\beta )^{-1}\,\!\]は、\[\beta\,\!\] の分散－共分散行列を与えます。

共変量または共変量の線形の組合せは、時間と共に単調に増加または減少しており、1つ以上の\[\beta _j^{\prime }s\]は無限大となります。

もし \[\lambda _0(t)\,\!\]は\[\nu\,\!\] の層でさまざまに変化すると、k番目の層にあるデータの数は\[n_k\,\!\]、 k = 1, ... , \[\nu\,\!\])で、\[n=\sum_{k=1}^\nu n_k\]を持ち、\[\hat \beta\,\!\]を取得するために(1)を最大化するのではなく、次の周辺尤度を最大化します。

\[L=\prod_{k=1}^\nu L_k\]

(4)

ここで\[L_k\,\!\]は、(1)で簡単なサンプルとして扱われるk番目の層にある\[n_k\,\!\]観測値に対する尤度への寄与となります。 共変量係数が層にまたがって一定であると結論付けするとき、異なるベースラインハザード関数\[\lambda _0(t)\,\!\]があります。

故障時間\[t_{(i)}\,\!\]と関連しているベースライン生存関数は次のように見積もられます。

\[exp(-\hat H(t_{(i)}))\] ,

ここで \[\hat H(t_{(i)})=\sum_{t(j)\leq t(i)}(\frac{d_i}{\sum_{l\in R(t_{(j)})}\exp (z_l^T\hat \beta +\omega _l)})\]

そして、\[d_i\,\!\] は、時間\[t_{(i)}\,\!\]における故障の数です。 I番目の観測値の残差は次式で計算されます。

\[r(t_l)=\hat H(t_l)\exp (-z_l^T\hat \beta +\omega _l)\]

ここで\[\hat H(t_l)=\hat H(t_{(i)}),t_{(i)}\leq t_l<t_{(i+1)}\]

逸脱は、\[-2^*\,\!\](logarithm of marginal likelihood)と定義されます。個々の共変量が十分であるかをテストする2つの方法があります。: ネストしたモデルの共変量間の差は、適切な\[\chi ^2\,\!\]の分布で比較されます。または、パラメータ推定の正規性がz検定を形作るために使われます。推定値を標準誤差で除算するか、帰無仮説下のモデルに対するスコア関数がz検定を形作るために使われます。