Rで何かをしたり、読書をするブログ

政府統計の総合窓口のデータや、OECDやUCIやのデータを使って、Rの練習をしています。ときどき、読書記録も載せています。

都道府県別の定期健康診断結果報告のデータの分析5 - 所見率を被説明変数にして重回帰分析をする。

www.crosshyou.info

の続きです。今回は、shokenritsu: 所見率を被説明変数にして重回帰分析をしてみたいと思います。

まずはじめに、per_jushin: 1事業所当たりの受診者数の数、log_place: 事業所の数の対数変換値、log_jushin: 受信者数の対数変換値を説明変数にして重回帰分析をしてみました。

p-value: 1.268e-05 とあります。これは、0.00001268 ですので、統計学的には有意なモデルです。ただ、各変数の p値は0.05よりも大きいです。R-squared が0.1978ですから、約20%くらいしか shokenritsu を説明できません。

このモデルに2乗項と交差項を加えてみます。

log_place の2乗項と log_jushin の2乗項の p値が0.05よりも小さくなりました。R-squareは0.2557と少し上昇しました。

step() 関数で不要な説明変数を削除してみます。

per_jushin: log_place: log_jushon の p値が 0.57626 ですから、これも削除できそうですね。他にも0.05以上のものがあるので、update() 関数で手動で削除していきます。

まだまだ削除できます。

I(per_jushin^2) を削除します。

これで、モデルの中の全ての変数の p値が0.05よりも小さな値になりました。R-squaredは0.2242です。このモデルで shokenritsu を 22% くらいは説明できるということですね。これだけ頑張っても 22% ですから、shokenritsu は 他の要因の影響が大きいということですね。

この最終的なモデル、lm_mod4v5 の予測値と実際の shokenritsu の散布図を描いてみます。

横軸がモデルの予測値、縦軸が実際の shokenritsu です。

今度は、このモデルに year と pref を加えます。

R-squared が 0.9654 に跳ね上がりました!都道府県と年と説明変数に加えると一気にモデルの説明力が増しました。

予測と実際の値の散布図を描いてみます。

さきほどの散布図と比べると、予測精度が改善されていることがわかります。

year と pref だけが説明変数のモデルも試してみます。

R-squared は 0.9642 ですので、ほとんど変わらないです。散布図を描きます。

散布図もほとんど変わらないですね。

今回は以上です。

次回は、

www.crosshyou.info

です。

 

はじめから読むには、

www.crosshyou.info

です。

今回のコードは以下になります。

#
# shokenritsuをper_jushin, log_place, log_jushinで回帰分析
lm_mod3 <- lm(shokenritsu ~ per_jushin + log_place + log_jushin, data = df)
summary(lm_mod3)
#
# modelに2乗項と交差項を加える
lm_mod4 <- lm(shokenritsu ~ per_jushin * log_place * log_jushin +
                I(per_jushin^2) + I(log_place^2) + I(log_jushin^2),
              data = df)
summary(lm_mod4)
#
# lm_mod4 から不要な変数を削除
lm_mod4v2 <- step(lm_mod4)
summary(lm_mod4v2)
#
# 手動で不要な変数を削除 1回目
lm_mod4v3 <- update(lm_mod4v2, . ~ . - per_jushin:log_place:log_jushin)
summary(lm_mod4v3)
#
# 手動で不要な変数を削除 2回目
lm_mod4v4 <- update(lm_mod4v3, . ~ . - per_jushin:log_jushin)
summary(lm_mod4v4)
#
# 手動で不要な変数を削除 3回目
lm_mod4v5 <- update(lm_mod4v4, . ~ . - I(per_jushin^2))
summary(lm_mod4v5)
#
# 実際のshokenritsuとlm_mod4v5の予測値の散布図
tibble(
  actual = df$shokenritsu,
  estimate = predict(lm_mod4v5)
) |> 
  ggplot(aes(x = estimate, y = actual)) +
  geom_point() +
  geom_abline(slope = 1, intercept = 0, color = "red") +
  theme_minimal()
#
# lm_mod4v5にyearとprefを追加する
lm_mod5 <- update(lm_mod4v5, ~ . + year + pref, data = df)
summary(lm_mod5)$r.squared
#
# 実際のshokenritsuとlm_mod5の予測値の散布図
tibble(
  actual = df$shokenritsu,
  estimate = predict(lm_mod5)
) |> 
  ggplot(aes(x = estimate, y = actual)) +
  geom_point() +
  geom_abline(slope = 1, intercept = 0, color = "red") +
  theme_minimal()
#
# shokenritsuをyearとprefだけで回帰分析
lm_mod6 <- lm(shokenritsu ~ year + pref, data = df)
summary(lm_mod6)$r.squared
#
# 実際のshokenritsuとlm_mod6の予測値の散布図
tibble(
  actual = df$shokenritsu,
  estimate = predict(lm_mod6)
) |> 
  ggplot(aes(x = estimate, y = actual)) +
  geom_point() +
  geom_abline(slope = 1, intercept = 0, color = "red") +
  theme_minimal()
#

(冒頭の画像は、Bing Image Creator で生成しました。プロンプトは、Natural landscape, there are white, red and pink Rhododendron transiens all over the field, under the blue sky, photo です。)