読者です 読者をやめる 読者になる 読者になる

tak0kadaの何でもノート

発声練習、生存確認用。

医学関連は 医学ノート

生存時間解析について -その2(定式化、その他)- 「Think Stats」第2章

統計 thinkstats python numpy

生存時間解析の定式化

生存時間解析は興味のあるイベントが起こるまで観察し、その時間について考察を加えるものである。

  • \(S\)、\(F\)、\(s\)、\(f\)

生存時間解析でまず最初に観察されるのが、ある時間\(t\)以上生存する確率を表す生存関数\(S(t)\)と、\(t\)より先に死ぬ確率を表す寿命分布関数\(F(t)\)である。 \[ \begin{align} S(t) &= &Pr(T > t)\\ F(t) &= &Pr(T < t) \end{align} \]

これらは \[ S(0) = 1, S(\infty) = 0, F(0) = 0, F(\infty) = 1, S(t) + F(t) =1 \] を満たす。

さらにこれらを時間微分した確率密度関数 \[ \begin{align} f(t) &= &\dot{F}(t)\\ s(t) &= &\dot{S}(t) \end{align} \] と定義すると、\(s\)、\(f\)は \[ s(t) + f(t) = 0 \] を満たす。

※ 今回の問題では、\(F\)がCDF、\(f\)がPDFに対応している。

  • \(\lambda\)、\(\Lambda\)

ある時間tでの生存者が死ぬ確率をハザード関数\(\lambda(t)\)とすると、 \[ \begin{align} \lambda(t) &= &\dfrac{Pr(t<T<t+dt | T>t)}{dt}\\ &= &\dfrac{S(t)-S(t+dt)}{dtS(t)}\\ &= &-\dfrac{\dot{S}(t)}{S(t)} \end{align} \] ハザード関数は常に正の値をとる。

これを積分したものである累積ハザード関数\(\Lambda\)は \[ \begin{align} \Lambda(t) &= &\int_{0}^{t}\lambda(u)du\\ &= &-\log{S(t)} \end{align} \] である。

条件付きPMF

ある時間\(t\)に生存していて、\(t+\Delta t\)までに死亡する確率を求める。これは\(\lambda\Delta t\)に相当する。

fig = plt.figure()
ax = fig.add_subplot(111)

ax.set_ylabel('conditional pmf')
ax.set_xlabel('days')

for i, arr in enumerate([lev, lev5fu, obs]):
    pmf_50, time_bins = np.histogram(arr, bins=70, range=(0,3500), density=True)
    S = 1 - np.cumsum(pmf_50*50)
    #delete 0 at the end of S (eg. [1, 0.9, ..., 0.001, 0, 0])
    S = np.delete(S, np.where(S==0))
    S = np.hstack((S,[0]))
    t = ((time_bins[:-1]+time_bins[1:])/2)[:len(S)-1]

    #condional pmf = hdt = -S'/Sdt = -dS/S
    dS = S[1:] - S[:-1]
    hdt = -dS/S[:-1]

    ax.plot(t, hdt, color=('b','g','r')[i])
    ax.legend(['Lev', 'Lev+5FU', 'Obs'], loc='best')

plt.show()

f:id:kuyata:20140430202408p:plain

時間が経過するに連れて分母の\(S\)が小さくなるので、条件付き確率は1に近づく。

相対危険度

現在のがん治療ではなかなか完全には治らないので「n年生存率がx%」のような言い方をする。ここでは5年生存率を考え、それぞれの治療について5年生存する相対危険度を求めてみる。

a = len(np.extract(lev>1825, lev))/len(lev)
b = len(np.extract(lev5fu>1825, lev5fu))/len(lev5fu)
c = len(np.extract(obs>1825, obs))/len(obs)

print("Lev+5FU/Obs relative risk="+str(b/c))
print("Lev+5FU/Lev relative risk="+str(b/a))

とすると

Lev+5FU/Obs relative risk=1.21104029605
Lev+5FU/Lev relative risk=1.1627487163

となって、Lev+5FUの方が5年生存する相対危険度が大きい(つまり効果がある)らしいと思われる。