www.crosshyou.info

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

毎月勤労統計調査の分析V2_2 - もっと多くのデータで基本統計量を算出

前回のデータ量が少なかったので、今回はもっと大きなデータで毎月勤労統計調査のデータを分析したいと思います。

こんな感じです。

f:id:cross_hyou:20180922174832j:plain

これがhead関数で始めの6行を表示したもの。

そして下図がsummary関数で要約統計量を表示したものです。

f:id:cross_hyou:20180922174949j:plain

行番号もMaxが4640です。つまりデータレコードが4640あるということです。

給与総額などでNAが360もあります。これは例えば、持ち帰り配達飲食サービスなどの業種で1000人以上の労働者がいる事業所が無い場合もあるからです。

前回は全産業トータルのデータだけでしたが、今回は細かい業種まで取込みましたので、例えば、給与総額の最小値は0円、最大値は163万0293円となっています。

給与総額0円は出勤日数0だとそうなっちゃうのかもしれないですね。

次回以降は、このデータでいろいろR言語の練習をしたいと思います。

 

毎月勤労統計調査の分析V2_1 - 基本統計量(summary関数、apply関数とsd関数, mean関数)

今回は、毎月勤労統計調査の分析をしたいと思います。

以前にもこの統計は分析したことがあるので、2回目です。

e-Stat(政府統計の総合窓口)のサイトを訪問したら、毎月勤労統計が新着でありました。

f:id:cross_hyou:20180922125428j:plain

クリックしてみました。

f:id:cross_hyou:20180922125441j:plain

「毎月勤労統計調査全国調査は、日本標準産業分類に基づく16大産業に属する常用労働者5人以上の事業所を対象に、賃金、労働時間及び雇用の変動を毎月把握する調査です。調査対象事業所は、常用労働者5人以上の約190万事業所(経済センサス-基礎調査)から抽出した約33,000事業所で、名目賃金(現金給与総額)や実質賃金、所定内及び所定外労働時間などがわかります。調査の結果は、景気動向を判断するための指標の一つとなっているほか、厚生労働政策や経済政策の基礎資料、企業の労働条件決定の際の参考資料として幅広く活用されています。」とのことです。

データのExcelファイルはこんなものでした。

f:id:cross_hyou:20180922130054j:plain

これをデータ分析しやすいように加工したCSVファイルがこちらです。

f:id:cross_hyou:20180922133208j:plain

規模、性別、形態を記号ではなくてわかりやすくしました。2018年7月のデータです。

このCSVファイルをread.csv関数でR言語に読込んで分析しましょう。

f:id:cross_hyou:20180922135001j:plain

head関数で始めの6行のデータを表示しています。summary関数で要約統計量を表示しています。

規模は、100-499人が5、1000人以上が5とそれぞれの規模が5個ずつデータを持っていいるようですね。table関数でどのような規模があるか見てみます。

f:id:cross_hyou:20180922135032j:plain

5人以上、30人以上、500人以上、1000人以上という「以上」で区別する分け方と、5-29人、30-99人、100-499人、500-999人というレンジで区別する分け方が混在しています。

性別は、男女計、男、女、NAと4種類です。NAというのは形態が一般とパートの分です。なので、性別と形態の組み合わせは、男女計&合計、男&合計、女&合計、NA&一般、NA&パートの5種類です。これと規模が8種類あるので、5 x 8 = 40のデータ行数ということですね。

現金給与総額を見ると、最小値は9万2196円、最大値は82万3642円です。

総労働時間は、最小値は80.3時間、最大値は173.6時間です。

summary関数は標準偏差を表示しないので、apply関数とsd関数で標準偏差を計算しましょう。

f:id:cross_hyou:20180922140825j:plain

apply関数とsd関数で計算した標準偏差をdatasdという箱に入れて、round関数で小数点以下2桁表示しました。na.rm=TRUEというのを付け加えてNAがあっても標準偏差を計算するようにしています。na.rm=TRUEを付けないとエラーになってしまいます。

標準偏差を平均値で割った変動係数も計算しましょう。

f:id:cross_hyou:20180922141256j:plain

apply関数とmean関数で平均値を計算してから、標準偏差を計算した平均値で割って変動係数を計算しています。

データのバラツキが一番小さいのは、出勤日数ですね。バラツキが一番大きいのは6月と7月の労働者数ですね。

 

貴金属流通統計調査の分析6 - 金の需要を線形回帰してみる。(lm関数)

今回は貴金属流通統計調査の金の需要データを使って、線形回帰の練習をしてみたいと思います。まずは、CSVファイルに保存しているデータをread.csv関数を使って読込みます。

f:id:cross_hyou:20180921140201j:plain

このようなデータです。電機機械の需要が一番多いですね。グラフにしてみましょう。plot関数を使います。

f:id:cross_hyou:20180921140441j:plain

f:id:cross_hyou:20180921140454j:plain

type="b"とするとこのように、グラフがプロット点とラインの両方を表示します。

まずは、時間の経過とともに値が増えている感じもしますので、時間ファクターで線形回帰してみます。lm関数です。

f:id:cross_hyou:20180921140916j:plain

まず、time <- c(1:17)で時間ファクターを作成して、lm関数を使っています。結果をmodel1という名前の箱に入れて、summary関数で結果を呼び出しています。

一番下の行のp-valueを見ると、値が0.05193なので有意ではないです。

このモデルに他の需要を加えてみましょう。どれを入れましょうか?相関係数をcor関数で出して、一番相関しているものを追加してみましょう。

f:id:cross_hyou:20180921141259j:plain

宝飾との相関係数が0.425で一番高いですね。宝飾を加えたモデルをみてみましょう。

f:id:cross_hyou:20180921141534j:plain

p-value = 0.1405 > 0.05 で有意ではないですね。宝飾を加えたモデルのほうがp-valueが大きくなってしまいました。

さらに、美術工芸品を加えたモデルはどうでしょうか?

f:id:cross_hyou:20180921141907j:plain

p-value = 0.08467 > 0.05 ですので有意ではないです。でも、model2よりはいいですね。

歯科医療を加えてみましょう。

f:id:cross_hyou:20180921142225j:plain

p-value = 0.1679 > 0.05 なので有意ではありません。

逆相関ですが、メッキよりも相関係数の絶対値は大きいメダルを加えてみます。

f:id:cross_hyou:20180921142645j:plain

p-value = 0.2197 > 0.05 です。有意ではないです。

最後はメッキも加えてみましょう。

f:id:cross_hyou:20180921143034j:plain

p-value = 0.3353 > 0.05 で有意ではないです。

ということで、金の電機機械の需要を他の用途の需要で線形回帰するのは難しいということがわかりました。

lm関数で簡単に線形回帰分析ができるのが凄いですね。

 

貴金属流通統計調査の分析5 - 金の各需要の相関関係を見る。(plot関数、cor関数、cor.test関数)

今回は、金の各需要の相関関係を見たいと思います。

まずは、CSVファイルに保存してあるデータをread.csv関数でR言語に読込みます。

f:id:cross_hyou:20180920112154j:plain

相関係数を出す関数はcor関数です。早速やってみましょう。

f:id:cross_hyou:20180920112637j:plain

gold_matrix <- gold[ , c(-1,-2)] で計算に必要のない年と月の列を削除しています。

gold_matrix <- as.matrix(gold_matrix) でgold_matrixのデータ型をデータ・フレームから行列(マトリックス)に変換しています。cor関数がデータ・フレームでは動かないからです。

そして、cor関数で相関を計算しています。ちょっと桁数が多くて見づらいので、round関数で値をまるめましょう。

f:id:cross_hyou:20180920113220j:plain

わかりやすくなりました。

相関係数の大きいのは、歯科医療とメッキの組み合わせが0.6、電機機械と宝飾の組み合わせが0.43ですね。

この二つの組み合わせの相関が有意かどうかをcor.test関数で確認しましょう。

f:id:cross_hyou:20180920113628j:plain

p-value = 0.01165 < 0.05 ですから、歯科医療とメッキには相関関係はありますね。正の相関です。

f:id:cross_hyou:20180920113840j:plain

p-value = 0.08888 > 0.05 ですから、電機機械と宝飾に相関関係があるとはいえません。

plot関数で2つ組み合わせの散布図を見てみましょう。

f:id:cross_hyou:20180920114314j:plain

par(mfrow=c(1,2)) で始めにグラフウィンドウを 1 x 2 に設定します。

そのあとにplot関数で散布図を書きます。

f:id:cross_hyou:20180920114452j:plain

このようにみると確かに歯科医療とメッキの組み合わせのほうが相関していることがわかりますね。

ちなみにplot関数をデータ・フレームに適用するとこうなります。

f:id:cross_hyou:20180920114926j:plain

f:id:cross_hyou:20180920114934j:plain

散布図のマトリックスになります。

 

 

貴金属流通統計調査の分析4 - 金需要の平成29年と平成30年で違いがあるか?(t.test関数、wilcox.test関数、binom.test関数)

今回は、金の需要が平成29年と平成30年で増加、もしくは減少しているのかを調べたいと思います。

まずは、CSVファイルに保存してあるデータをread.csv関数でR言語に読込みます。

f:id:cross_hyou:20180919203803j:plain

平成29年は1月から10月まで、平成30年は1月から7月までのデータがありました。季節によって金需要に差があるかもしれませんので、平成29年も平成30年も1月から7月までのデータで比較しましょう。

とりあえず、一番数量の多い電機機械のデータで調べてみましょう。

まずは、それぞれのベクトルを作成します。

f:id:cross_hyou:20180919204237j:plain

f:id:cross_hyou:20180919204248j:plain

こうして2つのベクトルが作成できました。

まずは2つの合計値を確認して、どっちが多いか少ないかを見ます。

f:id:cross_hyou:20180919204455j:plain

平成29年よりも平成30年のほうが多いですね。

まずは、平均値で検定します。ほんとうは2つのデータとも7個しかデータないので平均値は検定できないんじゃないかと思いますが、R言語の練習なのでかまわずやります。t検定のt.test関数です。

f:id:cross_hyou:20180919204911j:plain

p-value = 0.02067 < 0.05 なので有意です。つまり平成29年と平成30年では違いがある、ということですね。

次は、2つのデータの中央値(分布の位置)に違いがあるかを検定します。ウィルコクソン=マン・ホイットニー検定のwilcox.test関数です。

f:id:cross_hyou:20180919205349j:plain

p-value = 0.03125 < 0.05 なので有意です。つまり、平成29年と平成30年ではデータの分布に違いがあるということです。

最後に符号検定をしてみましょう。平成29年と平成30年を比較して、増加しているのか減少しているのかを調べます。

f:id:cross_hyou:20180919205619j:plain

平成30年のほうが多い月が6か月あって、平成29年のほうが多い月が1月ありました。

これで二項検定のbinom.test関数を使います。

f:id:cross_hyou:20180919205904j:plain

p-value = 0.125 > 0.05 なので有意とは言えないです。つまり平成30年のほうが平成29年よりも多いとは言い切れない、ということです。

平均値と標本の位置では違いがある、二項検定では違いが無いという結果になりました。

貴金属流通統計調査の分析3 - 金の需要の割合を円グラフで表示する。(pie関数)

今回は、貴金属流通統計調査のデータを使って、金の需要の割合を円グラフにしたいと思います。

まず、CSVファイルに保存してあるデータをread.csv関数でR言語に読込みます。

f:id:cross_hyou:20180919190625j:plain

このようなデータです。金の需要としては、電機・機械で使われたり、歯科医療で使われたり(いわゆる金歯ですね。)、メッキ、宝飾、美術工芸品、メダルなどがあります。

この数値はグラムを表しています。

まず、H29年1月からH30年7月までのデータを合計しましょう。apply関数とsum関数を使います。

f:id:cross_hyou:20180919191124j:plain

これらの合計値を計算しましょう。sum関数です。

f:id:cross_hyou:20180919191259j:plain

比率を計算します。

f:id:cross_hyou:20180919191637j:plain

これで円グラフが書けます。

f:id:cross_hyou:20180919192225j:plain

desという変数を作って、そこに「電気機械 61.78%」というような説明(description)を流し込みます。paste関数で文字と文字をつなぎ合わせます。sep=""でつなげる文字の間にスペースを入れないようにしています。

そして、pie関数で円グラフを作成します。clockwise=TRUEで12時方向からデータがスタートするように設定しています。

そしてできた円グラフがこちらです。

f:id:cross_hyou:20180919192523j:plain

というように、なりました。私が思っていた金の需要って宝飾とか美術品で使われるものでしたが、実際は電機・機械が6割もあるのですね。

 

 

 

貴金属流通統計調査の分析2 - 各需要の比率を計算する。

今回は、金の各需要の比率を計算しようと思います。

まずは、CSVファイルのデータをread.csv関数で読込みます。

f:id:cross_hyou:20180915153557j:plain

どのようなデータか全体表示します。

f:id:cross_hyou:20180915153706j:plain

元のデータファイルには、「その他」の項目もありましたが、このブログの分析では除外しています。なので、各需要の比率といっても実際の比率ではなく、あくまでもこのブログ内での比率になります。

年と月の列は計算にいらないので、年と月の無いデータフレームを作成しましょう。

f:id:cross_hyou:20180915154012j:plain

gold[ , c(-1,-2)] としていますが、gold[ , c(3:8)]としても同じですね。

各行の合計値を計算します。apply関数とsum関数ですね。

f:id:cross_hyou:20180915154333j:plain

この合計値の推移をグラフにしてみましょうか。plot関数を使います。

f:id:cross_hyou:20180915154625j:plain

f:id:cross_hyou:20180915154635j:plain

 

全体としては右肩上がりですね。金の需要は増えているということですね。

本題に戻りまして、比率を計算しましょう。

f:id:cross_hyou:20180915154948j:plain

こうして比率を見ると、電機機械が6割程度で、歯科医療と宝飾が2割弱、メッキが5パーセント前後、美術工芸品が1~2パーセントぐらい、メダルはほとんど無視していいぐらい、ということがわかります。

この比率の推移をグラフにしましょう。for構文とplot関数を使っていっぺんに作成します。

f:id:cross_hyou:20180915160016j:plain

for (i in c(1:6)) が i を1から6まで繰り返しますよ、という意味です。
{ と } で囲んだ部分、plot(goldP[ , i], main=colnames(goldP)[i], type="b" を繰り返す、っていうことですね。type="b" でグラフのタイプをプロットとラインの両方にしています。そうしてできたのが下の図です。

f:id:cross_hyou:20180915160356j:plain

宝飾の比率が高まっている感じですね。歯科医療や美術工芸品、メッキは低下気味、電機機械は横ばいという状態です。

今回は各需要の比率を計算しました。