- 事情があって、リスク解析に暴露されることになっている。
- なんだかんだでJapanSCOREなるものが出てきたので、それに付随するものをまとめてみる。
- こちら
- そもそもJapan SCOREとは、心臓外科における、リスク解析モデルで、日本中の心臓血管に関する手術のデータベースを用いて、まだ手術を受けていない患者さんの死亡率や、術後30日以内の合併症発症率などを割り出す手法であり、インフォームドコンセントの場でも使われることがある、信頼度の高いものだ。
- この文献についてはこちら。
-
https://www.ahajournals.org/doi/10.1161/CIRCULATIONAHA.107.756684
- なんだか小難しく見えるが、要素を分解しつつ、噛み砕きながらやってみよう。
- 関連する単語をいくつかpick upしていく。
- まずは、単純なものとして、二変数解析を取り上げる。
- Bivariate analysis
- wiki
- 基本的な意味あい
- ２つの変数についての解析で、それらの関係を学ぶためのもの
- 仮説生成に役立つ。例えば、片方の変数について、他方についての情報から推定することが出来る。
- 変数同士が依存関係にあるとき
- graphical methods
- 2変量が両方とも連続的であるなら、散布図
- 両方とも離散的なら、モザイク図
- 片方だけ連続的、他方が離散的なら、箱（ひげ）図
- 大事にしていることが何かを振り返る。
- まず、目的のものと、それに付随する何かの塊(data)がある。
- そして、目的のものと、付随するものとの関係を知りたいとする。
- その時に、どんな関係かを記述する手法は何か、という問題。これが回帰の本質。
- そして、その際に気にすることとしては、それぞれのものが、どんな性質を持っているかに応じて、使い分けることになる。
- 観点としては、離散的か連続的か、順序づけがあるか。それぞれのものについて考慮する。
- ここまで来たので、脚力がついただろうので、JapanSCOREについて見ていくとする。
- multiple logistic regression
- どの変数を選ぶかは、二変数解析により行う
- 実際、JapanSCOREを用いて、リスク解析をして、それをもとにデータを取るか否か（つまり手術をするか否か）を決めるので、それによるバイアスがどの程度なのか気になったので、ロジスティックモデルを使ってそれを試すおもちゃを作った。
- 普通に何も考えないで、ロジスティック回帰をおこなった場合と、時々ロジスティック回帰での係数を求めて、それによってデータを取る・捨てる操作をするのを繰り返すとどうなるかを見てみる。
import numpy as np
import math
import random
import matplotlib.pyplot as plt
np.random.seed(314)
beta_0 = -4
beta_1 = 1.5
beta_2 = 2.5
num_of_sample = 1000
def f(beta_0, beta_1, beta_2, x_1, x_2):
y_prob = 1/(1 + np.exp(-(beta_0 + beta_1 * x_1 + beta_2 * x_2)))
return y_prob
x_1 = np.random.random(num_of_sample)
x_2 = np.random.random(num_of_sample)
y_prob = f(beta_0, beta_1, beta_2, x_1, x_2)
y_true = np.random.binomial(1, y_prob)
#print(y_true)
fig1,ax1 = plt.subplots(figsize = (8,8))
from sklearn.linear_model import LogisticRegression
lr = LogisticRegression(penalty =None)
X1 = np.array([[x_1[i], x_2[i]] for i in range(len(y_true))])
Y1 = np.array(y_true)
lr.fit(X1, Y1)
print(lr.intercept_[0], *lr.coef_[0])
print(1)
for i in range(num_of_sample):
if y_true[i] == 0:
ax1.scatter(x_1[i], x_2[i], color = "blue", marker = "x")
else:
ax1.scatter(x_1[i], x_2[i], color = "red", marker = "o")
fig1.suptitle("control")
plt.savefig("normal_logistic_model")
plt.show()
#samplingの仕方を変える
beta_0 = -4
beta_1 = 1.5
beta_2 = 2.5
num_of_sample = 1000
num_of_annual_surgeon = 50
cutoff = 0.1
fig2,ax2 = plt.subplots(figsize = (8,8))
def f(beta_0, beta_1, beta_2, x_1, x_2):
y_prob = 1/(1 + np.exp(-(beta_0 + beta_1 * x_1 + beta_2 * x_2)))
return y_prob
x_1_all =
x_2_all =
y_all = []
beta_0_estimated = 0
beta_1_estimated = 0
beta_2_estimated = 0
for i in range(100):
x_1 = np.random.random(50)
x_2 = np.random.random(50)
y_prob = f(beta_0, beta_1, beta_2, x_1, x_2)
y_estimate = f(beta_0_estimated, beta_1_estimated, beta_2_estimated, x_1, x_2)
y_do_surgeon = np.array([y_prob[i] for i in range(50) if y_estimate[i] > cutoff])
y_true = np.random.binomial(1, y_do_surgeon)
for i in range(len(y_true)):
x_1_all.append(x_1[i])
x_2_all.append(x_2[i])
y_all.append(y_true[i])
lr2 = LogisticRegression()
X2 = np.array([[x_1_all[i], x_2_all[i]] for i in range(len(y_all))])
Y2 = np.array(y_all)
lr2.fit(X2, Y2)
beta_0_estimated = lr2.intercept_[0]
beta_1_estimated = lr2.coef_[0,0]
beta_2_estimated = lr2.coef_[0,1]
"""
if y_true[i] == 0:
ax2.scatter(x_1[i], x_2[i], color = "blue", marker = "x")
else:
ax2.scatter(x_1[i], x_2[i], color = "red", marker = "o")
"""
print(beta_0_estimated, beta_1_estimated, beta_2_estimated)
for i in range(len(y_all)):
if y_all[i] == 0:
ax2.scatter(x_1_all[i], x_2_all[i], color = "blue", marker = "x")
else:
ax2.scatter(x_1_all[i], x_2_all[i], color = "red", marker = "o")
fig2.suptitle("cutoff = 0.3")
plt.show()
plt.savefig("adaptive_logistic_model_cutoff=0.3")
かなり精度は悪いように思う。
(↑cutoff 0.3と書いているけれど、コード中では0.1となっています。すみません)
この辺りから、臨床の際の、選択の手法の多様性などが、いかに影響を与えているかが容易に想像できる。
内科から外科に回すか否か。TAVIをするか否か。CABGを選ぶか否か。
そもそも病院に来たか来ないか。
医療データのサンプル収集ほどこんがらがったものは無いだろう。