Where My Logistic Regression Went Wrong

Adjusting for Bias IV − Where My Logistic Regression Went Wrong
Keywords: confounding & collapsibility, effect measure, generalized linear model, language & writing
どうしてオッズ比が無限大になっちゃったの?
お父さん「あのさ、昨日のロジスティック回帰の表の続きをしない?コーヒー淹れたから」
私「えーお手やらわかに。コーヒーもありがとう。この表の話だよね」
| パラメータ | 自由度 | 推定値 | 標準誤差 | 95%信頼区間 | \(\chi^2\) | p値 |
|---|---|---|---|---|---|---|
| Intercept | 1 | -27.25 | 3.03 | -33.20~-21.30 | 80.55 | <0.01 |
| STOMA | 1 | 1.41 | 1.54 | -16.1~4.45 | 0.84 | 0.3588 |
| AGE | 1 | 1.47 | 1.48 | -1.42~4.36 | 0.99 | 0.3197 |
| SEX | 1 | 25.72 | 0.00 | 25.72~25.72 | ||
| INCOME | 1 | -0.00 | 0.00 | -0.00~0.00 | 0.58 | 0.4481 |
お父さん「そう、これね、間違いではないんだけど、まだよくできる。表を作成するときには、情報が自然に入ってくるように、読者の立場で考えなくちゃ。具体的な注意点はこんなところかな」
- 読者の立場で考える
- 表のフォーマットに従う
- 単位・表示桁数に気を配る
- 表タイトル・脚注を利用して補足情報を伝える
- 数え尽くしの原則を守る
お父さん「ロジスティック回帰からいろんな指標が出せるけど、表にはオッズ比、95%信頼区間、p値だけあれば普通はじゅうぶんだよ」
私「それはこの表にも載ってるよね」
お父さん「載っていません。この表の”推定値”は、みたところ回帰係数の値でしょ。それを換算してオッズ比のスケールに直さないといけない。換算するには、回帰係数\(\beta\)の指数をとって、\(OR=\exp(\beta)\)を求めればいい」
私「そうなのか。ストーマのオッズ比は\(\exp(1.41)=4.1\)ってことね。あれ?\(\exp(25.72)\)を計算してみると、パソコンに表示できないくらいに大きい数字になったんだけど。ほとんど無限大に近い数字だ。性別の影響が強いのかな。年収は\(\exp(-0.00)\)?なんか変だな」
お父さん「これはよくあるトラブルだよ」
私「そうなの?」
お父さん「復職しなかった患者さんはほとんど女性だったんじゃない?擬似完全分離(quasi-complete separation)っていうデータの配置上の問題が起きちゃってて、性別の回帰係数を、うまく推定できていないんだ。この先は少し長くなりそうだから、ちょっとたばこ吸わせてもらうね」
ロジスティック回帰の回帰係数は、最尤法(maximum likelihood)という汎用的な推定方法によって計算されます。これはある種の方程式を解くことに対応していますが、データから必ずしも解が求まるわけではありません。これが完全分離(complete separation)や擬似完全分離(quasi-complete separation)という問題です。これらの問題が起きると、回帰係数が求まらなかったり、不安定だったりするため、共変量の数の削減や、最尤法以外の推定方法の利用など、なんらかの根本的な対処が必要になります。
お父さん「(灰皿にたばこを押しつけながら)正しい計算ができていないわけだから、ロジスティック回帰の結果を示すのはやめた方がいい。性別だけじゃなくて年収も問題だね。単位が円のまま扱うと、桁が大きすぎて数字が読めない。単位を万円にして、年収が1万円増えたら復職率がどう変化するかを示すといいんじゃないかな。性別と年収の扱いは気になると思うけど、今はこれくらいにしよう。表の体裁を直してあげるよ」
| オッズ比 | 95%信頼区間 | p値 | |
|---|---|---|---|
| ストーマなし | (基準) | ||
| ストーマあり | 4.13 | 0.20~85.23 | 0.36 |
| 65歳未満 | (基準) | ||
| 65歳以上 | 4.35 | 0.24~78.60 | 0.32 |
| 女性 | (基準) | ||
| 男性 | ∞ | ∞~∞ | |
| 年収(万円) | 1.00 | 0.98~1.01 | 0.45 |
*手術後1年以内の復職をアウトカムとした多変量ロジスティック回帰
私「あれ?オッズ比4.13倍ってことは、復職率がだいたい4倍になったってこと?そんなことある?」
お父さん「違う違う。あのね、オッズ比はオッズの比をとったもの。オッズっていうのは、復職率を\(\pi\)で表すと、それを\(\pi/(1-\pi)\)と変換したものなんだ。\(\pi\)が0に近いときは、オッズ比は\(\pi\)の比に近い値になるんだけど、比の分子と分母の取り方で、向きが変わることに注意しないといけない。今回は、モデルの定義上、ストーマありの方が”復職できない”オッズが約4倍という意味だよ。ただし、2値アウトカムを”復職できない”に取ると、向きは逆になる」
私「わかりにくいのう」
お父さん「じゃあ、どこに気を付けて手直ししたか説明するね。ロジスティック回帰では、共変量はなにで、どう扱ったかを正確に伝えることが大切なんだ。データをみると、年齢は65歳未満か以上かで分類した2値データみたいだね。そういうときは、何歳で分類したか、どの分類を基準にしたかわからないと、結果を読み取れない。また、連続データでは、忘れがちだけど単位が必要。年収の単位を円から万円にして、表に追記しておいたよ」
論文で示される表には、決められたフォーマットがあります。表は、行(row)と列(column)から構成され、行と列は概念的には対称ですが、視覚的には役割が異なります。英語圏では、表は左から右に、行ごとに読むことが慣習になっています。そのため、このサイトに示されている表のように、縦線で区切ることはしません。
数字の意味をわかりやすく伝えるためには、どこにどのような情報を配置するかに気を配るべきです。さきほどの表では、列が指標の種類を表すというシンプルな方針を採用しています。そのため、数字がオッズ比、95%信頼区間、p値に対応することが一目でわかります。また、必要な情報から順に読み取れるように、各行または各列の順序にも気を配りたいところです。
ある程度大きな表であれば、情報はグループにわけられることが普通です。表のフォーマットも、グループが視覚的にわかりやすいように工夫できます。たとえば後に出てくる表1では、行と列でそれぞれグループ化がなされています。行をグループ化するときには字下げ(indent)を用います。性別の内訳を示すために字下げがなされていますよね。また、列をグループ化するときにはスパナ(spanner)を用います。“ストーマ保有者”、“ストーマ非保有者”といった列の見出しを、スパナヘッド(spanner head)といいます。
表では測定単位(unit)を明らかにすることが重要です。臨床検査値を扱う論文を投稿するときは、その雑誌の編集方針をみて、どのような単位系を採用しているか確認しましょう。
数字の表示桁数は、ある程度桁数やルールが統一されていなければ、印象を悪くします。統計ソフトウェアが出力する数字は桁数が大きいことがありますが、そのまま表に示す必要はありません。有効数字2桁または3桁か、もともとの測定値より桁数を1桁大きく示すのが基本です。
日本工業規格(Z 9041-1:1999)によると、平均の有効数字は、測定値より1~2桁大きく、標準偏差の有効数字は、最大でも3桁程度といわれています。つまり、むやみに数字の表示桁数を増やしても精確性が増すわけではないのです。年齢を例に挙げると、測定単位は「歳」なので、小数点1桁まで示すとよいでしょう。表1の年齢の平均・標準偏差では、そのようになっています。p値については、小数点以下2桁か3桁まで報告すれば、有意水準(0.05が多い)との大小関係が読み取ることができるので、必要十分でしょう。
私「表の作り方をご助言いただけるのはありがたいんだけどさ、肝心のロジスティック回帰はどうすればいいんだろ」
お父さん「擬似完全分離を避けるには、共変量の数を削るのがわかりやすい。今回はサンプルサイズが不足しているみたいだから、他の共変量はすべて使わないで、単にストーマ造設と復職率だけの関連を調べたらどうかな?1パラメータあたり5〜10イベント程度は欲しいといわれるからね。その場合も、通常の\(\chi^2\)検定は使わない方がいいよ。Fisherの正確検定とClopper-Pearson信頼区間っていう方法が推奨されている」
- ロジスティック回帰の共変量の数を減らす(目安としてパラメータの数のおよそ5倍のイベント数が必要)
- \(\chi^2\)検定ではなくFisher正確検定を用いる
- 2項分布に基づく信頼区間(Clopper-Pearson信頼区間)を用いる(田中2022)
足りなかったたくさんの情報
私「でも、ロジスティック回帰で交絡因子を調整しなさいっていったのはお父さんなんだよ」
お父さん「そうだね、それは重要な問題だよね。じゃあ、ストーマ保有者と非保有者で、年齢、性別、年収の分布がどれくらい違うかをまず確認しようか。前にも話したとおり、年収を尋ねたことでバイアスが生じていないかとか、年収を尋ねて欠測が多いんじゃないかとか、いろいろ気になるし。もし分布が揃ってたら、調整しなくても影響は小さいかもしれない」
私「あ、その表なら資料にあるよ」
お父さん「じゃあ、あとはストーマ造設と復職率の関連を表にするだけだね。資料には、以下の表1と表2を載せればいいよ」
表1. 20XX~20YY年にZZ病院で手術を受け、がんサバイバー復職率調査に回答した直腸がん患者80人の背景因子の記述
| ストーマ保有者(N=20) | ストーマ非保有者(N=60) | |
|---|---|---|
| 年齢(歳) | 57.0 (8.6) | 65.5 (10.4) |
| 65歳以上 | 10.0% | 51.7% |
| 性別 | ||
| 男性 | 35.0% | 53.3% |
| 女性 | 60.0% | 46.7% |
| 不明 | 5.0% | 0.0% |
| 年収(万円) | 508 (77) | 514 (101) |
*平均(SD)または%
表2. 直腸がん患者のうちストーマ保有者20人と非保有者60人の復職率の比較
| N | 1年以内に復職せず | 復職率 | 95%信頼区間* | p値† | |
|---|---|---|---|---|---|
| ストーマあり | 20 | 1 | 95.0% | 75.1~99.9% | |
| ストーマなし | 60 | 2 | 96.7% | 88.5~99.6% | 1.00 |
*Clopper-Pearson法
†Fisher正確検定
お父さん「集団のアウトカムを比較したいときは、その前に表1のようにそれぞれの集団の背景因子を記述するべきだよ。ストーマ保有者は20人しかいなくて、非保有者に比べて、平均年齢が低くて、女性が多いわけだ。やっぱり最終的にはロジスティック回帰で、背景因子の違いを調整しないといけないかもね。現時点での調査の回答率はどうだった?」
私「調査票は100人に送って、80人から回収できたから、回答率は80%かな」
お父さん「項目ごとの無回答(item non-response)についてはどう?」
私「そっか、調査票を回収できたかと、項目ごとの回答があったかは別だもんね。年齢は全員回答してくれたよ。性別は1人記入していなかった。項目ごとの回答数まで、表に示した方がよかった?」
お父さん「ううん、論文ではそこまで細かい数字は示さないことが多い。重要な変数が欠測してたりするなら別だけどね。たとえばさ、ロジスティック回帰を行うとき、性別不明はどう扱ったの?」
私「統計ソフトウェアのマニュアル調べたら、欠測は自動的に除外されるんだって」
お父さん「そうすると解析対象の人数が80人じゃなくて79人に減るんだから、表を作るときにそこを説明しておかないといけないんじゃない?一般に、ロジスティック回帰の結果を示す表だけだと、どんなデータかわかりにくいっていう問題がある。復職できなかったのは3人しかいないって、表2をみるまで気づかなかったよ?まだ1施設からのデータしか集まってないからでしょ。それじゃあロジスティック回帰はしたくてもできないよね」
私「確かにね。そのあたりの状況を、読者に伝えなきゃって気持ちが足りなかったみたい」
読者が表を読むときのことを想像してみてください。最初に目にするのは、表に示された結果自体ではなく、表タイトルの情報です。簡潔な表タイトルをよく目にしますが、これはもったいないことです。表1、表2の表タイトルには、調査が行われた状況、対象疾患、人数、統計解析の目的などが含まれています。表タイトルを決めるときは、可能な限り多くの情報を伝えるように工夫しましょう。
表には様々な数字が示されますから、どういった指標かわかりにくいときがしばしばあります。もしかしたら、複雑な統計手法から計算されたものかもしれません。データを解析するために用いた統計手法は、方法のセクションで説明することが基本ですが、読者への配慮として、表の脚注(footnote)を利用することもあります。表1でいえば、括弧内の数字の意味は自明ではないので、“平均(SD)または%”と脚注がついていますよね。表2の脚注には、どの統計手法を用いて信頼区間とp値を計算したかが説明されています。このように読者にとってなじみのない指標や表現には、補足説明が必要です。
脚注を参照するための記号には、a、b、cまたは*、†、‡、§、¶を用いるのが正式です。
統計解析の経験があるものなら誰でも、様々な理由で対象者が除外されるため、個々の解析で、人数が必ずしも同じにならないことを知っているものです。透明性を高めるために、統計解析ごとに、何人が対象となったか報告するべきです。
統計家の間では、数え尽くしの原則(principle of exhaustion)という言葉があります。分類変数の分類が複雑だったり、欠測データがあったりするときには、分類ひとつひとつの人数を合計し、それが全体の人数に合うかどうかを確認せよ、というのがこの原則の教えです。最初に作成した表では人数が示されておらず、どのようなデータが解析されたのかわかりません。一方、表1・表2では、何人が解析対象で、何人が復職したのかが記述されています。さらに、数え尽くすと全体80人や100%に一致することから、正しく集計されたことが読み取れます。
私「ちなみにお父さん、もし”臨床疑問と研究仮説”で例を挙げてくれたときの研究仮説3みたいに、研究の目的が復職率と関連する因子の探索だとしたら、何か変わることはある?」
お父さん「どんな特徴を持った患者が、復職が難しいのか、予測因子をみつけたいってことかな。予測因子の候補となるデータがあるなら、それをロジスティック回帰の共変量に入れることで、予測因子の探索を行うことができるよね」
私「うんうん」
お父さん「こういった探索的検討では、特定の因子に関心がないってことが第一の違いだよ。ストーマと復職率の関連を調べるときは、ストーマのオッズ比やp値を求めることがゴールでしょ。予測因子の探索では、予測因子の候補をできるだけたくさん集めて、どの因子をモデルに入れたら予測精度が高いか、逆にどの因子は要らないかを調べることがゴールになる。これを変数選択っていうんだ」
私「ゴールは違ってもさ、p値を使えばよくない?有意な因子をモデルに入れたらいいんでしょ」
お父さん「そういうやり方もあるけどね。予測2乗誤差、ROC曲線、AICといった予測精度専用の指標を使うことも多い。大事なのはね。予測因子の候補や変数選択の基準を、研究計画書などに書いておくこと。事前に決めておかないと、どうしても研究者の主観が入っているんじゃないかって批判されるからね」
ロジスティック回帰の起源として正しいのは、次の時代のうちどれでしょうか。
- 1900年以前
- 1900~30年
- 1931~60年
- 1961年以降
- 正解は1です
少なくとも、ベルギーの数学者Queteletが1838年に出版した書籍に、ロジット関数の式が示されています。ロジスティックという単語は、同じくベルギーの数学者Verhulstが1845年にProceeding of Belgian Royal Academyで、初めて用いられたとされています。
文献
次のエピソードとRスクリプト
このシリーズのエピソード
- Understanding Collapsibility of Effect Measures: Marginal vs Stratified
- From Risk to Logistic Regression
- Logit: How a Transformation Shapes an Effect
- Where My Logistic Regression Went Wrong
- Why Logistic Regression Fails in Small Samples
- Understanding Confounding in Effect Measures: Marginal vs Stratified
過去のシリーズ
用語集