crosshyou

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

国税庁の申告所得データの分析1 - 基本統計量

国税庁のホームページに長期時系列のデータがありました。

f:id:cross_hyou:20181215160655j:plain

この「申告所得税」のファイルをダウンロードしてみました。

f:id:cross_hyou:20181215160808j:plain

 

ファイルの4シート目の「所得種類別金額」のデータをR言語で分析してみようと思います。

read.csv関数でR言語に読込みやすいように、下図のようなCSVファイルにしました。

f:id:cross_hyou:20181215161051j:plain

 

それでは、read.csv関数でデータをR言語に読込み、summary関数で各変数の基本統計量をみてみましょう。ちなみに、金額はすべて百万円単位です。

f:id:cross_hyou:20181215161612j:plain

 

西暦の最小値は1949で最大値は2016なので1949年から2016年までのデータであることがわかります。NAを含む変数がいくつかあります、その他事業所得、利子所得、損益通算差額、分離長期譲渡所得、分離短期譲渡所得、株式譲渡所得などです。

変数がたくさんあるので、どれから分析したらいいでしょうか。。まずは「総計」からみてみましょうか。総計のヒストグラムをhist関数で描いてみます。

f:id:cross_hyou:20181215162231j:plain

f:id:cross_hyou:20181215162244j:plain

山が二つある感じですね。density関数とplot関数でカーネル密度グラフを描きます。

f:id:cross_hyou:20181215162602j:plain

f:id:cross_hyou:20181215162634j:plain

山が二つありますね。

小さい順に並び替えてグラフにしてみましょう。

f:id:cross_hyou:20181215162844j:plain

f:id:cross_hyou:20181215162856j:plain

右の2つのプロットが外れ値っぽいですね。

boxplot関数で箱ひげ図を描いてみましょう。

f:id:cross_hyou:20181215163059j:plain

f:id:cross_hyou:20181215163112j:plain

外れ値を表す丸がないので、2つの値は外れ値といえるほど外れてはいないようです。

このデータは時系列ですから、plot関数で時系列グラフを描いてみましょう。

f:id:cross_hyou:20181215163503j:plain

f:id:cross_hyou:20181215163516j:plain

1989年、90年のバブル期が一番金額が大きいですね。

2010年以降、回復しているようですが、バブル期に比べるとまだまだですね。

総計だけでsummary関数を実行してみましょう。

f:id:cross_hyou:20181215164006j:plain

 

最大値が、59兆1143億9800万円です。大きい順に並び替えて、何年に記録したのか確認しておきましょう。order関数を使います。

f:id:cross_hyou:20181215164420j:plain

 

1991年が最高の申告所得の年ですね。

申告所得の総計の動きはだいたいわかりました。各種類別の所得の動きをみてみましょう。

平均値を計算して、大きい順に表示してみましょう。apply関数とmean関数を使って平均値を出します。

f:id:cross_hyou:20181215164847j:plain

総計の次に大きいのは給与所得ですね。平均で9兆7299億4013万円です。

総計と給与所得の散布図をplot関数で表示してみましょう。

f:id:cross_hyou:20181215165715j:plain

f:id:cross_hyou:20181215165725j:plain

正の相関関係がありますね。

今回は以上です。

 

木質バイオマスエネルギー利用動向調査の分析3 - R言語で茨城、神奈川、千葉、佐賀とその他の都道府県の違いを確認する

 

www.crosshyou.info

 の続きです。

今回は事業所当りの木質バイオマス使用量を計算することからはじめましょう。

f:id:cross_hyou:20181215114257j:plain

 

全国の事業所の合計が1343事業所、全国の使用量合計が888万0772トンなので、1事業所当り6612トンの木質バイオマスを使用している計算になります。

都道府県ごとの値も計算しましょう。

f:id:cross_hyou:20181215114829j:plain

 

茨城の事業所は1事業所当り7万4988トンも使用しています。

少ないところはどこでしょうか?

f:id:cross_hyou:20181215115058j:plain

 

東京都がダントツで少ないですね。1事業所当り48トンです。

事業所当り使用量の各都道府県の値がどのように分布しているか、ヒストグラムで確認します。hist関数を使います。

f:id:cross_hyou:20181215115625j:plain

 

f:id:cross_hyou:20181215115635j:plain

 

1事業所当り2万トン以上使用する、規模の大きな都道府県が5つありますね。のこりは2万トン以下、大半は1万トン以下です。

plot関数とdensity関数で確率密度グラフで表示してみます。

f:id:cross_hyou:20181215120053j:plain

f:id:cross_hyou:20181215120107j:plain

 

小さい順に並べてグラフを書いてみましょう。

f:id:cross_hyou:20181215120609j:plain

 

f:id:cross_hyou:20181215120555j:plain

上位の4つ、茨城、神奈川、千葉、佐賀が外れ値っぽいですよね。boxplot関数で箱ひげ図を作成してみます。

f:id:cross_hyou:20181215121003j:plain

f:id:cross_hyou:20181215120911j:plain

外れ値を表す丸い図形が4つありますので、やっぱり茨城、神奈川、千葉、佐賀は外れ値ですね。

この4件とその他の43都道府県で事業所の構成比率や使用する材料の構成比率に違いがあるでしょうか?

まずは、4件だけのデータフレームを作りましょう。fourという名前にしましょう。

subset関数で作ります。

f:id:cross_hyou:20181215122514j:plain

この4件は木質ペレットは使用していないですね。

残りの43都道府県のデータフレームも作りましょう。othersという名前にします。

f:id:cross_hyou:20181215123410j:plain

このように、茨城県、神奈川県、千葉県、佐賀県が除外されているのが確認できます。

それでは4件の事業所の構成比をみてみましょう。colSums関数で3列目から5列目の合計値を計算します。

f:id:cross_hyou:20181215123937j:plain

発電機が8、ボイラーが9、両方が4となりました。

その他の都道府県はどうでしょうか?

f:id:cross_hyou:20181215124152j:plain

 

圧倒的にボイラーが多いですね、ボイラーの事業所は規模が小さいのでしょう。

両グループの事業所構成比率に統計的に有意な違いがあるかをカイ自乗検定で確認します。

まずは、そのためのマトリックスを作成します。matrix関数を使います。

f:id:cross_hyou:20181215124859j:plain

colnames関数で列名を付与して、rownames関数で行名を付与しました。

カイ自乗検定はchisq.test関数です。

f:id:cross_hyou:20181215125217j:plain

 

p-value = 5.42e-10 < 0.05 ですから、統計的に4県とその他では事業所の構成比率に差があるといえます。それでは、どの部分が有意に違っているでしょうか?

f:id:cross_hyou:20181215125945j:plain

 

chisq.test関数の結果を変数に格納して、$residualsでどの部分が有意なのかわかります。絶対値で1.96以上だと有意な違いとみます。こうしてみると、4県は発電機と両方が有意に多く、ボイラーが有意に少ないことがわかります。

同じように使用する材料の種類についても分析してみましょう。

使用量の合計値をcolSums関数で6列目から10列目を合計します。

f:id:cross_hyou:20181215130742j:plain

 

数値が大きくて違いがわかりにくいですね。パーセント表示にしてみましょう。

prop.table関数を使います。

f:id:cross_hyou:20181215131041j:plain

 

木質ペレットを使うかどうかと木粉の使用に差がある感じですね。それではカイ自乗検定をしてみましょう。

まずはマトリックスを作成します。

f:id:cross_hyou:20181215131904j:plain

chisq.test関数でカイ自乗検定を実行します。

f:id:cross_hyou:20181215131935j:plain

 

p-value < 2.2e-16 < 0.05ですから4県とその他の件の使用する材料の比率に違いがあると言えます。

どの部分が有意に違うのかを見てみましょう。

f:id:cross_hyou:20181215132047j:plain

 

すべての材料で有意に違いがありますね。4県を基準にして言うと、4県は木材チップ、木粉、その他が多く、木質ペレット、薪が少ないということです。

今回は以上です。

 

木質バイオマスエネルギー利用動向調査の分析2 - R言語で列の合計を算出する

 

www.crosshyou.info

 の続きです。

今回は事業所のタイプや使用材料のタイプごとの比率を算出してみましょう。

まずは、種類ごとの合計値をcolSums関数で算出してみます。

f:id:cross_hyou:20181213191040j:plain

総使用量が888万0772トンです。木材チップが773万4236トンで大半ですね。

割合を計算してみましょう。

f:id:cross_hyou:20181213191844j:plain

木材チップが87%と大半なのがわかります。

続いて事業所のタイプを見てみましょう。

f:id:cross_hyou:20181213192211j:plain

 

事業所総数が1343か所で、ボイラーが1175か所ですからボイラーが大半ですね。

ほとんどの事業所は木材チップをボイラーで燃やして熱を得ている、ということだと言えます。

割合も計算しておきましょう。

f:id:cross_hyou:20181213192535j:plain

ボイラーの事業所が87%ほどだとわかります。

今回は以上です。

 次回は

 

www.crosshyou.info

 

です。

 

木質バイオマスエネルギー利用動向調査の分析1 - R言語で基本統計量など。

今回からは、木質バイオマスエネルギー利用動向調査のデータをつかってR言語でのデータ分析の練習をしていきたいと思います。

政府統計の総合窓口(e-Stat)からデータを取得しました。

「本調査は、毎年、木質バイオマスのエネルギー利用動向を把握するため、木質バイオマスをエネルギー利用している事業所の概要、利用した設備の動向、公的補助の活用状況、利用した木質バイオマス量等の調査を行い、その動向について、全国、都道府県別に提供しています。」とのことです。

f:id:cross_hyou:20181212193129j:plain

いくつかファイルがありましたが、使用したファイルは、2016年のデータです。

f:id:cross_hyou:20181212193254j:plain

この都道府県別集計表の「木質バイオマスエネルギー利用機器の所有形態別事業所数」と「木質バイオマスの利用料」です。

この2つのファイルを合体させて以下のようなファイルにしました。

f:id:cross_hyou:20181212193542j:plain

都道府県別に、事業所数、発電機のある事業所数、ボイラーのある事業所数、両方あり事業所数、木材チップ使用量、木質ペレット使用量、薪使用量、木粉使用量、その他の使用量です。使用量はトン単位です。

それではこのファイルをread.csv関数で読込みましょう。全国の行が余計なので、この行は削除してしまいましょう。

f:id:cross_hyou:20181212195225j:plain

read.csv関数でCSVファイルのデータをR言語に読込み、df[-1, ]として最初の行、すなわち全国の行を削除しています。rownames関数で北海道の行が1になるように設定して、head関数ではじめの6行を表所しています。北海道の行が1となっていることがわかります。そして、summary関数で最小値、最大値、中央値、平均値などの基本的な統計量を算出しています。

summary関数では標準偏差が算出されませんので、apply関数とsd関数を組み合わせて各変数の標準偏差を出します。

f:id:cross_hyou:20181212200108j:plain

同じようにapply関数とmean関数を組み合わせて各変数の平均値を出します。

f:id:cross_hyou:20181212200555j:plain

 

こうして標準偏差と平均値が算出できましたから、変動係数(標準偏差 / 平均値)を計算してみます。

f:id:cross_hyou:20181212200959j:plain

データのバラツキが小さいのは事業所総数と木材チップです。そしてデータのバラツキが大きいのは木質ペレットです。

事業所総数のヒストグラムをhist関数で作成します。

f:id:cross_hyou:20181212202032j:plain

f:id:cross_hyou:20181212201915j:plain

 

150以上の事業種のある都道府県がありますね。大きい順のデータを並び替えみましょう。

f:id:cross_hyou:20181212202358j:plain

北海道が149か所、岩手が99か所、高知が79か所、宮崎が62か所、秋田県が55か所、岐阜県が52か所となっています。

事業所数が少ないところはどこでしょうか?東京かな?

f:id:cross_hyou:20181212203728j:plain

 

沖縄が1か所、佐賀が2か所、神奈川が3か所、東京が4か所、千葉が5か所、大阪が6か所となっています。

使用量についても調べてみましょう。

もともとのデータは各原料データのみで総使用量が無いので、作成しましょう。rowSums関数で木材チップ(6列目)からその他(10列目)の値を合計します。

f:id:cross_hyou:20181212205053j:plain

はい、こんな感じです。それでは総使用量で並び替えてみましょう。

f:id:cross_hyou:20181212205606j:plain

北海道が一番ではないんですね。茨城が1番で82万4872トンです。福島、静岡、北海道、宮崎、秋田と続きます。

小さい順ではどうでしょうか?

f:id:cross_hyou:20181212205908j:plain

東京がダントツで小さいですね。195トンしか使用していないです。埼玉、山梨、長崎、和歌山、香川と続きます。

総使用量の平均値や標準偏差、変動係数も算出しておきましょう。さきほどと同じくapply関数とsd関数、mean関数を使います。

f:id:cross_hyou:20181212210248j:plain

総使用量の平均値は18万8952トンで変動係数は0.06です。

次は、hist関数でヒストグラムを作図します。

f:id:cross_hyou:20181212210804j:plain

f:id:cross_hyou:20181212210833j:plain

一番小さいレンジが一番度数が多いですね。

事業所総数と総使用量の散布図をplot関数で作成しましょう。

f:id:cross_hyou:20181212211103j:plain

f:id:cross_hyou:20181212211119j:plain

プラスの相関がありそうですね。cor関数で相関係数を算出します。

f:id:cross_hyou:20181212211322j:plain

 

0.4097ですね。

今回は以上です。

 次回は

 

www.crosshyou.info

 

です。

読書記録 - 「抗生物質と人間 マイクロバイオームの危機」 山本太郎著 岩波新書

 

 抗生物質を使うと、人間の体内にある有益な微生物も死滅してしまうので、抗生物質の使用は極力控えなければならない。

病気は多くはある細菌の存在によって引き起こされるが、ある細菌が存在しなくなるとなる病気もある。

家畜は抗生物質を与えると成長がよくなるので、病気予防ではなく、成長促進の目的で抗生物質を投与している国がまだある、日本と米国はそうで、ヨーロッパは成長促進目的の投与を禁止している。

 

経済サンセスの事業所に関する集計データの分析4 - R言語で1事業所当りの従業員数を調べる

 

www.crosshyou.info

 に引き続き、経済サンセスのデータをR言語で分析します。

前回までは男女比率に注目していましたが、今回は1事業所当りの従業員数を調べてみましょう。

まずは1事業所当りの従業員数の算出に必要なデータを確認しましょう。

f:id:cross_hyou:20181211171429j:plain

 

subset関数で必要な変数だけを抽出しました。

総従業者数 / 事業所数を計算して1事業所当りの従業員数を見てみましょう。

f:id:cross_hyou:20181211172130j:plain

 

order関数を使って大きい順に並び替えてみましょう。

f:id:cross_hyou:20181211172642j:plain

 

関東が人数が多いですね。2016年が12.83人です。一番少ないのは2012年の新潟で9.1人ですね。

2016年のヒストグラムを描いてみます。hist関数です。

f:id:cross_hyou:20181211174533j:plain

 

f:id:cross_hyou:20181211174546j:plain

 

subset関数で2016年の事業所当り従業員数のみのデータフレームを作成し、hist関数でヒストグラムを作成しています。ヒストグラムの引数はデータフレームでは動かず、ベクトルでないと動かないので、empl16[ , 1]としてベクトルにしています。

同じように2012年もやってみましょう。

f:id:cross_hyou:20181211175509j:plain

 

f:id:cross_hyou:20181211175521j:plain

 

う~~ん、2016年と2012年でどっちがどうなのか、わかりにくいですね。

一つの画面でヒストグラムを描いてみます。

f:id:cross_hyou:20181211220841j:plain

 

f:id:cross_hyou:20181211220859j:plain

 

こうしてみると、2012年のほうが左側に分布しているようにみえます。

それぞれの平均値や中央値をsummary関数で確認します。

f:id:cross_hyou:20181211221211j:plain

最小値、第1分位値、中央値、平均値、第3分位値、最大値すべて2016年のほうが大きい値です。

それでは、この分布位置の違いが統計的に有意かどうかをウィルコクソン=マン・ホイットニー検定をします。今回は2012年と2016年で同じ地域で比べていますので、paired = TRUEを指定します。

f:id:cross_hyou:20181211221735j:plain

 

p-value = 0.0006104 < 0.05 ですから、2012年と比較して2016年のほうが、事業所当りの従業員数は増えているといえます。

今回は以上です。

 

経済サンセスの事業所に関する集計データの分析3 - R言語で各地域の男女比率をカイ二乗検定する

 

www.crosshyou.info

 の続きです。

前回の分析では、全地域を合計した男性従業者数と女性従業者数の比率が2012年と2016年では変化していて、カイ二乗検定の結果、その変化は統計的に有意なものであること、各地域別では松山都市圏のみ男女比率が増加(男性比率増加)していたこと、がわかりました。

今回は各地域ごとにカイ二乗検定をして男女比率が変化しているのが統計的に有意かどうかを調べてみたいと思います。

地域が14あるので、14の地域についてクロス表を作成してカイ二乗検定をすればいいです。一つ一つクロス表を手作業で作成するのではなく、for関数で自動化してみたいと思います。

まず、クロス表なデータを確認しておきましょう。

f:id:cross_hyou:20181211080409j:plain

 

男性従業者数は5列目、女性従業者数は6列目です。

札幌のクロス表を[行, 列]のインデックス座標形式で表現するとこうなります。

[1, 5]   [1, 6]
[15, 5] [15, 6]

札幌の次の仙台のクロス表はこうなります。

[2, 5]    [2, 6]
[16, 5] [16, 6]
札幌と仙台のクロス表の構造を比べると行の値が1増加、列の値は5、6で変化がないことがわかります。なのでこのクロス表の構造は

[i, 5]          [i, 6]
[i + 14, 5] [i + 14, 6]
i = 1 ~ 14
という構造になります。

これをR言語のfor関数で再現すればいいわけですね。

具体的にはこうなります。

f:id:cross_hyou:20181211083336j:plain

 

mという変数にクロス表を格納、mresultという変数にカイ二乗の結果を格納、regionnameという変数に格納して、それらをprint関数で出力します。

それでは結果をみてみましょう。

f:id:cross_hyou:20181211084150j:plain

 

まずは、札幌、仙台、関東です。赤く記しをつけたところがp値ですが3地域のp値は0.05以下なので統計的に有意な違いが2012年と2016年にある、ということです。

f:id:cross_hyou:20181211084424j:plain

新潟、静岡浜松、中京、近畿の4地域もp値が0.05以下ですので、有意です。

f:id:cross_hyou:20181211084726j:plain

 

岡山、広島、北九州福岡、熊本の各地域でもp値は0.05以下で有意です。

f:id:cross_hyou:20181211085137j:plain

 

宇都宮、松山、鹿児島の各地域のp値も0.05以下です。松山都市圏の男性比率の増加は有意な増加だったとわかります。

今回は以上です。

 次回は

 

www.crosshyou.info

 

です。