ロジスティック回帰分析

1
ロジスティック回帰分析
with the assistance of Mr. M. Torii
二値データの回帰分析法
2
規準変数が二値の場合の
回帰分析
•
•
•
•
「成功・失敗」を原因系の変数で予測
「発症・非発症」を原因系の変数で予測
「賛成・反対」を原因系の変数で予測
「賛成・どちらでもない・反対」を「賛成する・
賛成しない」の二値に落として原因系の
変数で予測
3
例：心疾患の発症
• 出展：丹後他(1996)ロジスティック回帰分析
原典(Truett et. al. 1967)
• n=2187, 男性
• 規準変数
– 冠状動脈性疾患の発症(12年後に発症したかどうか)
• 説明変数
–
–
–
–
–
–
–
年齢
コレステロール
血圧
相対体重
ヘモグロビン
喫煙
ECG所見
4
分析結果
説明変数
(リスクファクター)
X1 年齢
X2 コレステロール
X3 収縮期血圧
X4 相対体重
X5 ヘモグロビン
X6 喫煙
X7 ECG所見
相対体重：
喫煙：
ECG所見：
βの推定値
0.0708
0.0105
0.0166
0.0138
-0.0837
0.3610
1.0460
SE
0.0083
0.0016
0.0036
0.0051
0.0542
0.0587
0.2700
Wald 検定
8.5301
6.5625
4.6111
2.7059
-1.5443
6.1499
3.8741
体重÷年齢別身長別体重の中央値
0: never, 1: 1箱未満, 2: 1箱, 3: 1箱より多い
0: 正常, 1: 正常でない
5
普通に回帰分析しては
いけないのか
• してはいけない
– ０・１変数を連続変数で予測するというモデル
に無理がある
– y^=0.8, 1.5, -0.4のような予測値はどのように
解釈すればよいか不明
yi   0  1 xi1   2 xi 2  ei
xi1 , xi ：連続変数でも固定変
数でもよい
2
e：正規分布に従う連続
変数
i
 y：正規分布に従う連続
変数
i
6
では，どう考えるか
• 原因系変数が結果の生起確率P(Y=1)に
影響すると考えるのが自然
• P(Y=1)=a+bxはどうか？
– ダメ
– a+bxは区間[0,1]に収まらないことがある
– 0.5→0.6とするための努力と0.85 → 0.95と
するための努力には違いがある
7
では，どうするか
• そこで，生起（成功）確率を支配する実力と
いう潜在変数（心理学的連続体）があり，
それが正規分布すると仮定する
• さらに，その潜在変数が原因系の変数（説
明変数）から影響を受けることを想定する
8
成功する確率
実力
失敗する確率
原因系の変数が実力に影響する
9
実力と成功確率
成功確率：50%⇒60%
実力の増分：0.25
成功確率：85%⇒95%
実力の増分：0.60
10
正規分布のロジット近似
成功確率  P (Y  1)  緑の部分


実力
1
1
 y2
e 2 dy

2
1 e
これを「実力」につい
て解くと
P(Y  1)
log
 c  実力
1  P (Y  1)

1
 c実力
( c  1. 7 )
左辺をロジットという
．
ロジスティック回帰分
析は生起確率のロジットに回帰モデルを
想定したもの：
P(Y  1)
log
 a  bx
1  P (Y  1)
11
ロジスティック回帰モデル
一般にある現象の発生する確率（割合）ｐを、そ
の現象の生起を説明するために観測された変数
群 x  ( x1 , x2 ,  , xr ) で説明しようと考える
場合、 x  ( x1 , x2 ,  , xr ) という状態のも
とで現象が生起するという条件付き確率を p (x )
で表し、これを、
p(x)  Pr{生起 | x1 , x2 ,  , xr }  F ( x1 ,  , xr )
という関数Fを用いてモデル化する。
12
ロジスティック回帰モデル_2
つぎのFを用いてモデル化：
p(x)  Pr{生起 | x1 , x2 ,  , xr }
 F ( x1 , x2 ,  , xr )
1

1  exp  (  0  1 x1   2 x2       r xr )
ロジスティック関数：
1
exp( Z )
F (Z ) 

1  exp(  Z ) 1  exp( Z )
Z   0  1 x1   2 x2       r xr
13
ロジット（logit）
1
p ( x) 
1  exp  (  0  1 x1   2 x2       r xr )
ここから，式変形を行
うと
p ( x)
log
  0  1 x1   2 x2       r xr
1  p ( x)
上記左辺を p(x)のロジットという
14
オッズ(odds)
p ( x)
log
  0  1 x1   2 x2       r xr
1  p ( x)
ここから，式変形を行うと
p ( x)
 exp  0  1 x1   2 x2       r xr 
1  p ( x)
上記左辺を p(x)の1  p(x)に対するオッズという
15
オッズ比(odds ratio)
x1の odds比とは（偏回帰係数で
ある！？）
p ( x' )
 exp  0  1 ( x1  1)   2 x2       r xr 
1  p ( x' )
p ( x)
 exp  0  1 x1   2 x2       r xr 
1  p ( x)

p ( x' )
1  p ( x' )
 exp( 1 )
p ( x)
1  p ( x)
16
心疾患の例
説明変数
(リスクファクター)
X1 年齢
X2 コレステロール
X3 収縮期血圧
X4 相対体重
X5 ヘモグロビン
X6 喫煙
X7 ECG所見
相対体重：
喫煙：
ECG所見：
βの推定値
0.0708
0.0105
0.0166
0.0138
-0.0837
0.3610
1.0460
SE
0.0083
0.0016
0.0036
0.0051
0.0542
0.0587
0.2700
Wald 検定
8.5301
6.5625
4.6111
2.7059
-1.5443
6.1499
3.8741
オッズ比
オッズ比を求
exp(β)
めた際の単位
2.03
10
1.69
50
2.29
50
1.99
50
0.66
5
2.85
1
体重÷年齢別身長別体重の中央値
0: never, 1: 1箱未満, 2: 1箱, 3: 1箱より多い
0: 正常, 1: 正常でない
17
なぜオッズ比か
• オッズ(odds) とは比のこと
• オッズ比．．．比の比
• なぜ「比」だけではダメか
A薬
治癒
未治癒
PA
B薬
PB
1  PA 1  PB
PB
ではまずいのか？
PA
18
例
Ａ薬Ｂ薬
治癒 90
99
未治癒 10
1
Ａ薬Ｂ薬
治癒 50
55
未治癒 50
45
比 PB / PA
0.99/0.90=1.1
0.01/0.10=0.1
?
0.55/0.50=1.1
0.45/0.50=0.9
?
19
解説
• 薬の効きを治癒率の比と未治癒率の比
でみたものとが異なるのは矛盾
• 100名中治癒した割合は1割違うだけで
あるが，90→99と50→55とは評価は異
なるべきであろう
• では，未治癒率で見ればよいというこ
とになるかもしれないが，数値が治癒
率と未治癒率が入れ替わっているきは
同じ問題が起こる
20
オッズ比でみると
• 治癒率のオッズ比は未治癒率のオッズ比
の逆数
PA
1  PB
– 治癒率のオッズ比=2
未治癒率のオッズ比=0.5
1  PA
PB

PB
1  PA
1  PB
PA
• 1の近くでの変化は中庸での変化より高く
評価される
21
例
Ａ薬Ｂ薬
治癒 90
99
未治癒 10
1
Ａ薬Ｂ薬
治癒 50
55
未治癒 50
45
オッズ比
0.99
1  0.99  99  11
0.90
9
1  0.90
0.55
1  0.55  1.22  1.22
0.50
1
1  0.50
22
補足 -種々のモデル式• プロビット回帰モデル
1
 ( p(x))  Z ⇒標準正規分布関数
• complementary log-log回帰分析
log(  log( 1  p(x))  Z
⇒二重指標関数
• ロジスティック回帰分析
p ( x)
log
Z
1  p ( x)
⇒ロジスティック関数
23
1.0
p
二重指数関数
（double exponential
function）
ロジスティック関数
（logistic function）
0.5
0.0
Z
標準正規分布関数
（standardized normal distribution function）
p:確率値
Z:変数の線形な合成変数
24
近似について
• 二重指数関数とロジスティック関数は標準正規
分布関数の近似
• ロジスティックが一般的だが，これといった理由
はない
– オッズ比との相性のよさ
• どの近似を採用しても，データが存在する説明
変数ｘの範囲の中では違いは小さい
• しかし，外挿するときは注意が必要
– 感度分析．．．3種類の関数で推定してみて大きな差
がないことを確認する
25
具体的事例とSASによる分析
26
１９８６年NASAスペースシャトル
CHALLENGER号爆発事故
• 事故調査班は原因は「O-ring」という部品
の故障だと断定
• また、調査班は事故につながる重要な要
因として温度を取り上げている
• 過去のデータから、当時の温度から故障
率を予測するとどのような結果になるか？
27
過去23回のスペースシャトル打ち上げ時の温度と
「O-ring」故障数（全６個中）
故障数
温度
故障数
温度
故障数
温度
故障数
温度
2
53
1
58
0
73
0
76
0
70
1
70
0
67
0
70
0
78
0
81
0
75
0
76
1
57
1
63
0
68
0
68
1
70
0
72
2
75
0
68
0
79
0
66
0
69
28
「故障数」は正規分布ではない
図１．故障の頻度
20
頻度
15
10
5
0
0
1
2
3
4
故障数（全6個中）
5
6
29
温度と故障率の散布図
故障率(故障数/6)
図2　温度と故障率との散布図
0.5
0.4
0.3
0.2
0.1
0
20
40
60
温度（F）
80
100
30
SASプログラム- proc logisticOPTIONS NOCENTER PS=54 LS=90;
DATA d1;
INPUT num nf no temp @@;
CARDS;
1
6
11
16
21
2
0
2
0
0
6
6
6
6
6
53
78
75
72
70
2
7
12
17
22
0
1
0
0
0
6
6
6
6
6
66
57
79
76
73
3
8
13
18
23
0
0
1
0
0
6
6
6
6
6
68
67
58
81
76
4
9
14
19
1
0
0
1
6
6
6
6
70
69
67
63
;
PROC LOGISTIC DATA=d1;
MODEL nf/no = temp / SCALE=NONE COVB
PLRL LACKFIT;
OUTPUT OUT=d2 C=COOK;
PROC PRINT DATA=d2;
RUN;
5
10
15
20
0
1
0
0
6
6
6
6
75
70
70
67
31
Details
• MODEL nf/no = temp / SCALE=NONE PLRL LACKFIT;
OUTPUT OUT=d2 C=COOK;
–
–
–
–
従属変数に「故障数/全体の数」を指定
SCALE=NONE…適合度
PLRL…オッズ比とその区間推定
LACKFIT…いくつかのデータをまとめて，モデ
ルによる予測頻度とデータの頻度との比較
– C=COOK Cook統計量による回帰診断
• 分析に過度の影響があるobservationの同定
32
分析の吟味
33
SAS出力：適合度
ロジスティック関数と線型回帰モデルのよさを吟味
Deviance and Pearson Goodness-of-Fit Statistic
Pr >
Criterion
DF
Value Value/DF Chi-Square
Deviance
21
18.0863
0.8613
Pearson
21
29.9803
1.4276
Number of events/trials observations: 23
0.6435
0.0924
モデルの適合度を調べる
統計量=デビアンス（のｐ値）
大きいほど良い
34
SAS出力：偏回帰係数
Analysis of Maximum Likelihood Estimates
exp(-0.1156)
Parameter Standard Wald
Pr > Standardized
Odds
Variable DF Estimate Error Chi-Square Chi-Square Estimate Ratio
INTERCPT 1
5.0850 3.0525
2.7751
0.0957
.
.
temp
1 -0.1156 0.0470
6.0435
0.0140
-0.441494
0.891
回帰式
p^ (t )
1

1  exp( 5.085  0.1156 t )
Estimated Covariance Matrix
Variable
INTERCPT
INTERCPT
9.3176671947
TEMP
-0.142565536
TEMP
-0.142565536
0.002211241
35
SAS出力：オッズ比の区間推定
Profile Likelihood Confidence Limits
Variable
temp
Unit
1.0000
Odds
Ratio
0.891
Lower
0.809
Upper
0.970
36
SAS出力：予測の「よさ」をみる
Association of Predicted Probabilities and
Observed Responses
Concordant = 65.4%
Somers' D = 0.382
Discordant = 27.1%
Gamma
= 0.413
Tied
Tau-a
= 0.047
c
= 0.691
=
(1161 pairs)
7.5%
37
順位相関係数
NUM
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
NF
2
0
0
1
0
0
1
0
0
1
2
0
1
0
0
0
0
0
1
0
0
0
0
推定確率
0.261
0.073
0.059
0.047
0.027
0.019
0.182
0.065
0.053
0.047
0.027
0.017
0.165
0.065
0.047
0.038
0.024
0.014
0.100
0.065
0.047
0.034
0.024
温度
53
66
68
70
75
78
57
67
69
70
75
79
58
67
70
72
76
81
63
67
70
73
76
NUM'
1-1
1-2
1-3
1-4
1-5
1-6
2-1
2-2
2-3
2-4
2-5
2-6
…
…
NF'
1
1
0
0
0
0
0
0
0
0
0
0
…
…
推定確率
0.261
0.261
0.261
0.261
0.261
0.261
0.073
0.073
0.073
0.073
0.073
0.073
…
…
温度
53
53
53
53
53
53
66
66
66
66
66
66
…
…
NF’と推定確率の順位相関
係数をとったものが
associationの指標
38
定義式
N  23  6  138
t  1161
nc  1161  0.654  759
nd  1161  0.271  315
Somer' s D  (nc  nd ) / t
c  (nc  0.5(t  nc  nd ) / t
Goodman - Kruskal Gamma  (nc  nd ) /(nc  nd )
Kendall' s Tau  a  (nc  nd ) /( N ( N  1) / 2)
NUM
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
NF
2
0
0
1
0
0
1
0
0
1
2
0
1
0
0
0
0
0
1
0
0
0
0
NO
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
6
TEMP
53
66
68
70
75
78
57
67
69
70
75
79
58
67
70
72
76
81
63
67
70
73
76
COOK
0.29503
0.02790
0.02110
0.11088
0.01097
0.00799
0.00340
0.02395
0.01894
0.11088
1.41421
0.00712
0.00002
0.02395
0.01721
0.01443
0.00993
0.00555
0.02428
0.02395
0.01721
0.01322
0.00993
39
COOKの統計量
図2　温度と故障率との散布図
故障率(故障数/6)
OBS
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
0.5
0.4
0.3
0.2
0.1
0
20
40
60
温度（F）
80
100
• 当分析に対する影響度の大
きなobservation を同定する
• No.11 (t=75) は容疑者
• 分析者に都合のよいデータ
除去は，慎むべき
40
分析結果の利用
41
回帰式の利用
Challenger が爆発したとき(t=31)の故障確率は？
• 回帰式のモデル：
p^ (t )
1

1  exp( 5.085  0.1156 t )
。
• 温度が31 Fでの故障確率の点推定値
^
p(31)  0.818
• 6つの「O-ring」のうち少なくとも１つが故障
する確率
1  (1  p(31))  0.999963
^
6
42
95%信頼区間
( p^
L (t ),
p^
^
U
1
(t )) 
SE (   t )  Var ( )  2tCov( ,  )  t Var (  )
^
^
1  exp( 5.085  0.1156 t  1.96 SE (   t ))
^
^
^
^
2
pˆ (31)の区間推定 : (0.16 , 0.99 )
少なくとも一つ故障す
る確率の区間推定
下限値 : 1  (1  0.16 ) 6  0.65
上限値 : 1  (1  0.99 ) 6  1.00
^
43
一つの問題点
• t=31は，分析に使ったデータ範囲を越えて
いる
• これを外挿(extrapolation) という
• 外挿をした場合は，その結果が採用した関
数Fに大きく依存して変化することが少なく
ない．選んだ関数の理論的根拠が希薄な
場合はなおさら
44
対策はどうするのか？
• 他に考えられる関数を適用してみて、結果
がどの程度異なるかという感度分析をする
のがよい。
• この場合に考えられる候補としては
– プロビット回帰分析
– complementary log-log回帰分析
45
6個のうち少なくとも1つが
故障する確率
感度分析
表２．温度３１Fでの予測確率
Logistic
p(31)
0.82
probit
0.70
0.97
0.12-0.99
0.19-1.00
1.000
95% CI
0.16-0.99
Pr(少1|6)
1.000
0.999
0.65-1.00
0.54-1.00
95% CI
comp.log-log
0.71-1.00
46
まとめ
• ロジスティック回帰分析は従属変数が
二値変数の時に用いる
– ３件法のデータを二値変数として分析する
こともある
– 外挿の時は感度分析も忘れない
• プロビット回帰分析
• Comp.log-log回帰分析
• SAS，SPSS等で分析可能
47
文献
• 丹後・山岡・高木(1996)．ロジスティック
回帰分析．朝倉書店

Download Report