ロジスティック回帰 (1)

例題
リハビリを受けている入院中の患者に対してリハビリに関する満足度調査を実施した。1週間に割り当てられた単位数をもとにランダムにインタビューした。「満足」と回答する割合と1週間のリハビリ単位数には関係があるのでしょうか？

僕が作った架空のデータで練習しましょう。単位数は5単位毎に設定してます。例題の内容や結果については色々とご意見があると思いますが、統計学の練習のためのサンプルですのでご了承ください。

データの準備と要約

データセット log01 をdat に読み込みます（ファイルの読み込み方）

dat <- read.csv("log01.csv", header=T, fileEncoding = "UTF-8")

データの確認

dat

（満足＝満足と回答した人数）

y軸を「満足と回答する確率」、ｘ軸を「単位数」とします

g1 <- function(){
    plot(満足/患者数~単位数, data=dat, ylab="満足と回答する確率")
}
g1()

そのまま単回帰分析を実行して、回帰直線を挿入してみます

reg <- lm(満足/患者数~単位数, data=dat)
g1()
abline(reg) 

#lines(range(dat$単位数), reg$coef[1] + reg$coef[2]*range(dat$単位数)) と同じ直線

詳細は単回帰分析をご参照ください

患者が満足と答える確率は0～1となるべきなのですが、回帰直線では0と1を突き抜けてしまいます

この直線では満足度とリハビリ単位数の関係が十分に説明できません

Rを使用したロジスティック回帰分析

目的変数が連続変数（誤差が等分散正規分布）ではないので、一般化線形モデルと呼ばれています。特に目的変数が0または1の二値の場合は、二項分布に基づくロジスティック回帰モデルを使用します。Rの関数は 一般線形モデル＝lm(linear model) 、一般化線形モデル＝ glm(generalized linear model) となります・・・少しややこしいですね。

Rでは以下のように書きます

fit <- glm(
    cbind(満足, 患者数-満足) ~ 単位数, 
    family=binomial, 
    data=dat
)

結果のまえに、この式のcbindの中を確認してみましょう

cbind(dat$満足, dat$患者数-dat$満足)

１単位から9単位までの満足と非満足の数が出力されました

$\dfrac{1列目}{2列目}$ がオッズ比になっているのが分かります

この表から感度・特異度が算出されるので、ROC曲線を描くことができることも納得できますね

結果

summary(fit)

ロジスティック回帰分析の結果、リハ単位数が増えると「満足」と回答する確率 $P$ が有意に上昇することが分かりました

リハビリ単位数が1週間に5単位上がると、満足と回答する割合のオッズが有意に$e^{0.1185}\fallingdotseq1.13$倍に上がる（13%アップ）ことを意味しています（理由は「結果の解釈」に記載してます）

ロジスティック回帰の基準(reference) 変更

リハビリテーション統計学でも適用となるロジスティック解析を行うときには、オッズとオッズ比の基準 (reference) が非常に重要になります。

y2pt.com

2023.03.27

結果のみでよければ、ここまででOKです