医療統計学を学ぶ大学生のブログ

医療統計学、因果推論を専攻しています。R, SASユーザーです。

SASによる解析:LOGISTICプロシジャを用いた二値アウトカムに対するロジスティック回帰

大まかな内容としてはロジスティック回帰分析の概要と、それを実際にSASで実行するとしたらどんな感じのプログラムを組めばいいかを簡単にまとめていこうと思います。

また今回の内容とは関係ないのですが、最近SNSなどを経由してこのブログに関しての応援のコメントもいただくことも増え大変感謝しております。医療分野における解析方法の理論、実装面をできる限り最新のものまで追えるように引き続き頑張っていこうと思うので、今後ともよろしくお願いします。というわけで、ロジスティック回帰の話に入っていきます。

二項ロジスティック回帰の概要

ロジスティック回帰は基本的には通常の回帰分析の延長上の話です。ただ通常の回帰分析では目的変数が連続量であるのに対して、ロジスティック回帰では質的な変数が目的変数となります。この質的な目的変数が二値であれば二項ロジスティック回帰、3つ以上のカテゴリーを持つのであれば多項ロジスティック回帰、順序尺度であるならば順序ロジスティックなどとよりその手法は細分化していきます。

今回は二項ロジスティック回帰を取り扱うのですが、医療統計の分野においては、(生存、死亡)、(無発症、発症)、(非介入、介入)など様々な二値のアウトカムを考える場合が多々あり、そういった場合によくロジスティック回帰は用いられていますし、以前に紹介したIPWのような手法でも傾向スコアの算出のために使われていたりする手法です。

 

モデル式

まず次のように文字を定義します。
Y:目的変数
Xj:p個の説明変数(j=1, ・・・, p)
β0:切片項
βj:p個の偏回帰係数(j=1, ・・・, p)

 

重回帰モデルにおいては以下のような回帰式を考えます。

f:id:NorihiroSuzuki:20211117135239p:plain

 

これに対して二項ロジスティック回帰ではイベントの発生(介入)確率pをとした場合にpのロジットを取ったものを目的変数として置きます。すなわち、

f:id:NorihiroSuzuki:20211117135151p:plain

という回帰モデルを考えているわけです。この目的変数をpではなく、logit(p)としている点については次のように考えることで理解が可能かと思います。

まずpは確率ですので当然のことですが0~1の範囲をとります。ですが右辺は-∞~∞の範囲をとるため対応はしません。そこでpのオッズ、すなわちp/(1-p)をまず考えてみるとこれは0~∞の範囲を取ります。ここで対数関数は真数が0より大きい範囲において定義可能で、このオッズを自然対数をとったものは前述のように0~∞を取りますので、オッズの対数をとったものを考えることができます。すなわちpに関してロジットを取ったもの、対数オッズは対数関数のグラフからも分かるように-∞~∞の範囲をとり、右辺と左辺の範囲は対応をします。

 

偏回帰係数の解釈

回帰モデル式についてざっとわかったところで次に問題となるのが、推定された偏回帰係数の結果をどう解釈するのかというところかと思います。重回帰の時と同じく各偏回帰係数の値は、説明変数の値を1増加させたときに増加する目的変数の値を指しますので、ロジスティック回帰モデルにおける偏回帰係数は説明変数が1増加したときの対数オッズの増加量を意味します。

ですが、対数オッズの増加量と言われても多くの方がいまいちその意味合いがわからないため、通常は両辺のexp(指数)をとったものを考えます。すなわち対数をとることによって、

f:id:NorihiroSuzuki:20211117143847p:plain

となり、どの程度オッズが変動するかをみることが出来ます。

 

例えば説明変数が一つである以下のようなロジスティックモデルを考えてみます。

f:id:NorihiroSuzuki:20211117145318p:plain

ここでX1=0, 1の時のそれぞれのオッズを考えると、

f:id:NorihiroSuzuki:20211117145500p:plain

となります。log(a/b)=log(a)-log(b)とlogの引き算の形に変換することができるため、X=1をイベントとした場合のオッズ比は、

f:id:NorihiroSuzuki:20211117150111p:plain
f:id:NorihiroSuzuki:20211117150237p:plain

となることから、偏回帰係数のexpをとったものがオッズ比の推定量となります。

 

ちなみにオッズ比にはpが小さい時、リスクが小さい時(稀なイベント)にリスクの矜持となるため、医療分野における解析ではしばしばリスクとして捉えられています。書籍やサイトによってはオッズ比と言わずにリスクと書いているものもたまに目にします。オッズとリスクに関してはこちらをどうぞ。

norihirosuzuki.hatenablog.com

 

SASでの二項ロジスティックモデル

使用データ説明

二項ロジスティック回帰の説明にあたって使用するデータセットは、SASHELP内にあるBirth Weight Data(BWeight)です。これはアメリカの国立衛生統計センター(NCHS; Natinal Center from Health Statistics)によって提供された1997年のアメリカの18~45歳の母親によるシングルトン妊娠での出生した胎児の体重に関するデータとなっています。データセットに含まれる項目としては以下のような形になっています。(各項目の内容はラベルの通り)

 

f:id:NorihiroSuzuki:20211118020715p:plain

データのサイズは50,000です。

 

LOGSTICプロシジャでの実装

LOGSTICプロシジャでは次のように指定を行うことで二項ロジスティック回帰を行うことが出来ます。

 

f:id:NorihiroSuzuki:20211118110409p:plain

カテゴリー変数はダミー化を行う必要があり、ダミー化を行う必要のある変数はmodelステートメント、modelステートメントの両方に記載することで自動的に処理が行われます。

 

一応全部のオプションはこんな感じで大量にあります(詳しくは公式のHELPを)

f:id:NorihiroSuzuki:20211117152921p:plain

 

SAS公式のHELPはこちら。

documentation.sas.com

 

二項ロジスティック回帰を行うにあたってはまずデータセットを少し弄ります。Weightという変数には出生児の体重に関しての数値データが含まれているのですが、ここでは2500(g)未満を低体重時(=1)とする新たな変数tr1_Weightを作成します。全部で50,000あるオブザベーションは、値が0(2500g以上)、値が1(2500g未満)に下表のように分かれます。

f:id:NorihiroSuzuki:20211118102415p:plain

今回はこのtr1_Weightのロジットを取ったもの(関心のあるイベントは=1)を目的変数として二項ロジスティック回帰を行います。

 

コピペ用のコードは最後にまとめてありますが、以下のようなコーディングを行います。

f:id:NorihiroSuzuki:20211118102921p:plain

 

この実行結果はこんな感じ。

f:id:NorihiroSuzuki:20211118103845p:plain

オッズ比が仮にリスク比に近似できると考えると、例えば人種の違い(黒人、白人)によって低体重児の出生リスクは約2倍ほど異なることが伺えます。今回は各変数の主効果のみを考え、交互作用に関しては特に何も考慮はしていません。

 

コード

title 'Birth Weight Data';
proc contents data=sashelp.BWeight varnum;
   ods select position;
run;
data test;
set SASHELP.BWEIGHT;
tr1_Weight=ifn(Weight<2500, 1, 0);
run;
ods noproctitle;
ods graphics / imagemap=on;

proc logistic data=test;
class Black Married Boy MomSmoke Visit MomEdLevel / param=glm;
model tr1_Weight(event='1')=Black Married Boy MomSmoke Visit MomEdLevel MomAge CigsPerDay MomWtGain / link=logit technique=fisher;
run;