今回は、都道府県別の総人口・可住地面積・総生産額のデータを分析したいと思います。これらのデータで、都道府県別の人口密度・一人当りの生産額・面積当りの生産額を算出しようと思います。
データは政府統計の総合窓口、e-Statから取得しました。
47都道府県を選択し、
データ項目を選択し、
エクセルファイルにデータをダウンロードしました。
このCSVファイルをR言語のread.csv関数で読み込んで分析をしたいと思います。
read.csv関数でデータを読込み、na.omit関数でNA行を削除して、summary関数で各変数の統計値を表示しました。10年間のデータがあるようですね。何年度からあるでしょうか？table関数でみてみます。
2006年度から2015年度までの10年間ですね。
まず、人口密度を計算してみます。Pop / Ha で計算できます。
df$Pop / df$Ha で人口密度を計算して、名前属性を都道府県にしています。
2015年度の人口密度を小さい順に表示してみましょう。sort関数を使います。
北海道が一番人口密度が低く、東京が高いです。barplot関数で棒グラフを描いてみます。
棒グラフの右の三つ、神奈川県、大阪府、東京都は他と比べると段違いという感じですね。
一人当り生産額を計算しましょう。
一人当りの生産額も2015年度の値を小さい順に表示してみます。
奈良県が一番低く、東京都が高いです。boxplot関数で箱ひげ図を描いてみます。
上方の外れ値が２つあります。東京都と愛知県ですね。
面積当りの生産額を計算してみましょう。
2015年度の面積当りの生産額を小さい順に表示してみます。
北海道が一番低く、東京都が一番高いです。
10年間のデータを使って、hist関数でヒストグラムを描いてみます。
hist関数で、prob = TRUEにして相対度数のヒストグラムを描いて、lines(density(ProdHa), col = "red")で密度関数を付け足しています。750ぐらいのところ、300ぐらいのところ、220ぐらいのところ、130ぐらいのところに山があります。東京都、大阪府、神奈川県、愛知県でしょうね。
まとめると人口密度は
下位３地域は、北海道、秋田県、岩手県
上位３地域は、東京都、大阪府、神奈川県
一人当り生産額は
下位３地域は、奈良県、沖縄県、鳥取県
上位３地域は、東京都、愛知県、静岡県
面積当り生産額は
下位３地域は、北海道、秋田県、岩手県
上位３地域は、東京都、大阪府、神奈川県
でした。
一人当りの生産額を人口密度と年度で回帰してみましょう。
まず、年度をファクタでなく、数値型の変数にしましょう。
Yearのファクタ水準を1975から2017にして、ファクタ型を文字列にしてから数値型にしています。ファクタ型から直接as.numeric関数で数値型にすると、うまく年を変換できないので。
これで準備が整いましたので、lm関数で回帰モデルを作成してみましょう。
モデル全体のp-valueは2.2e-16より小さいですので有意ですが、各係数のp値は有意ではないですね。Year:Mitsudoは削除したモデルを考えます。
モデル全体のp値は、2.2e-16より小さいので有意です。modelとmodel2で有意な違いがあるかどうかanova関数で確認します。
p値が0.1669と0.05よりも大きいので、modelとmodel2には有意な違いはありません。よって、より単純なほうのmodel2を採用します。
model2を見ると、Yearの係数のp値は0.839とかなり大きいです。Yearを削除したモデルを考えます。
モデル全体のp値は2.2e-16なので有意です。model2とmodel3をanova関数で比較します。
p値は0.8387なのでmodel2とmodel3に有意な違いは無いです。
よって単純なmodel3を採用します。
model3のMitsudoの係数の符号は正の符号ですので、人口密度が上がるほど、一人当りの生産額は増える、ということですね。
plot関数で人口密度と一人当り生産額の散布図を描いて、abline関数でモデルの回帰式を追加しましょう。
このグラフを見ると、Mitsudoの2乗項を加えたほうがよさそうな気がします。
2乗項を加えたモデルを考えましょう。
モデル全体のp値は2.2e-16なので有意ですね。Mitsudoの2乗項の係数も、Mitsudoの係数も有意ですね。anova関数で、model3とmodel4を比較します。
p値が2.2e-16よりも小さいので、model3とmodel4では有意な違いがあります。model3の調製済みR2は、0.3818です。model4の調整済みR2は0.5131です。model4のほうが良いモデルですね。
それでは、MitsudoとperProの散布図にmodel4で推定された値を加えましょう。
predict関数で各xvの値に対してのyvを計算します。xvがMitsudo, yvがperProに相当します。
赤い直線がmodel3の回帰線で、青い曲線がmodel4の回帰線です。青い曲線のほうがフィットしていることがわかります。
今回は以上です。
