探索的データ分析 - EDA

課題

以下の指標の中から、一つを選択して、データの概要（description）を記録し、データを WDI で取得し、以下の分析をする。

1. 各年毎のデータの数の棒グラフ
2. 日本のデータの年の降順での表示
3. 経年変化を表す折れ線グラフ
   1. 日本
   2. 南部アフリカ関税同盟の５カ国
   3. 選択したいくつかの国
4. データが十分ある最近の年の値のヒストグラム
5. データが十分ある最近の年の値の10カ国の値の棒グラフ
   1. 値が大きい方から
   2. 値が小さい方から

それぞれについて考察（気づいたこと、疑問など）を記す

2023.1.27. 23:59 までに Moodle の演習の課題ボックスに提出したものについては、なるべく、早く見て、フィードバックを書きます。それ以降に提出されたものも見ますが、フィードバックは遅くなると思ってください。

データ

1. Government expenditure on education, total (% of GDP)：SE.XPD.TOTL.GD.ZS [Link] 変数名：ed_exp

2. School enrollment, primary (% gross)：SE.PRM.ENRR [Link] 変数名：primary

3. School enrollment, secondary (% gross)：SE.SEC.ENRR [Link] 変数名：secondary

4. School enrollment, tertiary (% gross)：SE.TER.ENRR [Link] 変数名：tertiary

5. Mortality rate, under-5 (per 1,000 live births)：SH.DYN.MORT [Link] 変数名：under5

6. Incidence of HIV (% of uninfected population ages 15-49)：SH.HIV.INCD.ZS [Link] 変数名：new_hiv

7. School enrollment, primary and secondary (gross), gender parity index (GPI)：SE.ENR.PRSC.FM.ZS [Link] 変数名：school_gpi

8. Ratio of female to male labor force participation rate (%) (modeled ILO estimate)：SL.TLF.CACT.FM.ZS [Link] 変数名：job_gpi

9. Unemployment, female (% of female labor force) (modeled ILO estimate)：SL.UEM.TOTL.FE.ZS [Link] 変数名：female_unemploy

10. Unemployment, male (% of male labor force) (modeled ILO estimate)：SL.UEM.TOTL.MA.ZS [Link] 変数名：male_unemploy

11. Net official development assistance and official aid received (current US$) DT.ODA.ALLD.CD [Link] 変数名：oda

例

貧困率の時に例として学んだ、1日　6.85 ドル以下で生活している人の割合のデータを例にします。

コードは難しいものもありますが、にた分析をするときは、どこを変えれば良いか、考えてください。大体は、データの名前（df_under_6.85）と、変数名（under_6.85）を変更することになります。

データ

• Poverty headcount ratio at $6.85 a day (2017 PPP) (% of population)：SI.POV.UMIC [Link]

• Poverty headcount ratio at $6.85 a day is the percentage of the population living on less than $6.85 a day at 2017 international prices.

データ情報

データの情報を書きます。リンクの情報のマーク、または、概要（Description）から、データの説明をみることができます。分析には、どのようなデータか理解することが必要ですから、必ずこの作業をしてください。

• データ名：貧困率（1日6.85ドル）

• データコード：SI.POV.UMIC

• 変数名：under_6.85

• 概要：1日6.85ドルの貧困人口比率は、2017年の国際価格で1日6.85ドル未満で生活している人口の割合

データの取得

準備

これまでとも同じように二つのパッケージを読み込み（load）ます。

library(tidyverse)
library(WDI)

WDI パッケージを使って、直接データをダウンロードし、変数名を、under_6.85 に指定。

この変数名は、今後使いますから、重要です。一応、例として、わたしが使った、変数名を書いてあります。他の変数名を使っても構いません。

df_under_6.85 <- WDI(indicator = c(under_6.85 = "SI.POV.UMIC"))

data ディレクトリに書き込んでおきます。

write_csv(df_under_6.85, "data/under_6.85.csv")

data ディレクトリから読み出します。Run All や、Run All Chunks Above をする時は、上の二つをスキップする設定にしてあります。

df_under_6.85 <- read_csv("data/under_6.85.csv")

Rows: 16758 Columns: 5── Column specification ─────────────────────────────────────────────────────────────────
Delimiter: ","
chr (3): country, iso2c, iso3c
dbl (2): year, under_6.85
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

データの確認

データ名で、中身を表示できます。head(df_under_6.85) とすると、6行表示されます。df_under_6.85 |> head(10) とすると、10行表示します。あとで使っています。

df_under_6.85

変数の情報を得ることができます。

str(df_under_6.85)

spc_tbl_ [16,758 × 5] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
 $ country   : chr [1:16758] "Africa Eastern and Southern" "Africa Eastern and Southern" "Africa Eastern and Southern" "Africa Eastern and Southern" ...
 $ iso2c     : chr [1:16758] "ZH" "ZH" "ZH" "ZH" ...
 $ iso3c     : chr [1:16758] "AFE" "AFE" "AFE" "AFE" ...
 $ year      : num [1:16758] 2022 2021 2020 2019 2018 ...
 $ under_6.85: num [1:16758] NA NA NA NA NA NA NA NA NA NA ...
 - attr(*, "spec")=
  .. cols(
  ..   country = col_character(),
  ..   iso2c = col_character(),
  ..   iso3c = col_character(),
  ..   year = col_double(),
  ..   under_6.85 = col_double()
  .. )
 - attr(*, "problems")=<externalptr> 

データには country のところに、国だけでなく、地域も入っているので、地域のリストを、iso2c で REGION に入れておきます。以前にも一度使いました。

REGION <- c("1A", "1W", "4E", "7E", "8S", "B8", "EU", "F1", "OE", "S1", 
"S2", "S3", "S4", "T2", "T3", "T4", "T5", "T6", "T7", "V1", "V2", 
"V3", "V4", "XC", "XD", "XE", "XF", "XG", "XH", "XI", "XJ", "XL", 
"XM", "XN", "XO", "XP", "XQ", "XT", "XU", "XY", "Z4", "Z7", "ZF", 
"ZG", "ZH", "ZI", "ZJ", "ZQ", "ZT")

地域名にはどのようなものがあるか見ておきましょう。あとで、どこを調べるか考える時も、地域名や、国名がわかっていると、便利です。

df_under_6.85 |> filter(iso2c %in% REGION) |> distinct(country, iso2c)

国名も表示しておきます。

df_under_6.85 |> filter(!(iso2c %in% REGION)) |> distinct(country, iso2c)

分析する国のリスト

南部アフリカ関税同盟 The Southern African Customs Union (SACU)

これは、みなさんも使ってください。

SOUTH_AFRICA_FIVE <- c("South Africa", "Namibia", "Eswatini", "Botswana", "Lesotho")

こちらは、自分で分析したい国のりすとを作ってください。

ラテンアメリカでジニ指数が大きい４カ国

CHOSEN_COUNTRIES <- c("Suriname", "Belize", "Brazil", "Colombia")

分析

1. 各年毎のデータの数の棒グラフ

まずは、NA の値を削除します。そのあとで、国の情報の数を数えたいので、!(iso2c %in% REGION) で、上でおいた、地域の iso2c を選択し、! は否定でしたから、地域ではないものを選択し、その数を棒グラフにしています。

df_under_6.85 |> drop_na(under_6.85) |> filter(!(iso2c %in% REGION)) |>
  ggplot(aes(year)) + geom_bar()

視覚化

2. 日本の貧困率（1日6.85ドル以下）

日本のデータを選択するには、filter(country == "Japan") とします。そのあとは、under_6.85 の値が、NA のものを削除し、year の descending order 降順に並べ替えています。

df_under_6.85 |> filter(country == "Japan") |> 
  drop_na(under_6.85) |> arrange(desc(year))

3. 経年変化

a. 日本

まず、日本のデータを選択、次に、under_6.85 が NA のものを削除し、x 軸に、year、y 軸に、under_6.85 の値をとり、折れ線グラフで表示します。

df_under_6.85 |> filter(country == "Japan") |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line()

必ず、気づいたこと、疑問を書いてください。それが、EDA です。

気づいたこと・疑問

• あまりにもデータが少ないので、評価はできない。

• 0.5% から 2.2% 程度、この基準で貧困者がいるということは、精度を高めながら、継続して、調査する必要があるように思われる。

b. 南部アフリカ関税同盟

選択するときに、日本ではなく、SOUTH_AFRICA_FIVE に入っている、5カ国のどれかに一致するものを選ぶのが、country %in% SOUTH_AFRICA_FIVE です。国を色で区別しています。実は、ggplot(aes(year, under_6.85, col = country)) + geom_line() でも同じ結果になります。

df_under_6.85 |> filter(country %in% SOUTH_AFRICA_FIVE) |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line(aes(col = country))

参考：平均的な値を曲線で表すことも可能です。loess を使うと滑らかな曲線で近似してくれます。

色は、国別、それ以外に、geom_smooth で近似曲線を加えています。Help から、geom_smooth などを確認してください。

df_under_6.85 |> filter(country %in% SOUTH_AFRICA_FIVE) |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line(aes(col = country)) +
  geom_smooth(formula = 'y~x', method = "loess", se = FALSE)

気づいたこと・疑問

• 貧困率を、6.85ドルでとると、非常に多くのひとたちが、貧困になっていることがわかる。

c. ラテンアメリカ４カ国

選択して国についての分析です。

df_under_6.85|> filter(country %in% CHOSEN_COUNTRIES) |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line(aes(col = country))

参考：平均的な値を曲線で表すことも可能です。loess を使うと滑らかな曲線で近似してくれます。

df_under_6.85 |> filter(country %in% CHOSEN_COUNTRIES) |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line(aes(col = country)) +
  geom_smooth(formula = 'y~x', method = "loess", se = FALSE)

分布

日本も、2008, 2010, 2013 年にしかデータがなく、南部アフリカもこのデータは少ないので、データの数から、2010年について見てみる。

binwidth または、bins の数を変更して、適切なものに調整します。

df_under_6.85 |> filter(year == 2010) |> filter(!(country %in% REGION))|>
  drop_na(under_6.85) |>
  ggplot(aes(under_6.85)) + geom_histogram(binwidth = 10)

参考：SACU の５カ国の値を縦線で書き込むには下のようにします。

このデータの値はあまりないので、難しいですが、日本のデータと、南部アフリカの両方にデータがある年を選んでいます。

df_under_6.85 |> filter(year == 2010) |> filter(country %in% SOUTH_AFRICA_FIVE) 

参考：日本とSACU の５カ国の値を縦線で書き込むには下のようにします。

少し難しいですが、JP に 2010 年の日本の値を入れてあります。

JP <- 0.5
SAF <- df_under_6.85 |> filter(year == 2010) |> filter(country %in% SOUTH_AFRICA_FIVE) |> pull(under_6.85)
df_under_6.85 |> filter(year == 2010) |> filter(!(country %in% REGION))|>
  drop_na(under_6.85) |>
  ggplot() + geom_histogram(aes(under_6.85), binwidth = 10) +
  geom_vline(xintercept = SAF, col = "red") + geom_vline(xintercept = JP, col = "blue") +labs(title = "貧困率（1日6.85ドル以下）", subtitle = "日本：青、SACU：赤")

データが十分ある最近の年の値の10カ国の値の棒グラフ

a. 値が大きい方から

arrange(desc(under_6.85)) |> head(10)　これで、降順にならべ、その上から 10 個を取り出しています。fct_reorder(country, under_6.85) は、アルファベット順ではなく、値の大きい順に表示するためのものです。

df_under_6.85 |> filter(year == 2019) |> drop_na(under_6.85) |> 
  filter(!(iso2c %in% REGION))|>
  arrange(desc(under_6.85)) |> head(10) |> 
  ggplot(aes(fct_reorder(country, under_6.85), under_6.85)) + geom_col() + 
  coord_flip() + labs(title = "Top 10 Countries", x = "country", y = "poverty rate (under $6.85 per day)")

b. 値が小さい方から

fct_rev(fct_reorder(country, under_6.85) の部分は、順序を上のものとは、逆にしています。

df_under_6.85 |> filter(year == 2019) |> drop_na(under_6.85) |> 
  filter(!(iso2c %in% REGION))|>
  arrange(under_6.85) |> head(10) |> 
  ggplot(aes(fct_rev(fct_reorder(country, under_6.85)), under_6.85)) + geom_col() + 
  coord_flip() + labs(title = "Lowest 10 Countries", x = "country", y = "poverty rate (under $6.85 per day)")

国の教育に関する支出

概要：国内総生産（GDP）に対する、国の教育に関する支出（Government expenditure on education, total (% of GDP)）のデータの分析を行う

データ

• Government expenditure on education, total (% of GDP)：SE.XPD.TOTL.GD.ZS [Link]

データ情報

• データ名：

• データコード：

• 変数名：ed_exp

• 概要：

データの取得

準備

library(tidyverse)
library(WDI)

WDI パッケージを使って、直接データをダウンロードし、変数名を、ed_exp に指定。

df_under_6.85 <- WDI(indicator = c(under_6.85 = "SI.POV.UMIC"))

write_csv(df_under_6.85, "data/under_6.85.csv")

df_under_6.85 <- read_csv("data/under_6.85.csv")

データの確認

df_under_6.85

str(df_under_6.85)

REGION <- c("1A", "1W", "4E", "7E", "8S", "B8", "EU", "F1", "OE", "S1", 
"S2", "S3", "S4", "T2", "T3", "T4", "T5", "T6", "T7", "V1", "V2", 
"V3", "V4", "XC", "XD", "XE", "XF", "XG", "XH", "XI", "XJ", "XL", 
"XM", "XN", "XO", "XP", "XQ", "XT", "XU", "XY", "Z4", "Z7", "ZF", 
"ZG", "ZH", "ZI", "ZJ", "ZQ", "ZT")

df_under_6.85 |> filter(iso2c %in% REGION) |> distinct(country, iso2c)

df_under_6.85 |> filter(!(iso2c %in% REGION)) |> distinct(country, iso2c)

分析する国のリスト

南部アフリカ関税同盟 The Southern African Customs Union (SACU)

SOUTH_AFRICA_FIVE <- c("South Africa", "Namibia", "Eswatini", "Botswana", "Lesotho")

いくつかの国を選択

CHOSEN_COUNTRIES <- c("Suriname", "Belize", "Brazil", "Colombia")

分析

1. 各年毎のデータの数の棒グラフ

df_under_6.85 |> drop_na(under_6.85) |> filter(!(iso2c %in% REGION)) |>
  ggplot(aes(year)) + geom_bar()

気づいたこと・疑問

視覚化

2. 日本

df_under_6.85 |> filter(country == "Japan") |> 
  drop_na(under_6.85) |> arrange(desc(year))

気づいたこと・疑問

3. 経年変化

a. 日本

df_under_6.85 |> filter(country == "Japan") |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line()

気づいたこと・疑問

b. 南部アフリカ関税同盟

df_under_6.85 |> filter(country %in% SOUTH_AFRICA_FIVE) |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line(aes(col = country))

気づいたこと・疑問

参考：平均的な値を曲線で表すことも可能です。loess を使うと滑らかな曲線で近似してくれます。

df_under_6.85 |> filter(country %in% SOUTH_AFRICA_FIVE) |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line(aes(col = country)) +
  geom_smooth(formula = 'y~x', method = "loess", se = FALSE)

気づいたこと・疑問

c. 選択した国・地域

df_under_6.85 |> filter(country %in% CHOSEN_COUNTRIES) |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line(aes(col = country))

気づいたこと・疑問

参考：平均的な値を曲線で表すことも可能です。loess を使うと滑らかな曲線で近似してくれます。

df_under_6.85 |> filter(country %in% CHOSEN_COUNTRIES) |> drop_na(under_6.85) |>
  ggplot(aes(year, under_6.85)) + geom_line(aes(col = country)) +
  geom_smooth(formula = 'y~x', method = "loess", se = FALSE)

気づいたこと・疑問

分布

データの数から、まずは、2020年について見てみる。

df_under_6.85 |> filter(year == 2010) |> filter(!(country %in% REGION))|>
  drop_na(under_6.85) |>
  ggplot(aes(under_6.85)) + geom_histogram(binwidth = 10)

気づいたこと・疑問

参考：SACU の５カ国の値を縦線で書き込むには下のようにします。

df_under_6.85 |> filter(year == 2010) |> filter(country %in% SOUTH_AFRICA_FIVE) 

気づいたこと・疑問

参考：日本とSACU の５カ国の値を縦線で書き込むには下のようにします。

JP <- 0.5
SAF <- df_under_6.85 |> filter(year == 2010) |> filter(country %in% SOUTH_AFRICA_FIVE) |> pull(under_6.85)
df_under_6.85 |> filter(year == 2010) |> filter(!(country %in% REGION))|>
  drop_na(under_6.85) |>
  ggplot() + geom_histogram(aes(under_6.85), binwidth = 1) +
  geom_vline(xintercept = SAF, col = "red") + geom_vline(xintercept = JP, col = "blue") +labs(title = "貧困率（1日6.85ドル以下）", subtitle = "日本：青、SACU：赤")

気づいたこと・疑問

データが十分ある最近の年の値の10カ国の値の棒グラフ

a. 値が大きい方から

df_under_6.85 |> filter(year == 2019) |> drop_na(under_6.85) |> 
  filter(!(iso2c %in% REGION))|>
  arrange(desc(under_6.85)) |> head(10) |> 
  ggplot(aes(fct_reorder(country, under_6.85), ed_exp)) + geom_col() + 
  coord_flip() + labs(title = "Top 10 Countries", x = "country", y = "poverty rate (under $6.85 per day)")

気づいたこと・疑問

b. 値が小さい方から

df_under_6.85 |> filter(year == 2019) |> drop_na(under_6.85) |> 
  filter(!(iso2c %in% REGION))|>
  arrange(under_6.85) |> head(10) |> 
  ggplot(aes(fct_rev(fct_reorder(country, under_6.85)), under_6.85)) + geom_col() + 
  coord_flip() + labs(title = "Lowest 10 Countries", x = "country", y = "poverty rate (under $6.85 per day))

気づいたこと・疑問