(Bing Image Creator で生成: プロンプト: Close up of pink and white Paeonia suffruticosa flowers, background is green grass field and blue sky, photo)
前回の回帰分析では、15歳~64歳人口割合、１人当たり最終エネルギー消費の２つの変数ともに1人当たり県民総所得と関連性があることがわかりました。今回は、2020年度のデータに2011年度のデータも加えて分析をしてみたいです。同じ観察対象の２つ以上の時間軸のデータなので、パネルデータ分析となります。
まず、2020年と2011年のデータフレームを作成します。
year というダミー変数を加えたモデルを回帰分析します。
year2011年度の係数が -494.559 ということは、他の条件が同じなら、2011年度の income は 2020年度の income よりも 50万円ほど少ない、ということですね。
summary() 関数でモデルサマリーを確認しましょう。
popratio, tokyo, energy, year2011年度、どの変数も p値は0.05 以下です。
plot() 関数でモデルのグラフを描きます。
左側の２つのグラフを見ると、Heteroskedasticity ではないようです。
Preush-Pagan 検定で Heteroskedasticity のチェックをします。
p-value は 0.5737 と 0.05 よりもかなり大きな値です。このモデルが Homoskedasticity である、という帰無仮説を棄却できません。つまり、このモデルが Heteroskedasticity ではないということです。
なので、 summary() 関数での p値は信用していいですね。
残差の2乗を回帰分析する方法でも確認します。
p-value: 0.5874 となってます。
実際の income の値と model4 の推計値をグラフにしてみます。
今回は２の期間のデータ、パネルデータを使って、年度のダミー変数を説明変数に加えて回帰分析しました。結果は、popratio の係数が 60.8 で energy の係数が 1.44 でした。
energy, tokyo, year が変わらなくて、popratio が 1 ポイント上昇すると、income はだいたい 6 万円増えます。
popratio, tokyo, year が変わらなくて, energy が 1 ギガ・ジュール上昇すると、 income は 1440 円増えます。
今回は以上です。
次回は、
です。
初めから読むには、
です。
今回のコードは以下になります。
#
# 2011年と2020年のデータフレームを作成
df_1120 <- df |>
filter(ycode == 2011 | ycode == 2020) |>
mutate(tokyo = if_else(pref == "東京都", 1, 0))
summary(df_1120)
#
# year を加えて回帰分析
model4 <- lm(income ~ popratio + tokyo + energy + year,
data = df_1120)
model4
#
# model4 のサマリー
summary(model4)
#
# model4 のグラフ
par(mfrow = c(2, 2))
plot(model4)
par(mfrow = c(1, 1))
#
# Breush-Pagan Test
lmtest::bptest(model4)
#
# model4 の残差との回帰分析
lm(I(model4$residuals^2) ~ popratio + tokyo + energy + year,
data = df_1120) |>
summary()
#
# 実際の income とmodel4 の予測値
tibble(actual_income = df_1120$income,
estimated_income = model4$fitted.values,
pref = df_1120$pref,
year = df_1120$year) |>
ggplot(aes(x = actual_income, y = estimated_income)) +
geom_point(aes(color = year)) +
geom_text(aes(label = pref), vjust = 1.1, hjust = 0.5, alpha = 0.2) +
geom_abline(intercept = 0, slope = 1) +
theme_bw()
#
