« 読了：Rothman(1990) 多重比較法は使うな | メイン | 読了：「あさひなぐ」「アヴァール戦記」「ひなまつり」「しろくまカフェ」「岳」 »
2012年3月13日 (火)
Mittlbock, M., & Schemper, M. (1996) Explained variation for logistic regression. Statistics in Medicine, 15, 1987-1997.
ロジスティック回帰分析のアウトプットには，なんか変な$R^2$がいろいろはいっているんだけど，あれっていったいなんなんだろうなあ，と前から不思議に思っていた。このたび関連する話題についてちょっと考える機会があったので，適当に論文を見繕って読んでみた。12種類の$R^2$指標を比較しました，という論文。そんなにあるのか。
12種類の$R^2$指標は，大きくいえば3グループにわかれる。ええと，個体$i$が持つ従属変数の値を$y_i (=\{0,1\})$とする。$y_i = 1$を仮に成功と呼ぶとして，全体における成功率を$\bar{p}$とする。また，独立変数の値を$x_i$とし、ロジスティック回帰モデルで推定された各個体の成功率を$\hat{p}_i$とする。
最初のグループは，$y_i$と$\hat{p}_i$の相関の二乗に基づく指標。
- $y$が二値であることは気にせずに、単純にPeasonの相関係数を二乗する。($r^2$)
- Spearmanの相関係数を二乗。($r^2_s$)
- Kendallの$\tau_a$を二乗。($\tau^2_a$)
- Kendallの$\tau_b$を二乗。($\tau^2_b$)
- Somersの$D$を二乗。($D^2_{\hat{p}y}$)
- Goodman-Kruskalの$\gamma$を二乗。($\gamma^2$)
順位相関係数を手当たり次第に集めてきて二乗しました、という感じですね。3番目以降はまあどれも似たような指標である(いずれも計算式の分子は同じ)。
第二のグループは，yの分散の縮減率に基づく指標。一般化していうと，conditionalな残差を表すなんらかの指標 $D(y_i|x_i)$の合計と、unconditionalな残差を表すなんらかの指標$D(y_i)$の合計を出して、比をとって1から引くタイプの指標である。
- yが二値であることなど気にせず、ふつうの回帰モデルと同じように，$D(y_i)=(y_i-\bar{p})^2$, $D(y_i|x_i)=(y_i-\hat{p}_i)^2$とする。($R^2_{SS}$)
- ロジスティック回帰モデルを信じれば、unconditionalな残差の分散の期待値は $\bar{p}(1-\bar{p})$、conditionalな残差の分散の期待値は$\hat{p}_i(1-\hat{p}_i)$なわけだから、これを$D(y_i), D(y_i|x_i)$として用いる。これをGiniのconcentration measureという。へー。($R^2_G$)
- アタリハズレだけに注目し、unconditionalな残差の絶対値$|y_i-\bar{p}|$と、conditionalな残差の絶対値$|y_i-\hat{p}_i|$を、なにかのカットオフ(たとえば0.5)よりも上か下かで二値変数に落としてしまい、これを$D(y_i)$, $D(y_i|x_i)$として用いる。モデルのせいでハズレ率がどれだけ減ったか、だけに注目するわけだ。ずいぶん荒っぽい話だなあと驚いたが、これはGoodman-Kruskalの$\lambda$と等価である由。ほー。 ($R^2_{CER}$)
- 残差指標をエントロピーとして捉える。$D(y_i)=-(y_i \log \bar{p} + (1-y_i) \log (1-\bar{p}))$, $D(y_i|x_i)=-(y_i \log \hat{p}_i + (1-y_i) \log (1-\hat{p}_i))$。前者の和は -(ナルモデルの対数尤度)、後者の和は -(モデルの対数尤度)となる(これを2倍したのが、ソフトの出力に出てくるところのモデルのデビアンスであろう)。で，これを引き算してカイ二乗検定するんじゃなくて，比をとって1から引く。こうして求めた$R^2$は、McFaddenのpseudo $R^2$というのと等価である由。($R^2_E$)
最後のグループは，モデルの尤度に基づく指標。
- さっきのは1-(対数尤度の比)だったが、今度は1-(尤度の比)^{2/n}とする。つまり、どれだけ尤度が上がったかを一人あたり幾何平均で出すわけだ。これを尤度比$R^2$、ないしCox-Snellの$R^2$という。($R^2_{LR}$)
- 仮に完全にフィットしているモデルでも，$R^2_{LR}$は1にならない。それじゃ使いにくかろうというので、$R^2_{LR}$をその最大値，すなわち1-(ナルモデルの尤度)^{2/n}でさらに割ってやる。実にムリヤリ感あふれるアプローチだが、これがあの，SPSSとかで出てくる，Negelkerkeの$R^2$なのだそうだ。へえええ！ ($R^2_{CU}$)
いま調べたら、SAS 9.22のproc logisticでは、$R^2_{LR}$が"RSquare"というラベルで、$R^2_{CU}$が"Max-rescaled RSquare"というラベルで出力されるらしい。
なんだかもうお腹一杯だが，論文のほうはここからが本番で，いろいろデータをつくっては12種類の$R^2$を求め，挙動を比較している。著者らいわく，
- 直観的にわかりやすいか
- ロジスティック回帰の性質とつじつまがあうか
- 0から1まで動くか
- ふつうの線形回帰モデルがうまく当てはまるようなデータでは線形回帰モデルの$R^2$と同じような値になるか
という4つの基準で検討すると，
- $R^2_E, R^2_{LR}, R^2_{CU}$は直観的に理解しにくい。
- $R^2_{CER}$や順位相関係数系の方法はロジスティック回帰とつじつまがあわない。
- $r^2_S$, $\tau$系, $R^2_{LR}$は完全にフィットしても1にならない。
- たいていの手法は線形回帰の$R^2$とずれる。
というわけで，消去法で結局$r^2, R^2_{SS}, R^2_{G}$が残ることになる。結局のところ，$y_i$が二値変数であることを無視しちゃった方がいいね，という，ちょっと奇妙な結論である。
勉強にはなったけど。。。ロジスティック回帰モデルのときも，慣れ親しんだ$R^2$のような奴が欲しいよ欲しいよ欲しいよ，というのがこの論文の前提になっているところがポイントだと思う。この前提そのものが，なんだかちょっとワガママなような気がしてならない。なんというかその，海外旅行先で味噌汁飲みたい，といっているように聞こえる。いや，それは飲みたいですけどね，私も。
それはまあいいや。ともかく，Negelkerkeの$R^2$というのが意外に無茶な発想で作られている，という点を学んだのが収穫であった。私の知る狭い範囲の話だが，Negelkerkeの指標はよくみかけるような気がする。SPSSが出力するせいかもしれない。
論文：データ解析(-2014) - 読了：Mittlbock & Schemper (1996) ロジスティック回帰の説明率指標を品定め