www.crosshyou.info

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

日銀の短観データの分析2 - R言語でANOVA。最近と先行きに違いはある。先行きのほうが悲観的。

 

www.crosshyou.info

 の続きです。今回もANOVAです。explanatoly variables(説明変数)をType(最近、先行きの2値をとるカテゴリカル変数)にresponse variable(反応変数)はValue(短観の値)です。

まずは、前回と同じくグラフで見てみます。

f:id:cross_hyou:20191009081425j:plain

f:id:cross_hyou:20191009081509j:plain

黒い点が最近、赤い点が先行きです。黒の水平線が最近の平均値、赤の水平線が先行きの平均値です。

tapply関数でそれぞれの平均値を求めてみます。

f:id:cross_hyou:20191009081956j:plain

先行きの平均値は、最近の半分以下ですね。

では、SSY、the sum of the squares of the differences between the y values and the overall meanを計算します。

f:id:cross_hyou:20191009085903j:plain

SSYは68298です。

SSE、the error sum of squaresを計算します。最近の平均値、先行きの平均値はさきほど計算していますから、これを利用しましょう。

f:id:cross_hyou:20191009090653j:plain

SSEは65572です。

SSA、the treatment sum of squaresはSSY - SSEで計算できます。

f:id:cross_hyou:20191009090948j:plain

SSAは2726です。SSA, SSE, SSYがわかりましたので、ANOVA表を作ってみます。

f:id:cross_hyou:20191009091802j:plain

D3のセルのところ、SSEのMean Squareを計算します。

f:id:cross_hyou:20191009091956j:plain

 

f:id:cross_hyou:20191009092047j:plain

E2のところがF ratioです。これは、2726/166です。

f:id:cross_hyou:20191009094221j:plain

f:id:cross_hyou:20191009094233j:plain

これでANOVA表が完成しました。このF ratioの16.38が大きいのか小さいのか、qf関数で判定します。


f:id:cross_hyou:20191009094426j:plain

3.865よりも大きなF ratioだとSSAには意味があるということですから最近の平均と先行きの平均には有意な違いがありますね。

aov関数とsummary関数でも同じようにANOVA分析をします。

f:id:cross_hyou:20191009094811j:plain

p値が5.24e-05と0.05以下ですから、最近と先行きでは、有意な差があります。先行きのほうが悲観的なのですね。

今回は以上です。

日銀の短観データの分析1 - R言語でANOVA。6月調査と9月調査に違いがあるとは言えない。

今回は先日発表された日銀の短観データを分析してみます。

f:id:cross_hyou:20191008205806j:plain

日銀のウェブサイトからダウンロードしたファイルは、上のようなものでした。

これを

f:id:cross_hyou:20191008205927j:plain

こんな感じにデータを組みなおして、read.csv関数でR言語で読み込み、分析してみます。

f:id:cross_hyou:20191008210957j:plain

summary関数でそれぞれの変数のおおよその様子をみてみます。

f:id:cross_hyou:20191008211132j:plain

Sectorは製造業、全産業、非製造業の3種類です。Indusはたくさんあります。Timeは調査時期で、Junは2019年6月調査、Sepは2019年9月調査です。Scaleは調査した企業の規模で大企業、中堅企業、中長期業の3種類です。Typeは短観の種類で最近の業況判断と先行きの判断の2種類です。Valueは短観の数値ですね。

Valueが反応変数で数値データ、その他が説明変数で全てカテゴリカルデータです。

この場合の分析は、ANOVA分析です。

とりあえず、短観の数値、Valueの箱ひげ図やヒストグラム、小さい順グラフを描いてみましょう。

f:id:cross_hyou:20191008212243j:plain

f:id:cross_hyou:20191008212255j:plain

hist関数でヒストグラムを描きます。短観は0より大きいと景気がよくて、0より小さいと景気が悪いとなりますから、ヒストグラムを見ると、0よりも大きい値のほうが多いですので、2019年6月と2019年9月はおおむね景気は良かったのではないでしょうか?赤の垂線が平均値、青の垂線が中央値です。

f:id:cross_hyou:20191008212839j:plain

f:id:cross_hyou:20191008212901j:plain

boxplot関数で箱ひげ図を描きます。下に2つ。上に3つの外れ値があります。

f:id:cross_hyou:20191008213240j:plain

f:id:cross_hyou:20191008213300j:plain

sort関数で小さい順に並べて、plot関数でグラフを書きます。

こうして三つのグラフを見ると、短観のデータは都道府県別の人口や県内総生産のデータと比べるとかたよりが少なく、平均値を中心にまとまっているように感じます。

それではANOVAをしてみましょう。ANOVAはAnalysis of Varianceの意味です。

まずは、説明変数をTimeにしましょう。JunとSepで値が短観の値は違うといえるのかどうかを調べてみます。いつものように参考図書は

 

Statistics: An Introduction Using R

Statistics: An Introduction Using R

 

 です。

まずは、グラフでみてみましょう。

f:id:cross_hyou:20191008220442j:plain

f:id:cross_hyou:20191008220453j:plain

黒の点がJun、つまり6月調査の点で赤の点がSep、9月調査の点です。黒い水平線が6月調査の平均値で、赤い水平線が9月調査の平均値です。2つの平均値の違いは、それぞれのデータのばらつきぐあいから見ると大した違いではなさそうですね。

まずは、SSY、sum of the square of the differences between y values and the overall meanを計算します。yというのは今回は短観の値ですね。

f:id:cross_hyou:20191008221142j:plain

SSYは68298.76と求められました。

次は、SSE、the error sum of squaresを計算します。

f:id:cross_hyou:20191008221929j:plain

SSEは67738.55です。
SSYとSSEがわかったので、SSA、the treatment sum of squaresを計算します。

f:id:cross_hyou:20191008222159j:plain

ここで短観のデータの数を確認しておきましょう。

f:id:cross_hyou:20191008223705j:plain

データの数は396個ですね。なので、ANOVA表はこうなります。

f:id:cross_hyou:20191008223841j:plain

SSEのMean squareは67338/394を計算します。そしてF ratioは560をSSEのMean squareで割った値です。

f:id:cross_hyou:20191008224210j:plain

なので、完成したANOVA表は、

f:id:cross_hyou:20191008225400j:plain



となります。

ANOVA表が完成したので、次は、qt関数を使ってF ratioの3.26という値が意味のある値なのかどうかを調べます。

f:id:cross_hyou:20191009095159j:plain

3.865が基準値ですので、それよりも小さい3.26というF ratioでは、Timeによる違いは無い、ということですね。

最後にaov関数とsummary関数を使って結果を確認しましょう。

f:id:cross_hyou:20191008225216j:plain

P値が0.07818と0.05よりも大きいので、JunとSepでValueの平均値に違いが無いという帰無仮説は棄却できません。今回は以上です。

 

 

都道府県別の老人福祉費と児童福祉費の分析6 - R言語で回帰分析の続き。一人当りの老人福祉費と児童福祉を計算した。

 

www.crosshyou.info

 の続きです。

児童福祉費と人口は関係ありましたね。児童福祉費と面積はどうでしょうか?

f:id:cross_hyou:20191005155022j:plain

p-valueは0.1507なので有意ではないですね。面積は児童福祉費には関係ないようです。

では、県内総生産を見てみましょう。

f:id:cross_hyou:20191005155300j:plain

p-valueが2.2e-16よりも小さいので有意なモデルです。I(avgGDP^2)はいらないようですね。削除します。

f:id:cross_hyou:20191005155535j:plain

p値が0.09358と0.05より大きいので、I(avgGDP^2)を削除しても大丈夫です。

モデルを見てみます。

f:id:cross_hyou:20191005155900j:plain

p-valueが2.2e-16より小さいので、有意です。

単回帰分析では人口とGDPが児童福祉費と関係があることがわかりました。

こんどは、coplot関数で児童福祉費が人口の関係がGDPの水準によって変わるかどうか見てみます。

f:id:cross_hyou:20191005160345j:plain

f:id:cross_hyou:20191005160356j:plain

よくわからないですね。対数のデータでやってみます。

f:id:cross_hyou:20191005160720j:plain

f:id:cross_hyou:20191005160730j:plain

GDPのレベルは児童福祉費と人口に大きな影響は及ぼしていないようです。

そもそも、老人福祉費にしても児童福祉にしても人口が多ければ多くなる、人口が少なければ少なくなりますよね。なので、都道府県別の一人当りの老人福祉費、児童福祉費を県内総生産と面積で回帰分析してみましょう。

まずは、一人当りの老人福祉費を計算します。

f:id:cross_hyou:20191005161438j:plain

f:id:cross_hyou:20191005161456j:plain

前に作っておいた自作のgraph関数で小さい順グラフ、ヒストグラム、箱ひげ図を描いています。箱ひげ図を見ると、外れ値がないですね。平均値は一人当り2万3640円です。最低は1万4190円、最大は3万2650円、中央値は2万4520円です。

sort関数で小さい順に表示しましょう。

f:id:cross_hyou:20191005161836j:plain

埼玉県が一番少ないです。老人よりも働き世代の人が多いイメージがあります。高知県が一番多いです。老人比率が高いのかな?

一人当りの児童福祉費も見てみましょう。

f:id:cross_hyou:20191005162147j:plain

f:id:cross_hyou:20191005162158j:plain

箱ひげ図を見ると、一人当りの児童福祉も外れ値はありません。最小値は6千289円、最大値は1万9533円、平均値は1万1408円、中央値は1万1142円です。

それでは、多い県、少ない県を見てみましょう。埼玉県は多いかな?高知県は少ないかな?

f:id:cross_hyou:20191005162513j:plain

一番少ないのは神奈川県、多いのは沖縄県です。埼玉県は7千854円と平均以下ですね。そして高知県は1万4169円と平均以上です。よく見ると神奈川県は老人福祉費では2番目に少なく、児童福祉費は一番少ないです。老人と児童には優しくない県なのか。。。他になにか恩恵があるのか?

一人当りの老人福祉費と児童福祉の散布図を描いてみます。

f:id:cross_hyou:20191005163137j:plain

f:id:cross_hyou:20191005163150j:plain

lm関数とablineで回帰直線も重ねてみました、右肩上がりの関係ですね。

今回は以上です。

 

都道府県別の老人福祉費と児童福祉費の分析5 - R言語で回帰分析をする。

 

www.crosshyou.info

 の続きです。

今回はR言語で回帰分析をします。老人福祉費または児童福祉を人口、面積、県内総生産の3つの変数で回帰分析します。事前の予想は、人口と県内総生産は関係あるけど、面積は関係ないんでは?と思います。

早速やってみます。lm関数を使います。

f:id:cross_hyou:20191005125818j:plain

まずは、各変数の交差項と2乗した項目を加えたmaximum modelからスタートします。

p-value は2.2e-16より小さいので有意なモデルです。Adjusted R-squaredは0.9858ですのであてはまりもいいですね。avgPop:avgArea:avgGDPの3変数の交差項はp値が0.05よりも大きいので削除します。

update関数を使います。

f:id:cross_hyou:20191005130251j:plain

anova関数でmodel1とmodel2を比較すると、p値が0.07087と0.05よりも大きいので、model1とmodel2では有意な違いはありません。従って、より単純なmodel2を採用します。model2をsummary関数で見てみます。

f:id:cross_hyou:20191005130532j:plain

I(avgArea^2)はp値が0.334315なので削除しましょう。

f:id:cross_hyou:20191005130857j:plain

p値が0.3343なので、model2とmodel3では有意な違いは無いとわかります。

model3を見てみます。

f:id:cross_hyou:20191005131120j:plain

avgAreaやavgGDPはp値が0.05ですが、avgAre:avgGDPなどの交差項のp値が0.05なので外せないです。残差プロットを見てみましょう。

f:id:cross_hyou:20191005131615j:plain

今度は、対数をとったデータ回帰分析してみましょう。

f:id:cross_hyou:20191005131945j:plain

さきほどと同じような手順で不要な変数を削除していきます。まずは、logPop:logArea:logGDPですね。

f:id:cross_hyou:20191005132207j:plain

anova関数でlmodel1とlmodel2を比較します。p値が0.5586ですから、二つのモデルに有意な違いはありません。なので、より単純なlmodel2をさらに詳しく調べます。

f:id:cross_hyou:20191005132413j:plain

I(logArea^2)を削除しましょう。

f:id:cross_hyou:20191005132618j:plain

p値は0.5492と0.05より大きいですから、lmodel2とlmodel3では有意な違いは無いです。lmodel3を見てみます。

f:id:cross_hyou:20191005132825j:plain

除外できる変数はなくなりましたね。残差プロットを見てみましょう。

f:id:cross_hyou:20191005133103j:plain

f:id:cross_hyou:20191005133140j:plain

老人福祉費は重回帰分析でしたので、児童福祉費は単回帰分析をしてみます。

まずは、人口で回帰してみましょう。

f:id:cross_hyou:20191005133710j:plain

avgPopのp値が0.202なので削除します。

f:id:cross_hyou:20191005133929j:plain

p値が0.202なので、avgPopは無くても大丈夫ですね。ChildPop2を見てみます。

f:id:cross_hyou:20191005134124j:plain

散布図の上に回帰線を描いてみましょう。

まずは、回帰線のxvを作成します。

f:id:cross_hyou:20191005135949j:plain

f:id:cross_hyou:20191005140005j:plain

今回は以上です。

 

都道府県別の老人福祉費と児童福祉費の分析4 -R言語で Variance(分散)を調べ、standard error(標準誤差)とconfidence interval(信頼区間)を求める。

 

www.crosshyou.info

 の続きです。

今回は、variance(分散)を調べます。上に画像のある、Michael J. CrawleyのStatistics An Introduction using Rを参考にします。

varianceはvar関数で簡単にわかります。

10年間平均の老人福祉の対数の分散は、

f:id:cross_hyou:20191002175641j:plain

0.07579028とわかります。

varianceは、

それぞれの観測値と平均値の差を2乗する

それを合計する

観測値の数-1で割る。

です。

やってみます。

f:id:cross_hyou:20191002180332j:plain

はい、0.07579028ですね。

10年間平均の児童福祉の平均値もvar関数で計算しましょう。

f:id:cross_hyou:20191002180523j:plain

0.07496646です。

対数でなくて実際の金額でも計算しましょう。

f:id:cross_hyou:20191002180726j:plain

2.282334e+15は小数点以下の数字が15個あるということです。

228京2334兆0000億0000万0千円^2です。

8.802981e+14は小数点以下が14個あるといことですので、

88京0298兆1000億0000万0千円^2です。

老人福祉費と児童福祉費では分散が228京円^2と88京円^2と大きく違います。二つの分散は違うと言えるでしょうか?var.test関数を使います。

f:id:cross_hyou:20191002181929j:plain

p-valueが0.001601と0.05よりも小さいので、ふたつのvarianceは違うと言えます。

varianceは、standard errors、confidence intervalsを計算するのに使います。

standard errorsは標準誤差と日本語では言うようです。

standard errorsは、

variance(分散)をデータの個数で割る

その平方根を計算する

で計算できます。

やってみます。

f:id:cross_hyou:20191002182726j:plain

standard errorは平均値と組み合わせて使うようです。

f:id:cross_hyou:20191002182928j:plain

10年間の老人福祉費の47都道府県の平均値は、

569億4643万7千円±69億6852万2千円(1 s.e., n = 47)

10年間の児童福祉費の47都道府県の平均値は、

276億9779万3千円±43億2778万8千円(1 s.e., n = 47)

と記述するようです。

confidence intervalは信頼区間と言うようです。

これは計算方法をここに書くのが難しいですが、t-valueにstandard errorを掛けた値だそうです。

t-valueはqt関数でできます。t-valueは95%だと観測数でかわりますが、だいたい2ぐらいです。

f:id:cross_hyou:20191002184130j:plain

こうして信頼区間がわかったらm

10年間の老人福祉費の47都道府県の平均値は、

569億4643万7千円±140億2690万7千円(95% CI, n = 47)

10年間の児童福祉費の47都道府県の平均値は

276億9779万3千円±87億1138万5千円(95% CI, n = 47)

と記述するそうです。

今回は以上です。

 

都道府県別の老人福祉費と児童福祉費の分析3 - Central Tendency(代表値)を調べる。算術平均値、中央値、幾何平均値、調和平均値。

 

www.crosshyou.info

 の続きです。

今回は、老人福祉費と児童福祉費のCentral Dendencyを調べようと思います。

いつものように Michael J. Crawleyの Statistics An introduction using R を参考にします。

 

Statistics: An Introduction Using R

Statistics: An Introduction Using R

 

 Central Tendencyって日本語だと代表値というのですね。

 まずは、arithmetic mean、算術平均値です。

f:id:cross_hyou:20190929090419j:plain

算術平均は学校で習う普通の平均値ですね。データの合計値をデータの個数で割ります。R言語ではmean関数です。老人福祉費の平均値は569億4643万7千円です。都道府県平均で570億円ぐらいです。

児童福祉費の算術平均を調べます。

f:id:cross_hyou:20190929090756j:plain

児童福祉費の算術平均は、276億9779万3千円です。老人福祉費の半分以下ですね。

次は、median, 中央値を調べます。データを小さい順に並べて、その真ん中の値ですね。都道府県は47都道府県ですから24番目の値です。

老人福祉費の中央値を調べます。

f:id:cross_hyou:20190929091134j:plain

老人福祉費の中央値は長崎県の399億5091万8千円です。median関数で中央値がわかります。

児童福祉費の中央値はどの県でしょうか?

f:id:cross_hyou:20190929091439j:plain

岐阜県が児童福祉費の中央値ですね。184億8356万0千円です。中央値でも児童福祉費は老人福祉費の半分以下ですね。

つぎは、geometric meanです。日本語だと幾何平均値というそうです。

f:id:cross_hyou:20190929092216j:plain

幾何平均は全てのデータを掛け合わせて、1/データの個数 の平方根をとります。

対数の平均値をexp関数で対数からもとの数に戻しても同じです。老人福祉費の幾何平均値は452億7913万7千円です。

児童福祉費の幾何平均値も求めます。

f:id:cross_hyou:20190929092636j:plain

児童福祉費の幾何平均値は、214億6529万0千円でした。

最後は、harmonic meanです。調和平均値というそうです。逆数の算術平均を求めて、その逆数だそうです。東京から大阪まで行きは時速100キロ、帰りは時速200キロのときの平均時速は?というときの平均です。

老人福祉費の調和平均を求めましょう。

f:id:cross_hyou:20190929093152j:plain

老人福祉費の調和平均値は、384億3977万2千円です。

児童福祉費の調和平均値を求めます。

f:id:cross_hyou:20190929093423j:plain

184億0644万3千円です。

今回は以上です。

 

都道府県別の老人福祉費と児童福祉の分析2 - グラフにしてデータの分布形状を確認する。

 

www.crosshyou.info

 の続きです。前回は都道府県別の10年間の平均値をtapply関数を使って算出しました。

せっかくなので、グラフにしてどのようなデータ分布なのか見てみます。

まずは、いっぺんに三つのグラフを描く関数を作ります。

f:id:cross_hyou:20190928122352j:plain

小さい順に並べたグラフ、ヒストグラム、箱ひげ図をいっぺんに描きます。

赤い線で平均値、青い線で中央値を表示します。それでは老人福祉費を見てみます。

f:id:cross_hyou:20190928122603j:plain

金額の大きい複数の都道府県があることがわかります。

次は児童福祉です。

f:id:cross_hyou:20190928122909j:plain

f:id:cross_hyou:20190928122920j:plain

社会福祉費と同じような形状ですね。

総人口を見てみます。

f:id:cross_hyou:20190928123130j:plain

f:id:cross_hyou:20190928123144j:plain

総人口も同じようにヒストグラムの裾野が右に広がっています。

総面積を見てみます。

f:id:cross_hyou:20190928123430j:plain

f:id:cross_hyou:20190928123443j:plain

面積は、北海道が桁違いに大きいことがわかります。

県内総生産はどうでしょうか?

f:id:cross_hyou:20190928123738j:plain

f:id:cross_hyou:20190928123752j:plain

県内総生産は東京都がダントツで、2番手グループに大阪、愛知、神奈川があります。

5つの変数がすべて分布の裾野が右側に広がっていますので、対数変換して分布の形状を左右対称に近づけましょう。

f:id:cross_hyou:20190928124426j:plain

f:id:cross_hyou:20190928124438j:plain

log関数で10を底数にして変換しました。すこし分布形状が左右対称に近づきました。

児童福祉費も定数変換します。

f:id:cross_hyou:20190928124734j:plain

f:id:cross_hyou:20190928124755j:plain

こちらも分布の峰が中央よりに移動しました。

総人口も対数変換します。

f:id:cross_hyou:20190928125051j:plain

f:id:cross_hyou:20190928125105j:plain

総面積を対数変換します。

f:id:cross_hyou:20190928125335j:plain

f:id:cross_hyou:20190928125346j:plain

対数をとっても北海道の面積はけた違いですね。

県内総生産の対数をとります。

f:id:cross_hyou:20190928125620j:plain

f:id:cross_hyou:20190928125633j:plain

北海道の面積ほどではないですが、東京の県内総生産も突出していますね。

今回はここまでです。

ここで、いままで作成してきたオブジェクトをls関数で確認しておきましょう。

f:id:cross_hyou:20190928125935j:plain

avgAreaは10年間平均の都道府県別の総面積

avgChildは10年間平均の都道府県別の児童福祉費

avgGDPは10年間平均の都道府県別の県内総生産

avgOldは10年年間平均の都道府県別の老人福祉費

avgPopは10年間平均の都道府県別の総人口

df1は生データ

df2は生データからNA行を削除したデータフレーム

graphは3つのグラフをいっぺんに描く関数

logAreaは10年間平均総面積を対数変換したもの

logChildは10年年間平均児童福祉費を対数変換したもの

logGDPは10年間平均県内総生産を対数変換したもの

logOldは10年間平均老人福祉費を対数変換したもの

logPopは10年間平均総人口を対数変換したもの

です。