生産者の米穀在庫等調査の分析1 - 都道府県別の基本統計量 北海道は別格だ。

政府統計の総合窓口(e-Stat)に新着データとして、「生産者の米穀在庫等調査」というデータがありました。

クリックしてみます。

本調査は、毎月、農家の米穀の在庫量等を調査し、農家1戸当りのうるち米及びもち米の供給量、消費量、販売量、在庫量等を全国、都道府県別に提供しています。

とのことです。

どういう種類のデータがあるでしょうか？

農家1戸当りのデータを都道府県別にまとめたものですね。ファイルを見てみます。

こういうデータですね。

これをR言語で読込みやすいようにCSVファイルに加工しました。

このCSVファイルをread.csv関数でR言語に読込み、summary関数で基本統計量を表示しましょう。

標準偏差を計算しましょう。sd関数です。

apply関数とsd関数で都道府県を除くすべての列にsd関数を適用しています。そして結果をround関数で整数表示にしています。

平均値も同じように計算します。mean関数です。

変動係数(CV)を計算します。標準偏差 / 平均値 です。

飯用の変動係数が一番小さい(0.19)ですね。まあ、農家1戸当りの飯用ですから、各県でそれほど大きくは変わらないですよね。それに比べると、販売量はCVが1.21ですから都道府県によって大きく変わる、ということですね。

販売量の箱ひげ図と飯用の箱ひげ図を比較してみましょう。boxplot関数です。

販売量のほうは、大きな外れ値が一つありますね。

各項目同士の相関係数を計算しましょう。cor関数です。

脱穀量と供給量の相関係数が1.000です。供給量 = 脱穀量 + 購入量 という関係ですが、脱穀量の平均値は6080、購入量の平均値は63と購入量は脱穀量の1％ほどですから供給量と脱穀量の相関が1.000になるのも納得です。

供給量と販売量の相関係数も1.000ですね。

それぞれのデータの分布の様子を見るために、データの小さい順に並び替えて棒グラフを描いてみましょう。棒グラフはbarplot関数、並び替えはsort関数、for関数で一度に10個のグラフを作ります。

はじめに、par(mfrow=c(2,5))で2行5列のグラフの箱を作っておいて、そのあとにfor関数で10個の棒グラフを作ります。

こうしてグラフを見ると、1つの都道府県が突出しているデータ項目が多いことがわかります。そしてその突出した都道府県とは。。。北海道ですね。。。確認してみましょう。order関数で並び替えて表示しましょう。

北海道は、飯用と無償譲渡量以外はダントツの1位だとわかります。

長崎は自分で食べる分(飯用)が多いんですね。

北海道がどれだけ突出しているかを見るために、各データを標準化して、標準化したデータの平均値を各都道府県ごとに計算してみましょう。

標準化した数値はZ値といいます。scale関数です。

scale(data[ ,-1]) で都道府県の列以外の数値データを標準化して、その結果をas.data.frame関数でデータフレームにしています。

apply関数とmean関数で行ごとのの平均値を算出してそれをAverageという列名でz_dataのデータフレームに追加しています。

cbind関数でdataの1列目、つまり都道府県とz_dataと結合しています。

colnames関数でz_dataの1列目の列名を"都道府県"と名付けています。

そして、order関数でAverageの大きい順に並び替えています。

こうしてみると、北海道がダントツだとわかりますね。

棒グラフにしてみましょう。

一番右の飛びぬけている棒が北海道ですね。北海道は別格ですね。

毎月勤労統計調査の分析V2_8 - 特掲産業と全データの比較(時給と残業時間比率)

今回は、毎月勤労統計調査の2018年7月のデータを使って、特掲産業だけに絞ったデータフレームを作成し、分析ごっこをしてみたいと思います。

「毎月勤労統計調査における記号の見方」には「特掲産業」に分類されている産業があります。

これらの産業だけに絞り込んだデータフレームを作成して分析をしようと思います。

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

業種コードの列で絞り込めばいいですね。

業種コードのところを見てみると、E091が40、E092が40などと絞り込まれていることがわかります。これまでの分析と同じように時給を計算してみましょう。

全データの時給と比較してみましょう。

平均値は特掲産業は2592円で全データは2661円

中央値は特掲産業は2321円で全データは2461円です。

特掲産業のほうが時給の水準は低い様子ですね。

ヒストグラムにして分布を比較してみましょう。hist関数です。

どちらも同じような分布形状ですね。

平均値に差があるか、t検定をしましょう。t.test関数を使います。

p-value = 0.1043 > 0.05 ですから、特掲産業の時給の平均値と全データの時給の平均値に有意な違いはないです。

時給に分布形状に違いがあるか、ウィルコクソン＝マン・ホイットニー検定をしてみましょう。wilcox.test関数を使います。

p-value = 0.008095 < 0.05 ですから、特掲産業の時給と全データの時給では分布形状に違いがあると言えます。特掲産業のほうが時給が低い方にシフトしていますね。

残業時間比率も計算して比較してみましょう。

所定外労働時間 / 総労働時間 x 100 で残業時間比率(%)が計算できます。

特掲産業の残業時間比率の平均値は、6.922%、中央値は6.670%

全データの残業時間比率の平均値は、7.201%、中央値は7.130%です。

全データのほうが残業時間比率は高いようですね。

ヒストグラムで分布を見てみましょう。hist関数です。

特掲産業のほうが残業時間比率が小さい区分が多いような印象です。

平均値に差があるか、t検定でみてみましょう。t.test関数です。

p-value = 0.03341 < 0.05 ですから、特掲産業の残業時間比率と全データの残業時間比率では、その2つの平均値には有意な差があります。特掲産業のほうがあんまり残業しないということですね。

分布の位置に違いがあるか検定しましょう。ウィルコクソン＝マン・ホイットニー検定をします。wilcox.test関数です。

p-value = 0.01024 < 0.05 ですから、特掲産業と全データでは、残業時間比率の分布位置に違いがあります。特掲産業のほうが残業時間比率が小さいほうに分布が多いということですね。

毎月勤労統計調査の分析V2_7 - 規模別、性別、形態別の時給を計算してみる。

今回は、規模別、性別、形態別の時給を計算してみましょう。

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

もう、面倒なので、このまま総労働時間と給与総額を使って時給を計算してしまいましょう。

最低時給が816円です。。。これって最低賃金に抵触しないのかな。。最高時給は8946円です。平均は2661円、中央値は2468円です。

ヒストグラムを見てみましょう。hist関数です。

箱ひげ図も書いてみましょう。boxplot関数です。

では、規模別の時給の平均値を計算してみましょう。tapply関数とmean関数です。

round関数で小数点以下を切り上げて、sort関数で小さい順に表示しました。

時給の一番低いのは規模が一番小さい5-29人で、時給の一番高いのは規模が一番大きい1000人以上です。棒グラフにしてみます。barplot関数です。

性別の時給を見てみましょう。tapply関数とmean関数です。

男の時給が一番高く、女の時給が一番低く、男女合計がその間というわかりやすい結果ですね。棒グラフにしましょう。barplot関数です。

勤務形態別の時給を計算します。tapply関数とmean関数です。

パートが低く、一般が高く、パートと一般の合計がその間というこれまたわかりやすい結果です。

棒グラフにしましょう。barplot関数です。

時給の高いのは規模は1000人以上の事業所、性別は男、勤務形態では一般です。

そこで、規模は1000人以上で性別は男、または規模は1000人以上で勤務形態は一般というデータフレームを作成して、その平均時給を計算してみましょう。

平均時給は4353円です。

その反対に規模は5-29人で、性別は女、または規模は5-29人で勤務形態はパートというデータフレームを作成して平均時給を計算してみましょう。

平均時給は1523円で高時給のグループの半分以下です。

この二つを棒グラフにして比較してみましょう。

ヒストグラムでも比較してみましょう。

par(mfrow=c(2,1))と始めにグラフエリアの配置を、2行1列と指定してから、hist関数を使うとこのように1つのグラフ画面の中に二つのヒストグラムを配置できます。

毎月勤労統計調査の分析V2_6 - 産業の大分類を性別や規模、勤務形態でもう少し細かく診てみる。

前回の分析では、産業の大分類ベースでは、時給とパート比率は逆相関の関係があることがわかりました。今回は産業の大分類ベースで企業規模別にしてみようと思います。

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

データの種類としては、行番号、業種、業種コード、規模、性別、形態、六月末人数、増加人数、減少人数、七月末人数、パート人数、出勤日数、総労働時間、所定内労働時間、所定外労働時間、給与総額、定期給与、所定内給与、超過労働給与、特別支払給与、の20種類です。

業種が調査産業計も大分類も中分類も特掲産業もごちゃ混ぜになっていました。なので、まずは中分類だけのデータセットを作りましょう。

このようなコマンドです。data$業種コード=="C" | data$業種コード=="D" と条件と条件を | でつなげるとOR検索になります。

規模をちょっと見てみましょう。table関数を使います。

規模をよく見ると、5人以上、30人以上、500人以上、1000人以上という「以上」で区分されているものと、5-29人、30-99人、100-499人、500-999人とレンジで区分されているものが混ざっています。なので、希望は5-29人、30-99人、100-499人、500-999人、1000人以上と人数が混ざらないものだけにしましょう。

性別を見てみましょう。table関数を使います。

男女合計を除きましょう。

形態の種類を見てみましょう。table関数です。

あ、パートと一般が無くなって合計だけになりました。パートと一般は区別して分析したいので、男女合計を除いたのが余計だったのですね。なので、男女合計を除くコマンドの前までのコマンドをもう一度走らせます。

これでいいでしょう。

それでは、七月末の人数が一番多いのはどのような業種、規模、性別、形態なのかorder関数で並び替えてみましょう。

卸売小売で規模は5-29人、性別は男女合計、勤務形態も一般とパートの合計が515万人で一番多いですね。

給与総額が多いのはどこでしょうか？

おお～～！建設業、1000人以上、男、合計が150万6745円です！七月はボーナスが出たからだと思いますが、すごいですね！

そうだ、人数の少ないセグメント、給与総額の少ないセグメントも確認しましょう。

order関数で小さい順に並び替えるのは面倒なので、tail関数を使いましょう。

あら、NAになってますね。。。na.omit関数でNAの列は削除してしまいましょう。

鉱業、採石、砂利採取が少ないですね。。。

給与総額はどうでしょうか？

宿泊飲食サービスなどのサービス業の給与総額が少ないのですね。。

毎月勤労統計調査の分析V2_5 - 時給換算で一番時給の良い業界は？

今回は前回作成した毎月勤労統計調査2018年7月データの大分類データで時給換算で一番時給の良い業界はどこなのか調べたいと思います。

まずは、前回の結果、産業分類大分類別の給与総額ランキングです。

学術研究、専門・技術サービスが一番給与総額が高かったですね。

時給換算したらどうなるでしょうか？

総労働時間も一緒に表示してみましょう。

データフレーム[ , c(”列名１", "列名２")]

のように表示したい列名を選択して表示することができます。

一番給与総額の低かった宿泊飲食サービスは、総労働時間も99.2時間で一番少ないです。

それでは、時給を計算してみましょう。

こうなりました。

では、時給の高い順に並び替えて表示しましょう。

こうなりました。時給でも学術研究、専門・技術サービスが一番高いです。時給3476円です。そして、宿泊飲食サービスが時給換算でも一番低く、1400円です。

宿泊飲食サービスは、パートの比率が高いからだと思うのですが、どうでしょうか？

パート比率を計算してみましょう。あ、その前にデータの全項目を見てみましょう。

summary関数です。

"パート人数" という項目と"七月末人数"という項目でパート比率が計算できますね。

やっぱり、宿泊飲食サービスのパート比率は77.64%ということで全産業分類の中で一番高いです。時給とパート比率の散布図を見てみましょう。plot関数です。

X軸が時給でY軸がパート比率の散布図です。パート比率が高いほど、時給が低いという逆相関ですね。cor関数で相関係数を計算しましょう。

相関係数は、-0.832とかなりの逆相関です。

この相関係数が統計的に有意かどうか、検定します。cor.test関数です。

p-value = 0.0001191 < 0.05 ですから、統計的に有意です。つまり、業界の時給はパート比率と関係があるということです。

毎月勤労統計調査の分析V2_4 - 大分類の業種での人数や給与総額を見てみる。(R言語でデータ整理)

今回は毎月勤労統計の2018年7月のデータを使って、産業ごとの値を見たいと思います。前回の分析では集約された産業分類も個々の産業分類も一緒になっていたので、いまひとつ納得感のある分析ではありませんでした。

まずは、毎月勤労統計がどのように産業分類されているのか確認しましょう。

「毎月勤労統計調査における記号の見方」というPDFファイルを見たら、大分類、中分類、特掲産業という3つのくくりがありました。

大分類はこうなっています。

TLは調査産業計ですから、実際の大分類は、Cの鉱業、採石業、砂利採取業からRのサービス業(他に分類されないもの)の16種類です。今回はこの16種類の産業分類で調べてみましょう。

ついでに規模はどのようになっているかというと、このようになっています。

形態もついでですから見ておきましょう。

規模や性別、就業形態です。

今回は大分類の16産業分類の規模は5人以上、性別は男女計で絞り込んで人数や給与総額などを見てみましょう。

まずは、CSVファイルに保存してあるデータをread.csvファイルでR言語に読込みます。head関数とsummary関数でどんなデータセットなのか確認しましょう。

業種コード、規模、性別、形態の4つの列で条件を絞り込めばいいのですね。

まず、業種コードが、C,D,E,F,G,H,I,O,P,Q,Rだけのデータフレームを作成します。

data$業種コード=="C" | data$業種コード=="D" 。。。と条件を | でつなげるとOR検索になりますね。summary関数でみてみましょう。

業種コードのところが、C : 40 、D : 40 などとなっていますね。

規模は5人以上だけに絞りましょう。

規模のところを見ると5人以上は75で、その他は0になっていますね。

性別を男女合計だけにしましょう。

性別のところが男女合計は45で、女、男は0になりました。形態を合計だけにしましょう。

これで条件が絞り込まれましたので、大分類の16業種、規模は5人以上、性別は男女合計、形態は合計になｔれいるはずです。全体を表示して確認しましょう。

最後の仕上げで、行番号(1列目)、業種コード(3列目)、規模(4列目)、性別(5列目)、形態(6列目)を削除してしまいましょう。

きれいになりましたが、一番左端、行番号が気に入らないです。これも直しましょう。

では、七月末人数の大きい順に並び替えてみましょう。

卸売小売が937万人で一番多く、製造業は808万人、宿泊飲食サービスは441万人と続きます。

barplot関数で棒グラフにしてみましょう。barplot関数です。

卸売小売と製造業がツートップという感じですね。

給与総額も並びかえましょう。

学術研究、専門・技術サービスが56万1677円でトップです。宿泊飲食サービスが13マン8864円で最下位です。

これも棒グラフにしましょう。

人数とは違った形ですね。こちらのそうが違いがあまりない感じです。

今回はここまでです。

毎月勤労統計調査の分析V2_3 - 人数の多い業種は？(order関数で並び替え)

引き続き、R言語で毎月勤労統計調査のデータを分析していきたいと思います。

データをread.csv関数で読込み、summary関数で最大値、最小値、平均値、中央値などを確認します。

とりあえず、七月末人数の多い順に並び替えてみましょうかね。order関数を使います。

「規模は5人以上、性別は男女合計、形態は合計」にデータを絞り込んでから、order関数で並び替えましょう。

data2 <- data[data$規模=="5人以上" & data$性別=="男女合計" & data$形態=="合計", ]

というコマンドで新しくdata2という「規模は5人以上かつ性別は男女合計かつ形態は合計」というデータだけのデータフレームを作成しています。

そして、data2[order(data2$七月末人数, decreasing=TRUE), ] というコマンドで七月末人数で大きい順に並びかえ、head関数で始めの20行を表示しています。このとき、c(2,4,5,6,10)と表示する列を2列目、4列目、5列目、6列目、10列目だけに指定しています。head(データフレーム, 20)としているので、始めの20のデータが表示されます。head関数は何も書かずに、head(データフレーム)だけだと始めの6行が表示されます。

一番上の調査産業計を見ると人数は約5000万人です。これだけの人数が日本では働いているということですね。

実際の業種ごとで見ると、1番目が、卸売小売で937万人ですね。2番目が製造業で808万人、3番目が763万人です。次が小売で609万人です。卸売小売は937万人で小売が609万人ですから、卸売は937-609=328万人かな？と思ったら、もうちょっと下のほうに卸売328万人ってありました。こうしてみると、業種が個別業種と集約業種が混ざっているのですね。個別の業種だけに絞りこむ必要があるようです。

個別業種だけに絞り込むのは次回以降にやってみたいと思います。

今回はここまでです。