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

No.11549　logistic regressionのcrude odds ratioについて　　【taipapa】　2009/12/18(Fri) 01:57

お世話になります．
ある手術80例の結果に影響する因子をlogistic regression modelで解析しようとしています．初歩的な質問で恐縮ですが，crude odds ratioについてご教示ください．

library(epicalc) [Rのepicalcというパッケージを使用してます]
glm.hogehoge <- glm(Y ~ x1+x2+x3+x4+....x10, binomial)
logistic.display(glm.hogehoge, crude.p.value=TRUE)
[logistic.displayはepicalcに入っている関数です]
Logistic regression predicting Y : positive vs negative 
 
                        crude OR(95%CI)      crude P value
x1                       1.03 (0.96,1.12)     0.396        
                                                                       
x2                       0.76 (0.13,4.27)     0.754        
                                                                       
x3                       1.5 (0.5,4.49)       0.464      

x4                       6.95 (0.86,56.23)    0.069

このあとに，asjusted odds ratioが続くわけですが，このcrude P valueは，glm(Y ~ x4, binomial)とした時のWald testのp値と一致するようです．しかし，尤度比検定のp値とは一致しません．私の質問は，crude odds ratioとしてプレゼンするときには，このように，単変量でやったときのWald testのp値を言えばよいのでしょうか？また，尤度比検定のp値はcrude p valueとしては使えないのでしょうか？というのは，ある変数が，単変量解析の時にcrude p valueとして，Wald testでは0.069なのに，尤度比検定では0.02と有意になっているからです．このような場合の解釈をご教示いただければ幸いです．

No.11552　Re: logistic regressionのcrude odds ratioについて　　【青木繁伸】　2009/12/18(Fri) 09:24

> ある変数が，単変量解析の時にcrude p valueとして，Wald testでは0.069なのに，尤度比検定では0.02と有意になっているからです．このような場合の解釈をご教示いただければ幸いです．

単変量の結果は単変量の結果，多変量の結果は多変量の結果であり，両者が一致する必要もないし，一致する必然性もありません。別物なのです。

No.11553　Re: logistic regressionのcrude odds ratioについて　　【taipapa】　2009/12/18(Fri) 16:56

青木先生，ご指摘有り難うございます．
私の書き方が良くなかったようです．ポイントを絞って書きますと，
glm(Y ~ x, binomial)としたときのcrude odds ratioのWald testのp値と尤度比検定のp値が異なる場合は，どちらを採用すべきでしょうか？　今回の解析の中で，Wald testで有意でないのに，尤度比検定で有意と表示されたものがあったので，お聞きしています．あくまで，単変量解析内部での話です．尤度比検定の方が良いと読んだことがあるのですが．．．

No.11554　Re: logistic regressionのcrude odds ratioについて　　【青木繁伸】　2009/12/18(Fri) 18:26

> Wald testで有意でないのに，尤度比検定で有意

尤度比検定ってどういう風にやるんですか？

> library(epicalc)
> d <- data.frame(y=factor(sample(0:1, 20, replace=TRUE)), x=factor(sample(0:1, 20, replace=TRUE)))
> ans <- glm(y~x, d, family=binomial)
> summary(ans)

Call:
glm(formula = y ~ x, family = binomial, data = d)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.1010  -1.1010  -0.9005   1.2558   1.4823  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)
(Intercept)  -0.6931     0.7071  -0.980    0.327
x1            0.5108     0.9309   0.549    0.583

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 26.920  on 19  degrees of freedom
Residual deviance: 26.615  on 18  degrees of freedom
AIC: 30.615

Number of Fisher Scoring iterations: 4

> logistic.display(ans, crude.p.value=TRUE)

Logistic regression predicting y : 1 vs 0 
 
          OR(95%CI)          P(Wald's test) P(LR-test)
x: 1 vs 0 1.67 (0.27,10.33)  0.583          0.581        ★★ ここの2つのP値の違いのこと？
                                                      
Log-likelihood = -13.3077
No. of observations = 20
AIC value = 30.6155

> anova(ans, test="Chisq")
Analysis of Deviance Table

Model: binomial, link: logit

Response: y

Terms added sequentially (first to last)

     Df Deviance Resid. Df Resid. Dev P(>|Chi|)
NULL                    19     26.921          
x     1  0.30501        18     26.616    0.5808

> dput(d) # 使ったデータフレーム
structure(list(y = structure(c(1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
1L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 1L), .Label = c("0", 
"1"), class = "factor"), x = structure(c(2L, 2L, 1L, 1L, 2L, 
1L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 2L), .Label = c("0", 
"1"), class = "factor")), .Names = c("y", "x"), row.names = c(NA, 
-20L), class = "data.frame")

No.11555　Re: logistic regressionのcrude odds ratioについて　　【青木繁伸】　2009/12/18(Fri) 18:39

オッズ比の信頼区間や検定は，何通りかありますし，anova にしろ Wald にしろ，漸近近似ですし。fisher.test で出てくるオッズ比の信頼区間は広いし（P値は大きめに出るし）。
No. 11554 の d だと，以下のようになるわけで。
> table(d)
   x
y   0 1
  0 6 6
  1 3 5
> fisher.test(table(d))

	Fisher's Exact Test for Count Data

data:  table(d) 
p-value = 0.6699
alternative hypothesis: true odds ratio is not equal to 1 
95 percent confidence interval:
  0.1969305 15.6493164 
sample estimates:
odds ratio 
  1.624460 

No.11556　Re: logistic regressionのcrude odds ratioについて　　【taipapa】　2009/12/19(Sat) 00:31

青木先生，何度もありがとうございます．

> logistic.display(ans, crude.p.value=TRUE)

> Logistic regression predicting y : 1 vs 0

> OR(95%CI) P(Wald's test) P(LR-test)
> x: 1 vs 0 1.67 (0.27,10.33) 0.583 0.581 ★★ ここの2つのP値の違いのこと？

その通りです．LR-testはlikelihood ratio testと考え尤度比検定かと．そして，ご指摘のとおり，anova(ans, test="Chisq")の結果と等しいです．
また，No. 11555で仰るように，odds ratioや信頼区間はfisher.testの値とも異なります．

以前出来合いのソフトでやっているときは，出てきた数字をそのまま何も考えずに使っていたのですが，Rを使うようになってからは，色々な方法が使え，情報もネットにいっぱいあって良いのですが，このcrude odds ratioのような単純そうなことでも，結局，種々の方法の中のどれが一番適切なのか難しくて，私のような統計の素人にはどれを使用すべきか迷ってしまいます．

つまりは，No. 11554でいうと，

OR(95%CI) P(Wald's test) P(LR-test)
x: 1 vs 0 1.67 (0.27,10.33) 0.583 0.581

0.583 を使うべきか0.581を使うべきかと言うことです．この場合は，どっちも有意ではないので，どっちでも良いのでしょうが，私のデータでは，Wald testでは0.069なのに，P(LR-test)では0.02と有意になっているので困っております．まぁ，値としてはたいした差ではないかもしれませんが．はたして，どちらを使うべきなのでしょう？おそらく大部分の論文はWald testの結果を使用していると推測しているのですが．（大部分のソフト（SPSSなど）はデフォルトはWaldをつかっているようです）

医学論文の場合は，まず，単変量解析で粗オッズ比と信頼区間とp値を見せてから，多変量解析で，model reductionをやって最適モデルで同じことをするのが常套と思うのですが，みなさん，単変量解析の時の検定には上記の何を使われているのでしょうか？

No.11557　Re: logistic regressionのcrude odds ratioについて　　【taipapa】　2009/12/19(Sat) 00:35

> 医学論文の場合は，まず，単変量解析で粗オッズ比と信頼区間とp値を見せてから，多変量解析で，model reductionをやって最適モデルで同じことをするのが常套と思うのですが，みなさん，単変量解析の時の検定には上記の何を使われているのでしょうか？

多変量解析の時は，調整オッズ比と調整信頼区間とそのp値ですね．訂正しておきます．

No.11558　Re: logistic regressionのcrude odds ratioについて　　【知ったかぶり】　2009/12/19(Sat) 14:39

Waldではダメというわけではないけれど,尤度比検定の方が良い,とされているようです.Wald testが広く使われているのは,計算が簡単だからでしょう.

http://www.ats.ucla.edu/stat/mult_pkg/faq/general/nested_tests.htm
http://en.wikipedia.org/wiki/Wald_test

No.11559　Re: logistic regressionのcrude odds ratioについて　　【taipapa】　2009/12/19(Sat) 22:13

「知ったかぶり」さん，情報ありがとうございます．やはり，尤度比検定が良いような記述ではありますね．ただ，最終的には，調整オッズ比が重要でしょうから，とりあえず，epicalcのlogistic.displayが出力するWald testの値を使う方が楽ですね．