GES 001 演習5

準備

1. （自分のPCまたは教室のPCに）ログイン

2. ウェッブ・ブラウザー（Google Chrome など）を起動

   • Moodle の GES001 経済と経済学のサイトから、このスライドのページを表示（リンク「Rでデータサイエンス」https://ds-sl.github.io/intro2r/ges001/）

3. （別のタブまたは ウィンドウで）PositCloud にログイン[Posit.cloud]

   • アカウントのない人はサイン・アップ [共有プロジェクト] から、Save a Permanent Copy）

   • RStudio を自分のコンピュータにインストールしている人は起動

4. リンクの右上の Raw ボタンの右の Copy a raw file からコピーして演習用 R Markdown ファイルを作成（あとで再度解説します）[Rmd]

ファイル

ファイルを作成するときの注意

1. Posit Cloud のときは、まず、Login し、intro2rj のプロジェクトに入り、 ファイルから、ges001 のフォルダーを選択して、それを選択して、移動します。この中に、data フォルダが作成されていることを確認し、このges001 に、新しく作成したファイルを保存します。新しく作成したファイルの入っているフォルダーの中に、data フォルダがあることが大切です。そこに、データを書き込みます。

2. RStudio の場合には、まずは、あたらしい　Project を作成します。File > New Project から作成します。すでに、作成してある場合は、それを、Open Project や、Recent Project から開きます。その中に、新しいファイルを作成します。作成したフォルダーに、data フォルダがあることを確認してください。新しく作成したファイルの入っているフォルダーの中に、data フォルダがあることが大切です。そこに、データを書き込みます。

---

ファイルを提出するときの注意

RStudio の場合には、自分の PC に作成したファイルがありますから、問題ないと思いますが、Posit Cloud で作成した場合には、提出したいファイルの左にあるチェックボックスをチェックします。Files の 右端にある、ギアマークの Export を押すと、ダウンロードできます。それを提出してください。末尾が、nb.html となっているものを提出していただくのがよいですが、よくわからないときは、nb.html ファイルと、Rmd ファイルと両方提出してください、

第5週

01/18(TH)　南部アフリカ諸国の貧困と不平等に対する対策

　　　　　　COVID-19が貧困に与えた影響 　　　　

• 機会の不平等と成果の不平等が相互に関係し、長期にわたって持続している。

• 人生の初期に出現する不平等を最小化するための政策が重要である。

• 要因：教育、労働市場、失業率、男女間格差、土地所有権、などなど　　

01/23(TU)　Rでデータサイエンス5　 [Main]

演習 1月24日（火）

内容

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

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

• School enrollment, primary (% gross)：SE.PRM.ENRR [Link]

• School enrollment, secondary (% gross)：SE.SEC.ENRR [Link]

• School enrollment, tertiary (% gross)：SE.TER.ENRR [Link]

---

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

• Incidence of HIV (% of uninfected population ages 15-49)：SH.HIV.INCD.ZS [Link]

• School enrollment, primary and secondary (gross), gender parity index (GPI)：SE.ENR.PRSC.FM.ZS [Link]

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

• Unemployment, female (% of female labor force) (modeled ILO estimate)：SL.UEM.TOTL.FE.ZS [Link]

• Unemployment, male (% of male labor force) (modeled ILO estimate)：SL.UEM.TOTL.MA.ZS [Link]

• Net official development assistance and official aid received (current US$) DT.ODA.ALLD.CD [Link]

---

準備

library(tidyverse)
library(WDI)

---

データの読み込み（importing）

df_poverty_inequality <- WDI(
  indicator = c(gini = "SI.POV.GINI",
                under_6.85 = "SI.POV.UMIC",
                ed_exp = "SE.XPD.TOTL.GD.ZS",
                primary = "SE.PRM.ENRR",
                secondary = "SE.SEC.ENRR",
                tertiary = "SE.TER.ENRR",
                under5 = "SH.DYN.MORT",
                new_hiv = "SH.HIV.INCD.ZS",
                school_gpi = "SE.ENR.PRSC.FM.ZS",
                job_gpi = "SL.TLF.CACT.FM.ZS",
                female_unemploy = "SL.UEM.TOTL.FE.ZS",
                male_unemploy = "SL.UEM.TOTL.FE.ZS",
                oda = "DT.ODA.ALLD.CD"), extra = TRUE)

---

保存と読み込み

２回目からは、data から読み込めるようにしておきます。

最初の１回目は、かならず実行してください。

write_csv(df_poverty_inequality, "data/poverty_inequality.csv")

df_poverty_inequality <- read_csv("data/poverty_inequality.csv")

Rows: 16758 Columns: 24── Column specification ─────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr   (7): country, iso2c, iso3c, region, capital, income, lending
dbl  (15): year, gini, under_6.85, ed_exp, primary, secondary, tertiary, under5, new_hiv, sch...
lgl   (1): status
date  (1): lastupdated
ℹ 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.

---

データを見てみよう (viewing)

df_poverty_inequality 

---

変数の選択（selecting）

df_pov_ineq <- df_poverty_inequality |> 
  select(country, iso2c, year, gini:oda, region, income, lending)
df_pov_ineq

---

変数同士の相関関係: NA ではないところのみ選択

df_pov_ineq |> drop_na(gini:oda) |> select(gini:oda)

---

cor(cars$speed,cars$dist)

[1] 0.8068949

cars |> ggplot(aes(speed, dist)) + geom_point() + 
  geom_smooth(formula = 'y~x',method = "lm", se=FALSE)

相関係数：直線の傾きが正なら正、負なら負、直線に近い程、1 または-1 に近い

---

変数相互の相関係数を一度に求める

df_pov_ineq |> drop_na(gini:oda) |> select(gini:oda) |> cor() |> 
  round(digits = 2) |> as.data.frame()

---

---

df_pov_ineq |> drop_na(gini:oda) |>
  ggplot(aes(gini, primary)) + geom_point(aes(col = income)) + 
  geom_smooth(formula = 'y~x', method = "lm", se = FALSE) +
  labs(title = "cor(gini,primary) = 0.34")

---

df_pov_ineq |> drop_na(gini:oda) |>
  ggplot(aes(gini, tertiary)) + geom_point(aes(col = income)) + 
  geom_smooth(formula = 'y~x', method = "lm", se = FALSE) +
  labs(title = "cor(gini,tertiary) = -0.29")

---

df_pov_ineq |> drop_na(gini:oda) |>
  ggplot(aes(under_6.85, tertiary)) + geom_point(aes(col = income)) + 
  geom_smooth(formula = 'y~x', method = "lm", se = FALSE) +
  labs(title = "cor(under_6.85,tertiary) = -0.77")

---

df_pov_ineq |> drop_na(gini:oda) |>
  ggplot(aes(under_6.85, under5)) + geom_point(aes(col = income)) + 
  geom_smooth(formula = 'y~x', method = "lm", se = FALSE) +
  labs(title = "cor(under_6.85,under5) = 0.71")

ファイルリンクと課題

基本的には、PositCloud（https://posit.cloud/）を使って実習

• 探索的データ分析（EDA） -
  • 練習と一つ目の課題（w5eda.Rmd） [リンク], [Rmd]
  • 二つ目以降の課題：w5eda1.Rmd [リンク], [Rmd]
• 課題：2023.1.27. 23:59 までに Moodle の演習の課題ボックスに提出したものについては、なるべく、早く見て、フィードバックを書きます。それ以降に提出されたものも見ますが、フィードバックは遅くなると思ってください。

R Notebook のソースファイル（Rmd）を取得方法

• Rmd のリンクをクリックし Raw の横の Copy a raw file からコピー、新規 RMarkdown ファイルを PositCloud または RStudio のプロジェクト内に作成しペースト

課題

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

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

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

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

---

提出前の確認

• Preview で確認。

• Web Browser で、w5_c123456.nb.html など、R Notebook を見て確認。

• もし、問題があれば、Run ボタンの右の三角から、Run All を選択し、エラーがでないか確認。

• 最初にもどる。

---

途中でのエラー

• 入力したときには、例を参照して、スペルなどを確認してください。全角になっていると問題がおきます。() がペアでマッチしているか、確認してください。
• 引用符が入っていなかったり、== のところが、= だったり、いろいろな可能性があります。Error message を読むこともたいせつです。エラーがでた、Code Chunk と、Error message を、ChatGPT や、Google Bard, Google Search に入れると、解決方法を教えてくれることもあります。
• File not found の、エラーがでたときには、上から順に、Run （Code Chunk の右上の三角印を押して実行）してみてください。または、エラーが出たところに、カーソルを置き、上の、Run ボタンの右の三角から、Run All Chunks Above を選択すると、そこまでのすべての　Code Chunk を実行してくれます。
• 上の方法でうまくいかないときは、data フォルダに、データ（***.csv）が入っているかを確認、なければ、data フォルダがあることを確認して、最初のデータ読み込みのところを実行してみてください。
• 実行できていても、結果が見えないこともあります。そのときは、Code Chunk の下にある、山二つの記号を押してみてください。これは、結果を表示、非表示にします。それが原因で隠れている場合があります。

---

データ（全体で11：最初のリスト参照）

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

例：国の教育に関する支出

スライドで見ている方は、RNotebook ファイルで見てください [リンク]

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

データ

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

データ情報

• データ名：

• データコード：

• 変数名：

• 概要：

---

データの取得

準備

library(tidyverse)
library(WDI)

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

df_ed_exp <- WDI(indicator = c(ed_exp = "SE.XPD.TOTL.GD.ZS"))

write_csv(df_ed_exp, "data/ed_exp.csv")

df_ed_exp <- read_csv("data/ed_exp.csv")

Rows: 16758 Columns: 5── Column specification ─────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr (3): country, iso2c, iso3c
dbl (2): year, ed_exp
ℹ 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.

---

データの確認

df_ed_exp

str(df_ed_exp)

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 ...
 $ ed_exp : num [1:16758] 3.91 4.63 4.35 4.54 4.74 ...
 - attr(*, "spec")=
  .. cols(
  ..   country = col_character(),
  ..   iso2c = col_character(),
  ..   iso3c = col_character(),
  ..   year = col_double(),
  ..   ed_exp = col_double()
  .. )
 - attr(*, "problems")=<externalptr> 

---

国と地域

country には、国と地域両方が入っています。地域の iso2c は以下のものです。

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_ed_exp |> filter(iso2c %in% REGION) |> distinct(country, iso2c)

df_ed_exp |> filter(!(iso2c %in% REGION)) |> distinct(country, iso2c)

---

分析する国のリスト

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

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

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

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

分析

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

df_ed_exp |> drop_na(ed_exp) |> filter(!(iso2c %in% REGION)) |>
  ggplot(aes(year)) + geom_bar()

視覚化

2. 日本の教育費（% of GDP）

df_ed_exp |> filter(country == "Japan") |> 
  drop_na(ed_exp) |> arrange(desc(year))

3. 経年変化

a. 日本

df_ed_exp |> filter(country == "Japan") |> drop_na(ed_exp) |>
  ggplot(aes(year, ed_exp)) + geom_line()

気づいたこと・疑問

• 1970年代の急激な上昇、1990年ごろの急激な現象は、何が原因なのだろう。

• 2014年ごろから減少、2018年ごろから増加、2020年から2021年は減少。

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

df_ed_exp |> filter(country %in% SOUTH_AFRICA_FIVE) |> drop_na(ed_exp) |>
  ggplot(aes(year, ed_exp)) + geom_line(aes(col = country))

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

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

気づいたこと・疑問

• 平均で見ると、上昇してきており、7% 程度という大きな割合になっている。

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

df_ed_exp |> filter(country %in% CHOSEN_GINI_COUNTRIES) |> drop_na(ed_exp) |>
  ggplot(aes(year, ed_exp)) + geom_line(aes(col = country))

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

df_ed_exp |> filter(country %in% CHOSEN_GINI_COUNTRIES) |> drop_na(ed_exp) |>
  ggplot(aes(year, ed_exp)) + geom_line(aes(col = country)) +
  geom_smooth(formula = 'y~x', method = "loess", se = FALSE)

---

分布

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

df_ed_exp |> filter(year == 2020) |> filter(!(country %in% REGION))|>
  drop_na(ed_exp) |>
  ggplot(aes(ed_exp)) + geom_histogram(binwidth = 1)

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

df_ed_exp |> filter(year == 2020) |> filter(country %in% SOUTH_AFRICA_FIVE) 

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

JP <- 3.416981
SAF <- df_ed_exp |> filter(year == 2020) |> filter(country %in% SOUTH_AFRICA_FIVE) |> pull(ed_exp)
df_ed_exp |> filter(year == 2020) |> filter(!(country %in% REGION))|>
  drop_na(ed_exp) |>
  ggplot() + geom_histogram(aes(ed_exp), binwidth = 1) +
  geom_vline(xintercept = SAF, col = "red") + geom_vline(xintercept = JP, col = "blue") +labs(title = "2020年の教育費の対GDP百分率", subtitle = "日本：青、SACU：赤")

---

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

a. 値が大きい方から

df_ed_exp |> filter(year == 2020) |> drop_na(ed_exp) |> 
  filter(!(iso2c %in% REGION))|>
  arrange(desc(ed_exp)) |> head(10) |> 
  ggplot(aes(fct_reorder(country, ed_exp), ed_exp)) + geom_col() + 
  coord_flip() + labs(title = "Top 10 Countries", x = "country", y = "Government expenditure on education, total (% of GDP)")

b. 値が小さい方から

df_ed_exp |> filter(year == 2020) |> drop_na(ed_exp) |> 
  filter(!(iso2c %in% REGION))|>
  arrange(ed_exp) |> head(10) |> 
  ggplot(aes(fct_rev(fct_reorder(country, ed_exp)), ed_exp)) + geom_col() + 
  coord_flip() + labs(title = "Lowest 10 Countries", x = "country", y = "Government expenditure on education, total (% of GDP)")

---

演習の内容と課題

基本的には、PositCloud（https://posit.cloud/）を使って実習

• 探索的データ分析（EDA） -
  • 練習と一つ目の課題（w5eda.Rmd） [リンク], [Rmd]
  • 二つ目以降の課題：w5eda1.Rmd [リンク], [Rmd]
• 課題：2023.1.27. 23:59 までに Moodle の演習の課題ボックスに提出したものについては、なるべく、早く見て、フィードバックを書きます。それ以降に提出されたものも見ますが、フィードバックは遅くなると思ってください。

参考文献

1. 「みんなのデータサイエンス - Data Science for All」[はじめてのデータサイエンス]

   • 導入として、GDP（国内総生産）のデータを使って説明しています。

2. Posit Recipes（旧 Posit Primers）: The Basics 対話型の演習サイトの最初 [Link]

3. Posit Cheat Sheet. 早見表です。印刷して使うために、PDF も提供しています。[Site Link]

4. DataCamp Cheat Sheet: Tidyverse for Biginners. データサイエンスの教育をしている会社の早見表の一つです。基本が簡単にまとまっています。[Link]