Rでロジスティック回帰分析
Rでロジスティック回帰分析をしたので、おおまかな流れをメモ。
無駄にオブジェクトや命令文が増えているが、自分にとっての分かりやすさを重視で。
分析
分析自体はデフォルトの関数(glm)を用いるが、追加でpsclパッケージを用いた。psclに含まれるpR2()で、対数尤度やMacFaddenの疑似R二乗値を計算してくれる。
他にもDesignパッケージのlrm関数というのもモデルの検定なんかに便利らしいので、今後また調べるつもり。
library(pscl) #今回必要なパッケージ #関数本体。dfはデータフレーム、yは0か1の従属変数、x1〜xnが独立変数。 ans <- glm(y ~ x1+x2+xn, family=binomial(link="logit"), df) #切片のみモデル ans0 <- glm(y ~ 1, family=binomial(link="logit"), df) summary(ans) #結果を表示 exp(ans$coef) #オッズ比 pR2(ans) #尤度比統計(詳細はヘルプファイル) anova(ans0, ans, test="Chisq") #尤度比検定(-2対数尤度でカイ二乗検定)
データを整えてファイル出力
ついでに結果をまとめるためのデータ操作も書いておく。
library(Hmisc) #ここで使うパッケージ ##データの整形 (temp <- summary(ans)$coe) #結果のスコアをデータフレームに (usb <- temp[,1:2]) #偏回帰係数と標準誤差 (or <- data.frame(exp(ans$coef))) #オッズ比 (table <- data.frame(usb, or)) #表の統合 (table <- round(table, 3)) #表を小数点以下3桁で揃える (logL <- round(pR2(ans), 3)) #尤度比統計も桁を揃える summary(ans) #有意確率など、スコアの確認用に出力
ここまでやれば、あとはExcelなりLaTeXなりで使えるように出力して、表を整える。サンプルサイズや有意確率も忘れずに。
