

Rでロジスティック回帰を行った際のresidual devianceを取得したいと考えています。  
glm.fit関数を用いたところ、以下のような完全分離の警告が出ました。

[code]

    glm.fit: fitted probabilities numerically 0 or 1 occurred 
    
[/code]

そこで、[こちら](https://www.bananarian.net/entry/logistic)を参考にlogistf関数を利用してみました。

[code]

    cont <- logistf.control(maxit=1000,maxstep = 100)
    plcont <- logistpl.control(maxit=1000)
    res <- logistf(y ~ x,firth = TRUE, control = cont,plcontrol = plcont)
    
[/code]

summary(res)の中身は以下のとおりです。

[code]

    > summary(res)
    logistf(formula = dat[, j] ~ Dat_res[, i], control = cont, plcontrol = plcont, firth = TRUE)
    
    Model fitted by Penalized ML
    Coefficients:
                       coef   se(coef)  lower 0.95 upper 0.95    Chisq        p
    (Intercept)   0.1375199 0.11492692 -0.08696857  0.3642894 1.439322 0.230249
    Dat_res[, i] -0.7204668 0.05044381 -0.82565105 -0.6269034      Inf 0.000000
                 method
    (Intercept)       2
    Dat_res[, i]      2
    
    Method: 1-Wald, 2-Profile penalized log-likelihood, 3-None
    
    Likelihood ratio test=895.6586 on 1 df, p=0, n=1000
    Wald test = 203.9916 on 1 df, p = 0
    
[/code]

この結果から、[こちらのサイト](https://www.slideshare.net/logics-of-
blue/2-4deviance)でも説明されているdeviance残差(residual
deviance)を計算しようと考えていますが、導出する方法がよく分かっていません。  
理解不足で申し訳ありませんが、ご教授よろしくお願いします。


CC BY-SA 4.0
2021-08-26T07:25:36.690
81058
2021-10-20T06:02:40.667
47922
question
r
Rでロジスティック回帰を行ったときのdeviance残差の計算方法


一般化線形モデル (GLM) の deviance (逸脱度) の定義は「-2 x モデルの対数尤度」です
(別の定義もありますがそのスライドや質問の文脈からこれを意図していると判断します), この場合, 通常 residual deviance
は推定したモデルと飽和モデル (フルモデル) の逸脱度の差のことを指します.

この場合, `logisf` のオブジェクトの `loglik` には罰則推定モデルとフルモデルの対数尤度が含まれているので, `-2 *
(res$loglik[1] - res$loglik[2])` で計算できることになります.

そしてこれはデフォルト設定での尤度比検定の検定統計量と等価なので, `summary()` で表示したときの `Likelihood ratio
test=...` の値と同じになるはずです.

おせっかいかも知れませんが, ネット上の記述は多くの場合断片的なので, 体系的に説明している教科書を1冊は読んでおいたほうがいいと思います.
例えば日本語で, かつ逸脱度について言及がある本の1つに, 久保拓哉『データ解析のための統計モデリング入門』 (2012, 岩波書店) があります.

ただし, 免責事項というか, 以下の点に注意してください.

  1. 私は **logistf** パッケージを今回初めて知ったので, 適切に動作するかは確認していません (説明書通りの方法を提示しただけです)
  2. そちらが何を目的にしているのかはこちらにはわからないため, このアプローチが適切なのかまでは保証できません (i.e., 罰則付き推定と罰則なし推定は別物ですし, 完全分離がどのような理由で発生しているかは私には知りえません)

**回答から時間が経ってしまいましたが, 誤解を招きそうなので以下追記**

質問文とコメントを読み返して気づきましたが, 私もあなたも residual deviance と deviance residuals
を混同して書いています. 私の方は「deviance残差」と表記してる箇所もありますが, 意図していたのは全て residual deviance です.
モデルに対して定義されるのがこの residual deviance で, **-2*(対数尤度 - フルモデルの対数尤度)** と定義され,
これは上記通りの方法で数値を得られます.

しかし, 参考にされたスライドでは, 上記の residual deviance を単に "deviance" と呼んでいます. deviance を
**-2*(対数尤度)** とするか **-2*(対数尤度 - フルモデルの対数尤度)** と定義するかは文献によってまちまちなので紛らわしいです. R の
glm のヘルプや上記の参考書では, 参考になさったスライドとは異なり, 前者を deviance とし, 後者を residual deviance
と呼んでいます. (以降はこの区別のため, 元スライドと違う呼び方になってることに注意してください).

一方で観測ごとに計算する deviance residuals / deviance 残差 は, 上記の (residual) deviance
の平方根に符号をつけたものを指します. deviance = **-2*(対数尤度)** と residula deviance = **-2*(対数尤度
- フルモデルの対数尤度)** のどちらを使うかもまた, 文献によってまちまちなので紛らわしいです.
(参照したスライドも式を明示していないので私も読み違えていたようです) つまり,

  1. deviance の定義には複数あり, しかもどれが主流か決まってない
  2. **(モデルに対する) residual deviance と (観測点に対する) deviance residuals 名前が似ているが全く別物である**

という二重の紛らわしい箇所があることに注意してください (後者は読み違えないように "unit --" とか "casewise ---"
とかを頭につけて区別することがあります).

そちらはタイトルでは「deviance残差」と書いて本文では「residual deviance」「deviance残差(residual
deviance)」などと書いています. 私は, 表記にゆらぎがあるものの, モデルに対する residual deviance
を意図していると判断してして回答しましたが, 実際に欲しいのが residual deviance ではなくこちらの (casewise な)
deviance residuals であるなら, glm では `residual(<glm オブジェクト>)` で得られます.
ただしスライドの定義とは異なり, R では **-2*(対数尤度)** を deviance と定義しており, この関数で計算されるものもこれに準じています
(GLM で使われる residuals の定義は他にもいくつもありますが, どれを使うべきかは目的によると思います). 一方で logistf には
residuals を計算する機能はありません. よって観測点ごとの対数尤度か事後対数確率から計算することになります.
例えば以下に計算例が挙げられています. これは尤度を基にした場合で, かつ deviance は **-2*(対数尤度)** と定義しているようです.

<https://github.com/georgheinze/logistf/issues/37>

なお, コメントに書いたように, 線型回帰モデルの残差とGLMのdeviance 残差は必ずしも同じふるまいをするわけではないため,
常に残差と同じように扱えると思い込まないでほしい, という点は変わりません.


CC BY-SA 4.0
2021-08-29T15:16:41.863
81103
2021-10-20T06:02:40.667
2021-10-20T06:02:40.667
40575
40575
81058
answer