「統計学関連なんでもあり」の過去ログ--- 041

No.06449　ロジスティック回帰における、定数項が「ないこと」の意義　　【まだまだ修業中】　2008/04/30(Wed) 18:32

はじめまして，SPSS利用者です。初めて質問させていただきます。

ロジスティック回帰において，個体の検査陽性／陰性を目的変数として，それに関連する疫学要因の検討をしております。調査手法は横断研究です。

ステップワイズ減少法で，オプションをデフォルトで実行すると，1変数と定数項（ともにp<0.01で有意）が残されました。しかしながら，適合度検定(Hosmer-Lemeshow)の有意確率が#Null!となってしまい，計算されません。いろいろと条件（投入，除去，Exp(B)の信頼区間，など）を変更しても状況は変わりませんでした。

そこで，オプションの「モデルに定数項を含む（S）」のチェックをはずしてみたところ，生物学的にも納得のいく変数が残り，適合度もまずまず，という結果が得られました。

この意義を調べているのですが，ロジスティック含め，回帰分析の定数項に関しては，「積極的に解釈をする必要はない」という記述がある一方，横断研究やコホート研究では「すべての曝露要因がない場合の有病確率（のロジット）」として解釈可能，というものを読みました。

しかしながら，「定数項を含まない」モデルの意義そのものについての記述を見つけることができません。私の今回の分析の目的は検査陽性／陰性に関連する要因の検討ですから，モデルに定数項を含もうが含むまいが気にする必要がない，と考えてよろしいでしょうか。

本サイトの過去ログを検索したところ，関連すると思われるものとして，「[5394へのレス(5413)] Re: Rで集計」にて，「実質的には，定数項はfactorの最初のカテゴリーに振り分けられ，その他のカテゴリーの回帰係数に調整が加えられます」とのコメントを拝見いたしました。
これを私なりに解釈しますと，「定数項を含めないで作成されたモデルは，当初あった（はずの）定数項を各要因に振り分けることで調整されたものである。したがって，今回の定数項がないモデルに残された要因については妥当なものとして考察可能」ということでよろしいでしょうか。

よろしくお願いいたします。
長文失礼しました。

No.06450　Re: ロジスティック回帰における，定数項が「ないこと」の意義　　【青木繁伸】　2008/04/30(Wed) 22:42

独立変数にカテゴリー変数が含まれる場合には，定数項はどれか一つの（普通は最初の）カテゴリー変数のベースラインとしてモデルに含まれる。以下では，変数 f がカテゴリー変数（R でいう factor）で，2つのカテゴリー a, b を持ち，モデルでは fa, fb と表される。x1, x2 は連続変数である。
以下の2つのモデルでは fb に対する係数が異なるように思うかも知れないが，(fa,fb) が (0,1) の場合と (1,0) の場合（(0,0) や (1,1) はない）を考えると両方の結果は同じであることがわかる。
> glm(y ~ f+x1+x2, data, family=binomial)

Coefficients:
(Intercept)           fb           x1           x2  
  -5.767555     0.252387     0.008084     0.011471  

> glm(y ~ f+x1+x2-1, data, family=binomial)

Coefficients:
       fa         fb         x1         x2  
-5.767555  -5.515168   0.008084   0.011471
しかし，カテゴリー変数を含まない場合には，このようなことはない（定数項の振り分け先がない）。場合によっては，かなり異なる（符号や有意性までも異なるような）結果になることがある。
> glm(y ~ x1+x2, data, family=binomial)

Coefficients:
(Intercept)           x1           x2  
  -5.645581     0.008297     0.011386  

> glm(y ~ x1+x2-1, data, family=binomial)

Coefficients:
       x1         x2  
-0.028487   0.008353
以下のようなテストデータで色々試すとよいかも。
data <- structure(list(x1 = c(118, 150, 110, 135, 130, 130, 110, 
152, 120, 136, 112, 125, 128, 176, 152, 140, 155, 150, 
150, 140, 135, 104, 155, 142, 124, 125, 140, 150, 115, 
130, 140, 135, 160, 130, 130, 120, 150, 136, 130, 158, 
124, 136, 120, 124, 124, 124, 130, 136, 128, 114, 140, 
138, 145, 128, 145, 130, 135, 110, 156, 135, 136, 128, 
126, 145, 156, 126, 130, 138, 136, 126, 135, 122, 180, 
132, 108, 120, 118, 122, 152, 118, 114, 148, 135, 124, 
125, 124, 138, 124, 118, 132, 126, 148, 146, 132, 120, 
106, 114, 124), x2 = c(227, 310, 212, 198, 265, 227, 
205, 284, 171, 178, 258, 255, 194, 263, 238, 338, 212, 
212, 162, 285, 162, 237, 270, 178, 254, 224, 216, 305, 
198, 227, 191, 200, 216, 270, 250, 227, 285, 178, 198, 
243, 220, 285, 460, 237, 191, 162, 250, 225, 292, 242, 
232, 225, 264, 134, 232, 213, 317, 200, 255, 175, 223, 
188, 167, 173, 217, 255, 184, 217, 263, 180, 232, 234, 
171, 180, 148, 192, 263, 173, 232, 188, 223, 255, 213, 
225, 202, 223, 263, 144, 188, 171, 227, 278, 196, 300, 
192, 150, 179, 134), y = c(1, 1, 1, 1, 1, 1, 0, 1, 
0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0), f = structure(c(1, 
2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 
2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 
2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 
2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 
2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 
2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 
2), .abe = c("a", "b"), cass = "factor")), .Names = c("x1", 
"x2", "y", "f"), row.names = c(NA, -98), class = "data.frame")
glm(y ~ x1+x2, data, family=binomial)
glm(y ~ x1+x2-1, data, family=binomial)
glm(y ~ f+x1+x2, data, family=binomial)
glm(y ~ f+x1+x2-1, data, family=binomial)
なぜ，「適合度検定(Hosmer-Lemeshow)の有意確率が#Null!となってしまい，計算されません。」ということになるのか，理由はわかりませんが，何かデータに特殊な状況が生じるのでしょうね。

No.06453　Re: ロジスティック回帰における，定数項が「ないこと」の意義　　【まだまだ修業中】　2008/05/01(Thu) 10:41

ご丁寧な解説，またサンプルデータまでご提示いただき，ありがとうございました。
じっくり読ませていただきます。

No.06455　Re: ロジスティック回帰における，定数項が「ないこと」の意義　　【まだまだ修業中】　2008/05/01(Thu) 11:29

ご説明を読み，確かに同じというのが分かりました。
これについては，式の展開のしようで常識的にわかってくるものでしょうか。

あるいは，こういったことを解説した成書ないし論文がありましたらご教示いただけるとありがたいです。これまでデータをソフトウェアについつい「流し込んでいた」私ですが，今後は理論的な背景もしっかりつかんでおかなければと痛感いたしました。

さて，説明不足でしたので，早速ですが追加させてください。

最初のモデルは検査陽性／陰性に関連する要因として，1つのカテゴリ変数A（なし／あり）と定数項が残されました。これの適合度が悪かったため，定数項を含まないと設定したところ，現在のモデルではカテゴリ変数A（なし／あり）とカテゴリ変数B（なし／あり）が残された，という状況でした。したがって連続変数は（もともとのデータ的にも，カテゴリ化しておいたため）含まれておりません。
なお，各カテゴリ変数（要因）のベースラインとなっているのは「なし」であり，算出されるオッズ比は要因が「なし」の個体に対して「あり」の個体はどれくらい（何倍くらい）陽性になりやすいか，ということを表現しているモデルであると理解しています。また，これらカテゴリ変数AおよびBの間に有意な交互作用は認められませんでした（p＝0.6程度）。

さてこの場合，話が戻って恐縮ですが，モデルに残されているカテゴリ変数AおよびBのそれぞれのベースラインのカテゴリ（この場合「なし」）に，調整された「定数項の振り分け」がなされている，と考えて差し支えないでしょうか。

よろしくお願いいたします。

No.06456　Re: ロジスティック回帰における，定数項が「ないこと」の意義　　【青木繁伸】　2008/05/01(Thu) 13:27

No. 6450 のデータ定義とプログラム例において，list が ist，class が cass になっておりました。編集ミスでした。訂正します。

No. 6450 のデータ定義とプログラムの後に，以下を加えて実験してみればわかるでしょう。
data$x1 と data$x2 で，それぞれ平均値より大きな数値を hi, 小さな数値を lo という factor に変換して data$A, data$B という変数にしています。
data$A <- factor(ifelse(data$x1>132, "hi", "lo"), levels=c("lo", "hi"))
data$B <- factor(ifelse(data$x2>223, "hi", "lo"), levels=c("lo", "hi"))
glm(y ~ A+B, data, family=binomial)
glm(y ~ A+B-1, data, family=binomial)
実行結果は，
> glm(y ~ A+B, data, family=binomial)

Coefficients:
(Intercept)          Ahi          Bhi  
    -2.7317      -0.1074       1.4774

> glm(y ~ A+B-1, data, family=binomial)

Coefficients:
   Alo     Ahi     Bhi  
-2.732  -2.839   1.477
つまり，定数項を持つモデルの定数項 -2.7317 は，定数項を持たないモデルの最初のカテゴリー変数のベースライン Alo に対する係数 -2.732 になっているのです。 Blo は2つめのカテゴリー変数のベースラインなので，定数項を持たないモデルにおいても現れない。両方のモデルの Ahi の係数が違うように見えるが本当は同じになるよというのは No. 6450 に書いたとおり。

> こういったことを解説した成書ないし論文がありましたら

このようなことは実は，当たり前の単純なことなので，わざわざ解説したりしたものはないのではないかなあ？数式で示すこともできるだろうけど，数値データで明らかにわかることなので，頭をひねって証明する必要もないのでは？

No.06458　Re: ロジスティック回帰における，定数項が「ないこと」の意義　　【まだまだ修業中】　2008/05/01(Thu) 15:58

>このようなことは実は，当たり前の単純なことなので，わざわざ解説したりしたものはないのではないかなあ？

お恥ずかしい限りです。

いずれにしても，大変よくわかりました。
今回の私の疑問は解消されました。本当にありがとうございました。

（と同時に，分析過程をいかに理解していないで使っていたかがわかりちょっと恐ろしくもなりました。）

今後ともよろしくお願いいたします。