crosshyou

主にクロス表(分割表)分析をしようかなと思いはじめましたが、あまりクロス表の分析はできず。R言語の練習ブログになっています。

都道府県別の建築物リフォーム・リニューアル調査のデータ分析7 - リフォーム・リニューアルの金額を県民1人当たり所得で回帰分析してみる。

UnsplashRedd Fが撮影した写真 

www.crosshyou.info

の続きです。

今回は、l_total: リフォーム・リニューアルの金額の対数変換値値をl_shotoku23: 平成23年基準の県民1人当たり所得の対数変換値で回帰分析してみます。1人当たり所得の所得が大きいほど、リフォーム・リニューアルの金額が大きいかどうかを調べてみます。

lm()関数で回帰分析してみます。

summary()関数でmodel1を表示してみます。

l_total = -23.283 + 3.8645(l_shotoku23) + u  という結果です。

これは、shotokuが1%増えれば、totalが3.8645%増える、と解釈できます。

グラフで視覚化してみます。

散布図を見ると、確かに右肩上がりになっています。

このモデルに、yearを加えるとどうなるでしょうか?

update()関数でmodel1にyearを加えてみました。

summary()関数で表示してみます。

l_shotoku23の係数が4.0380になりました。
yearを入れたこのモデルは、l_shotokuの傾きは一緒で、2016年度、2017年度、2018年度と切片が違うモデルです。

これもグラフにしてみたいと思います。

これは、moderndiveパッケージのgeom_parallel_slopes()関数を使ってみます。

2016年度が一番上に回帰直線があります。

このmodel1とmodel2では統計的に有意な違いがあるのかどうか、anova()関数で確認します。

p値が0.0004271となっているので、model1とmodel2は統計的に有意な違いがあります。つまり、model1にyearを加えた意味は有った、ということです。

今度は、年度ごとに傾きが違うモデルを考えてみます。

l_shotoku23 * year と + ではなくて、* を使うと切片だけでなくて、傾きも違うモデルになります。

サマリーを表示してみます。

l_shotoku23の係数が、3.8667になっています。これは、2016年度の傾きです。

2017年度の傾きは、これにl_shotoku23:year2017年度の0.2228を足して、

3.8667+0.2228=4.0895になります。

2018年度の傾きは、3.8667+0.2995=4.1662となります。

傾きが年度が増えるごとに急になっています。

これもグラフにしてみます。geom_smotth()関数で視覚化できます。

傾きの違いは、あまり感じられないですね。

anova()関数で、model2とmodel3を比較してみます。

p値が0.9413ですから、model2とmodel3で統計的に有意な違いは無いということです。つまり、傾きを年度ごとに変える意味は無いですね。

こんどは、model2にprefを加えてみます。都道府県ごとにも切片を変えるということです。

summary()関数でmodel4をみてみます。

中略

l_shotoku23の係数が-0.09525とほとんど0に近いマイナスになりました。p値も0.957340となっているので、prefを考慮すると、l_shotoku23はl_totalとは無関係ということですね。

anova()関数で、model2とmodel4を比較してみます。

p値は2.2e-16とゼロといってよいぐらい小さい値です。moel2にprefを加えてmodel4を作成した意味はあった、ということですね。

今回の結論としては、都道府県の違いを考慮しないと、県民1人当たり所得が多いほどリフォーム・リニューアルの金額も多い。都道府県のごとに見ると、県民1人当たり所得とリフォーム・リニューアルの金額は無関係、ということでした。

今回は以上です。

初めから読むには、

www.crosshyou.info

です。

都道府県別の建築物リフォーム・リニューアル調査のデータ分析6 - Rのcaretパッケージでratioを予測するモデルを作る。その2

UnsplashLeonid Antsiferovが撮影した写真 

www.crosshyou.info

の続きです。引き続き、caretパッケージでratioを予測するモデルを作成していきます。

まずは、linear regression modelで、もう少し複雑なモデルにしてみます。

l_totalの2乗項と、l_total x year を追加してみました。

このlm_fit2と名付けたオブジェクトを表示してみます。

names()関数でどういう要素があるか見てみます。

たくさんありますね。finalModelが最終的なモデルですかね。みてみます。

このような推定結果でした。

では、このモデルを使って、df_testのデータを使い、ratioを予測してみます。

predict()関数で予測できます。

これを評価します。

この評価結果を保存しておきます。

複雑な式のほうが、予測精度は低いですね。

この調子で、他のモデルも試していきます。まずは、Random Forestモデルです。

method = "rf"にするだけです。

rf_fitをみてみます。

mtryというのがパラメータで、mtry = 2 が最終的に採用されたようです。

finalModelを表示してみます。

見方がわからないですが、500回の決定木モデルを試したってことですかね。。

予測と評価をします。

この結果も保存しておきます。

いままでの3つのモデルの中では、Random Forestが一番予測精度が高いですね。

次は、gam, Generalized Addaptive Model with using Splinesを試してみます。

gam_fitを表示してみます。

finalModelを表示してみます。

GVC scoreってなんでしょうね。わからないです。

予測して、評価します。

この評価結果も保存しておきます。

round()関数で小数点以下4桁までの表示にしました。

最後にSupport Vector Machines with Linear Kernelのモデルをやってみます。

svm_fitオブジェクトを表示してみます。

final modelを表示してみます。

パラメータには、spsilionとcost C というのがあるようです。

予測をして結果を評価します。

全部の予測結果を表示してみます。

RMSEの小ささで評価すると、Random Forestが一番予測精度が高く、2番目が単純なlinear regression modelでした。

今回は以上です。

次回は、

www.crosshyou.info

です。

 

初めから読むには、

www.crosshyou.info

です。

読書記録 - 「資本主義の方程式 経済停滞と格差拡大の謎を解く」 小野善康 著 中公新書

日本の経済が停滞し、格差拡大しているのは何故か、ということを「資産選好」という観点から解説しています。

「資産選好」というのは、人々が(たとえ将来使う予定がなくても)金融資産をもっともっと保有したい、という執着心のことです。

日銀がいくら異次元の金融緩和をしても、人々は増えたお金を消費にまわすことなく、金融資産として保有するだけなので景気は回復しない、ということです。

都道府県別の建築物リフォーム・リニューアル調査のデータ分析5 - Rのcaretパッケージでratioを予測するモデルを作る。その1

UnsplashDiego PHが撮影した写真 

www.crosshyou.info

前回は、l_totalを2016年度、2017年度、2018年度で回帰分析してみました。理論ベースでもbootstrapでのシミュレーションベースでも、年度によってl_totalの平均値には違いがあることが確認できました。

今回は、caretパッケージのツールを使って、ratioをl_total, year, prefで予測するモデルをいくつか作ってみましょう。

The caret Package (topepo.github.io) にcaretパッケージの解説があります。

まずは、caret パッケージを読み込みます。

ratio, l_total, year, prefだけのデータフレームも作っておきます。

次は、df_selectedをトレーニング用のデータとテスト用のデータに分割します。

そして、トレーニング用のデータフレーム、df_trainとテスト用のデータフレーム、df_testを作ります。

両方のサマリーを表示してみましょう。

ratioの平均値はdf_trainは0.5002, df_testは0.5169でした。

t.test関数で二つのratioの平均値に統計的な違いがないかどうかみてみます。

p-valueは0.6711ですので、統計的に違いはないです。

l_totalもみてみます。

p-valueは0.3008なのでこちらも統計的に有意な違いはないですね。

それでは、caretによる予測モデル作成をしていきましょう。

まずは、基本として一番簡単な線形モデル、ratio = l_total + year + pref を作成して評価してみます。

lm_fit1のモデルをみてみましょう。

RMSE, Rsquared, MAEがモデルがどのくらいratioの値を予測しているかを表す評価項目ですね。RMSEとMAEは小さいほうがよくて、Rsquaredは大きいほうがいいですね。

summary()関数で表示してみます。

中略

それでは、lm_fit1のモデルで、df_trainのデータを使い、ratioを予測してみます。

predict()関数を使います。

予測データがlm_pred1として保存されました。

これを、postReample()関数で評価します。

lm_eval1として評価を保存しておきました。次回以降でいろいろな予測モデルを試してみましょう。

今回は以上です。

次回は、

www.crosshyou.info

です。

 

初めから読むには、

www.crosshyou.info

です。

 

都道府県別の建築物リフォーム・リニューアル調査のデータ分析4 - inferパッケージでBootstrapでの回帰分析

UnsplashNick Fewingsが撮影した写真 

www.crosshyou.info

の続きです。前回までのデータの分析で、年度によって、l_totalの値が違うっぽいことがわかりました。今回はこれが統計定期に有意な違いなのかどうかを調べてみます。

inferパッケージで、シミュレーションベースで調べてみます。

はじめに、年度ごとにl_totalの値が違うかどうかをグラフでみてみましょう。

geom_boxplot()関数で箱ひげ図を描きました。真ん中の黒い水平線が中央値です。

group_by()関数とsummarize()関数とmean()関数で年度ごとのl_totalの平均値も確認しておきます。

2016年度は7.75, 2017年度は7.41, 2018年度は7.43です。

それでは、inferパッケージを読み込みます。

Full infer Pipeline Examples • infer の

MULTIPLE EXPLANATORY VARIABLES の例を参考にします。

まずは、サンプルからのobserved fitを取得します。

2016年度が基準で、2017年度は-0.244, 2018年度は-0.320小さいということです。

7.51(2017年度の平均) - 7.75(2016年度の平均) = -0.24

7.43(2018年度の平均) - 7.75(2016年度の平均) = -0.32

という結果と同じですね。

これは、理論ベースの lm()関数でも同じです。

続いて、Bootstrapをして1000回繰り返します。

この結果を視覚化してみます。

steelblue色の垂直線が0の位置です。2017年度の推定値も2018年度の推定値も分布は明らかに0よりも小さいエリアにあります。

この値をlm()関数の値と比較しますと、

ほとんど同じ値ですね。

95%の信頼区間を求めてみます。

inferパッケージでbootstrapをしてほうは、get_confidence_interval()関数で、

lm()関数での回帰分析のほう(理論的なほう)は、confint()関数で身体区間が求まります。

bootstrapのほうが範囲が広いですね。

2018年度の信頼区間もみてみます。

2018年度もbootstrapでの信頼区間のほうが範囲が広いです。

p値を求めてみます。

bootstrapのほうは0、lm()関数のほうは1.148e-06です。どちらも0とみなしていいくらいの小さい値です。

inferパッケージでbootstrapで回帰分析した場合は、get_p_value()関数でp値がわかります。

2018年度のp値も同じように求めてみます。

bootstrapでのp値は0,lm()関数でのp値は9.22e-10です。どちらも0とみなしていいくらいの小ささですね。

l_totalの値は年度によって値が違うようです。

今回は以上です。

次回は

www.crosshyou.info

です。

初めから読むには、

www.crosshyou.info

です。

都道府県別の建築物リフォーム・リニューアル調査のデータ分析3 - geom_point()関数で散布図を描く。top_n()関数で上位/下位を抽出する。

UnsplashZoe Schaefferが撮影した写真 

www.crosshyou.info

の続きです。前回はヒストグラムで一つの変数の分布をグラフにしてみました。

今回は、散布図で二つの変数の関係性をグラフにしてみます。

まずは、l_shotoku23とl_totalです。

l_shotoku23が大きいと、l_totalも大きいという正の相関関係が見られます。

次は、l_shotoku23とratioをみてみます。

あんまり相関はない感じです。

l_totalとratioをみてみます。

これもあまり相関は無いようです。

cor()関数で相関係数をみてみます。

l_shotoku23とl_totalの相関係数は、0.65, l_shotoku23とratioは-0.117, l_totalとratioは-0.067となっています。

l_shotoku23のtop5をみてみます。

top5といっても2016年度、2017年度、2018年度の3つの年度があるので、上位15の観測を表示しています。東京都、愛知県、栃木県、静岡県、富山県、滋賀県、茨城県が表示されました。top_n()関数で上位15の観測だけにして、arrange(desc())関数で大きい順に表示しています。

bottom5も表示してみます。

沖縄県、鳥取県、宮崎県、鹿児島県、奈良県、青森県が表示されました。top_n()関数で、top_n(-15, l_shotoku23)のようにマイナスの数字を入れると、小さいほうの観測をとってきます。

l_totalのtop5をみてみます。

東京都、神奈川県、愛知県、大阪府、北海道、埼玉県がランクインしました。

bottom5をみてみます。

高知県、沖縄県、奈良県、島根県、佐賀県、大分県、鳥取県、和歌山県、山梨県が表示されました。

ratioのtop5を見てみます。

島根県、岡山県、岩手県の3県がratioが1.0以上のときがありました。この3県と宮崎県、福井県、大分県、兵庫県、埼玉県、神奈川県、鳥取県、和歌山県、福岡県が表示されました。

ratioのbottom5をみてみます。

徳島県、沖縄県、佐賀県、山梨県、栃木県、静岡県、群馬県、福井県、茨城県、福島県、香川県が表示されました。これらの県は、住宅のリフォーム・リニューアルの金額よりも、非住宅のリフォーム・リニューアルの金額のほうがうんと多いということですね。

今回は以上です。

次回は、

www.crosshyou.info

 

初めから読むには、

www.crosshyou.info

です。

都道府県別の建築物リフォーム・リニューアル調査のデータ分析2 - データをヒストグラムにして分布を可視化する。

UnsplashMarek Piwnickiが撮影した写真 

www.crosshyou.info

の続きです。前回はデータをRに読み込ませるところまでやりました。

そして、shotoku17: 平成17年基準の県民1人当たり所得のデータはすべてNAであることがわかりました。なので、shotoku17は削除しておきます。

今回はデータをグラフにして、どんな様子のデータなのかをみてみます。

はじめは、total: リフォーム・リニューアルの合計金額ですね。

ggplotパッケージを使ってヒストグラムにしてみます。

geom_histogram()でヒストグラムにしました。さらに、facet_wrap()で年度ごとにしてみました。右に外れ値のように大きな値の観測があります。

totalは、対数変換したほうが正規分布に近くなりそうですね。

こんどは、このl_totalでヒストグラムにしてみます。

こちらのほうが分布の形状がわかりやすいですね。2016年度、2017年度、2018年度とだんだんと左側に分布が移動しているような気がします。

次は、house: 住宅のリフォーム・リニューアルの金額です。これも対数変換してからヒストグラムにしてみます。

住宅はトータルほどには右側の外れ値は目立たないですね。こちらもだんだんと左に分布が動いているように思われます。

次は、nohouse: 住宅以外のリフォーム・リニューアルの金額のヒストグラムです。前の二つと同じく対数変換値を作成してヒストグラムにしてみます。

住宅以外のリフォーム・リニューアルの金額の分布も左方向に移動、つまり金額が減少しているように見えます。

続いて、shotoku23: 平成23年基準の1人当たり所得をヒストグラムにしてみます。これも対数変換値を作成してからヒストグラムにしてみます。

1人当たり所得は対数変換値を作成しても右側に外れ値がありますね。

houseとnohouseの比率を作成して、それのヒストグラムもみてみます。

house/nohouseの比率をratioと名付けました。ヒストグラムを見ると、1以上の観測値もあります。住宅のリフォーム・リニューアルの金額ののほうが多い都道府県があるということですね。
今回は以上です。

次回は、

www.crosshyou.info

です。

初めから読むには、

www.crosshyou.info

 

です。