1 例題

バーセル指数(barthtot: 移乗動作や着替え、食事といった日常生活動作の自立度)の高さ/低さが、どの程度の長さの介護時間の長さ(c12hour)と結びついているのかを調べたい。ただし、性別(c161sex)を統制変数とする。この測定結果に基づいて、バーセル指数の限界効果(marginal effect)をプロットする。限界効果とは、他の変数を一定(通常は平均値で固定)とした場合、注目する独立変数の値が1単位増加した場合の従属変数に与える効果のことである。

1.1 データを読み込む

# パッケージを読み込む
library(sjPlot)

# データを読み込む
data(efc)

1.2 モデルをあてはめる

fit1 <- lm(c12hour ~ barthtot + factor(c161sex), data = efc)

# 結果を見てみる
summary(fit1)
## 
## Call:
## lm(formula = c12hour ~ barthtot + factor(c161sex), data = efc)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -88.519 -22.919  -9.257   6.466 155.019 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      93.88865    4.52285  20.759   <2e-16 ***
## barthtot         -0.85539    0.05018 -17.046   <2e-16 ***
## factor(c161sex)2  4.63045    3.47836   1.331    0.183    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 43.95 on 879 degrees of freedom
##   ( 26 個の観測値が欠損のため削除されました )
## Multiple R-squared:  0.2512, Adjusted R-squared:  0.2495 
## F-statistic: 147.5 on 2 and 879 DF,  p-value: < 2.2e-16

1.3 プロット

# ライブラリを読み込む
library(ggplot2)
library(sjPlot)

# プロット
theme_set(theme_bw())  # ggplot2の白黒テンプレートを使用(指定しない場合にはデフォルトのテンプレートが使用される)

plot_model(
  model = fit1,                # あてはめた結果のobject  
  type = "pred",               # 回帰分析の予測値に基づくプロットを指定
  terms = "barthtot",      # x軸にどの変数を指定するか
  title = "バーセル時間から予測される介護時間", 
  axis.title = c("バーセル指数", "平均介護時間")) #x軸名、y軸名の順に

参考文献 Lüdecke, Daniel. 2023. “Plotting Marginal Effects of Regression Models”. https://strengejacke.github.io/sjPlot/articles/plot_marginal_effects.html