杜 依濛
長崎大学　経済学部

この練習ファイルは、Rの基本的なコーディング操作に基づいた授業内演習用の資料です。 演習課題では、システムデータおよびオリジナルデータを使用します。 コーディングを始める前に、初期設定を行い、必要なパッケージをインストールしてください。

  1. 注意事項
  1. コーディング練習

1 RとRstudioのインストール（初期設定）

参考資料

重要：まず参考資料を確認してください。インストール前に必要な作業がある。

1.1 Rのインストール

ダウンロード

インストール手順：

  • インストールパスに日本語やスペースを含めないよう確認、インストール前の作業を必ず確認してください。
  • ダウンロードしたインストーラーをダブルクリックして実行
  • 画面の指示に従って進める

1.2 RStudioのインストール

重要：Rを先にインストールしてから、RStudioをインストールしてください。

ダウンロードページ

1.3 RStudioの環境構築

「Tools」→「Global Options」

CRANミラーの設定

「ツール」→ 「グローバルオプション」 → 「パッケージ」 → 「CRANミラー」

「Tools」→「Global Options」→「Packages」→「Management」→「Change」(Primary CRAN repository)

文字化けを防ぐためのテキストエンコーディング設定

「デフォルトのテキストエンコーディング（Default Text Encoding）」 を 「UTF-8」 に変更

「Tools」→「Global Options」→「Code」→「Saving」→「Change」→「UTF-8」

デフォルトの作業ディレクトリの設定（必要に応じて）

RStudio を起動するたびに、毎回同じフォルダを作業ディレクトリにしたい場合は、次のように設定

「Tools」→「Global Options」→「General」→「Browse」(Default Working Directory)

1.4 パッケージのインストール

よく使うパッケージのインストール

  • パッケージをインストール：install.packages(“パッケージ名XXX”)
  • パッケージを読み込んでRで使えるようになる：library(パッケージ名XXX)
  • 以下のコードをRスクリプトにコピペして実行してください
# tidyverse：データ操作や可視化、読み込みなどの必須パッケージの集合体。主な関数例：filter(), select(), mutate(), summarize(), arrange(), group_by(), ggplot()など。

# "dependencies = TRUE"は関連するすべてのパッケージをインストールする意味
install.packages("tidyverse", dependencies = TRUE)

# tinytex：R MarkdownでPDFを出力するためのLaTeX環境を簡単に構築
install.packages("tinytex", dependencies = TRUE)

# tinytexをインストール（小さなエラーは無視してOK、失敗してもスキップ可能）
tinytex::install_tinytex()

# rticles：学術論文用のR Markdownテンプレート集
install.packages("rticles", dependencies = TRUE)

# readxl：Excelファイル（.xls/.xlsx）読み込み, 関数：read_excel()
install.packages("readxl", dependencies = TRUE)

# writexl：Excelファイルを書き出し, 関数：write_xlsx()
install.packages("writexl", dependencies = TRUE)

# psych：多変量解析用関数
install.packages("psych", dependencies = TRUE)

# stargazer: モデルを自動で読み取ってテーブルを出力, 関数：stargazer(m1, m2, type = "latex")
install.packages("stargazer", dependencies = TRUE)

# modelsummary: モデルを自動で読み取ってテーブルを出力, 関数：modelsummary(list(m1, m2), output = "latex")
install.packages("modelsummary", dependencies = TRUE)

# wooldridge：データセットコレクションパッケージ wooldridgeを使う
install.packages("wooldridge", dependencies = TRUE)


# よく使われる他の便利なパッケージ集
install.packages(c(
  "sf", "hrbrthemes", "patchwork", "manipulateWidget",
  "ggthemes", "tidyquant", "rvest", "DT", "basetheme",
  "pacman", "corrplot", "formattable", "paletteer", 
  "reprex"), 
  dependencies = TRUE)
  • インストール済みパッケージ数の確認
library(tidyverse)
installed.packages() %>% 
  as_tibble() %>% 
  count() %>% 
  pull() %>% 
  paste("このパソコンには", ., "個のRパッケージがインストールされている！") %>% 
  message()

1.5 Rスクリプト

  1. 「File」 → 「New File」 → 「R Script」

2. スクリプトを書いてみましょう

library(tidyverse)
library(readxl)
library(writexl)
library(psych)
library(stargazer)
library(modelsummary)
library(ggplot2)
library(wooldridge)

data(wage1)
head(wage1)

ggplot(data=wage1, aes(x=wage)) +
  geom_histogram()

  1. スクリプトファイルを名前付けて保存

任意のディレクトリに保存しても構わないが、後から見つけやすくするために、自分で作成した「プロジェクトフォルダ（英語名）」内に保存するのがおすすめ

ショートカット：
Windows: Ctrl + S
Mac:     Cmd + S

1.6 その他の設定 (必要に応じて)

Rprofileの活用

Rprofile を設定する目的：Rprofile は RStudio起動時に自動で読み込まれる設定ファイル、よく使うパッケージやオプションをここに書いておくと、毎回手動で設定する手間が省ける

# usethisパッケージを使って .Rprofile を開く
usethis::edit_r_profile()
  • R Profileテンプレートの設定例

以下の内容を .Rprofileファイルに記述したら設定完了、Rstudio起動時に設定を自動ロードし、必要に応じてあとから自由に編集できる

#　パッケージがインストール済みであれば一括呼び出す
pacman::p_load(tidyverse, psych, ggplot2, dplyr, writexl, readxl)

# 科学的記数法を使わず、通常の数値で表示
options(scipen = 1)

# ggplot2のカラーパレットをviridisに設定（色覚バリアフリー対応）
options(
  ggplot2.continuous.colour = "viridis",
  ggplot2.continuous.fill = "viridis"
)

フォント設定（日本語 + 英語）

  1. お気に入りのフォントを毎回便利に使えるように覚えやすい名前をつける
  2. 自分のPCにあるフォントをRで使う
  3. .Rprofileに書いて起動時に自動設定されるようにするのがおすすめ
  • Windowsユーザー
windowsFonts(
  <覚えやすい名前> = windowsFont("<Windowsにインストールされた実フォント名>")
)

# 設定例：
windowsFonts(
  enfont = "Avenir"   # 英語フォント（実フォント名に置き換えて）
　jpfont = "YuGothic" # 日本語フォント（実フォント名に置き換えて）
)
  • Macユーザー
# 設定例：
enfont = "Avenir" # 英語フォント（実フォント名に置き換えて）
jpfont = "YuGothic" # 日本語フォント（実フォント名に置き換えて）

2 Rを使った練習：データ処理

2.1 関連パッケージ

library(wooldridge)
library(tidyverse)
library(readxl)
library(writexl)
library(psych)
library(ggplot2)

2.2 Rのデータ構造

データにある１つの要素としての変数のタイプ

  • <-、->はそれぞれ左向き代入と右向き代入を意味する
# 数値型（Numeric）
a <- 3
class(a)
## [1] "numeric"
b <- 4
a+b
## [1] 7
3 -> a
4 -> b
a+b
## [1] 7
# 文字列型（Character）
b <- "2"
b
## [1] "2"
class(b)
## [1] "character"
# 因子型（Factor）
c <- c( "Bad", "Good", "Excellent" )
c
## [1] "Bad"       "Good"      "Excellent"
class(c)
## [1] "character"
# 他のタイプ（日付など）
library(lubridate)
d <- ymd(Sys.Date())
d
## [1] "2025-10-06"
class(d)
## [1] "Date"

Rにあるデータの構造

  • ベクトル(Vector)はRの最も基本的なデータ構造で、1種類のデータを持つ1次元データ構造
# 数値型（Numeric）
v <- c(1, 2, 3)
v
## [1] 1 2 3
class(v)
## [1] "numeric"
# 文字列型（Character）
vc <- as.character(v)
vc
## [1] "1" "2" "3"
class(vc)
## [1] "character"
# 日付型
vd <- ymd(c("2020-01-01", "2020-01-02", "2020-01-03"))
vd
## [1] "2020-01-01" "2020-01-02" "2020-01-03"
class(vd)
## [1] "Date"
  • ベクトル(Vector)の操作例
a <- c(1, 2, 3)
b <- c(2, 3, 4, 5, 6, 7)

# 関数や計算に`( )`を使い、データの抽出に`[ ]`を使う場合が多い
# bの3番目の要素を取り出す
b[ 3 ]
## [1] 4
# bの2番目から4番目を取り出す
b[ 2:4 ]
## [1] 3 4 5
# bの1, 3, 5番目の要素を取り出す
b[c( 1, 3, 5 )]
## [1] 2 4 6
# cを昇順に並び替え
c <- c(2, 3, 1)
sort(c)
## [1] 1 2 3
# ベクトルaの平均値を計算
mean(a)
## [1] 2
  • 行列(Matrix)は、1種類のデータを持つ2次元データ構造
#  1:20　→ 1から20までの整数を生成
#  matrix(1:20, c(5, 4), ...) →　20個の数を 5行×4列 の行列に並べる。デフォルトでは 列方向（column-wise） にデータが埋められる。
#  dimnames = list(...) → 行と列の名前をつける。list(行名ベクトル, 列名ベクトル) の形。行名は "r1"〜"r5", 列名は "c1"〜"c4"。

m <- matrix(1:20, c(5, 4), 
            dimnames = list(c("r1", "r2", "r3", "r4", "r5"), 
                            c("c1", "c2", "c3", "c4"))) 
m
##    c1 c2 c3 c4
## r1  1  6 11 16
## r2  2  7 12 17
## r3  3  8 13 18
## r4  4  9 14 19
## r5  5 10 15 20
# 1行2列の値を抽出
m[1, 2]
## [1] 6
# "r1"行"c2"列の値を抽出
m["r1", "c2"]
## [1] 6
# 第1〜３行の3列目の値を抽出
m[1:3, 3]
## r1 r2 r3 
## 11 12 13
# 第1行・第3行の3列目の値を抽出
m[c(1, 3), 3]
## r1 r3 
## 11 13
  • 配列(Array)は、1種類データを持つ多次元のデータ構造
arr <- array(1:24, c(2, 3, 4), 
      dimnames = list(c("A1", "A2"),　　　　　# 第1軸（行）
                      c("B1", "B2", "B3"),　　# 第2軸（列）
                      c("C1", "C2", "C3", "C4"))) # 第3軸（層）
arr
## , , C1
## 
##    B1 B2 B3
## A1  1  3  5
## A2  2  4  6
## 
## , , C2
## 
##    B1 B2 B3
## A1  7  9 11
## A2  8 10 12
## 
## , , C3
## 
##    B1 B2 B3
## A1 13 15 17
## A2 14 16 18
## 
## , , C4
## 
##    B1 B2 B3
## A1 19 21 23
## A2 20 22 24
# 第3軸の1番目の層を抽出
arr[,,1]　
##    B1 B2 B3
## A1  1  3  5
## A2  2  4  6
class(arr[,,1])
## [1] "matrix" "array"
  • データフレーム(data.frame)はテーブル形式のデータ集合を扱い、複数種類データを持つ2次元のデータ構造 (リストは各種のデータ型を格納する複雑なデータ構造で、その扱い方は演習に含まれていないため、紹介を省略する)

2.3 練習課題

練習：行列mの “r1”行・“r3”行の3列目の値をどのように抽出するか?

m <- matrix(1:20, c(5, 4), 
            dimnames = list(c("r1", "r2", "r3", "r4", "r5"), 
                            c("c1", "c2", "c3", "c4"))) 
m
##    c1 c2 c3 c4
## r1  1  6 11 16
## r2  2  7 12 17
## r3  3  8 13 18
## r4  4  9 14 19
## r5  5 10 15 20
# "r1"行・"r3"行の3列目値を抽出
m[c("r1", "r3"), 3]
## r1 r3 
## 11 13
# "r1"行・"r3"行の"c3"列目値を抽出
m[c("r1", "r3"), "c3"]
## r1 r3 
## 11 13

2.4 データフレームを確認

  • データパッケージ wooldridge に含まれるデータセット wage1 を読み込む
# パッケージ名: wooldridge、インストールされてない場合まずインストールする
install.packages('wooldridge') 
# `wooldridge` パッケージのデータ一覧を表示
data(package = "wooldridge")

# パッケージ `wooldridge` の中のデータセット `wage1`と`intdef`を読み込む
library(wooldridge)
data(wage1)

class(wage1)
## [1] "data.frame"
# データの最初の数行だけをざっと見る
head(wage1)
##   wage educ exper tenure nonwhite female married numdep smsa northcen south
## 1 3.10   11     2      0        0      1       0      2    1        0     0
## 2 3.24   12    22      2        0      1       1      3    1        0     0
## 3 3.00   11     2      0        0      0       0      2    0        0     0
## 4 6.00    8    44     28        0      0       1      0    1        0     0
## 5 5.30   12     7      2        0      0       1      1    0        0     0
## 6 8.75   16     9      8        0      0       1      0    1        0     0
##   west construc ndurman trcommpu trade services profserv profocc clerocc
## 1    1        0       0        0     0        0        0       0       0
## 2    1        0       0        0     0        1        0       0       0
## 3    1        0       0        0     1        0        0       0       0
## 4    1        0       0        0     0        0        0       0       1
## 5    1        0       0        0     0        0        0       0       0
## 6    1        0       0        0     0        0        1       1       0
##   servocc    lwage expersq tenursq
## 1       0 1.131402       4       0
## 2       1 1.175573     484       4
## 3       0 1.098612       4       0
## 4       0 1.791759    1936     784
## 5       0 1.667707      49       4
## 6       0 2.169054      81      64
data(intdef) 
# データの最初の数行だけをざっと見る
head(intdef)
##   year   i3  inf  rec  out        def i3_1 inf_1      def_1        ci3 cinf
## 1 1948 1.04  8.1 16.2 11.6 -4.6000004   NA    NA         NA         NA   NA
## 2 1949 1.10 -1.2 14.5 14.3 -0.1999998 1.04   8.1 -4.6000004 0.06000006 -9.3
## 3 1950 1.22  1.3 14.4 15.6  1.2000008 1.10  -1.2 -0.1999998 0.12000000  2.5
## 4 1951 1.55  7.9 16.1 14.2 -1.9000006 1.22   1.3  1.2000008 0.32999992  6.6
## 5 1952 1.77  1.9 19.0 19.4  0.3999996 1.55   7.9 -1.9000006 0.22000003 -6.0
## 6 1953 1.93  0.8 18.7 20.4  1.6999989 1.77   1.9  0.3999996 0.15999997 -1.1
##        cdef y77
## 1        NA   0
## 2  4.400001   0
## 3  1.400001   0
## 4 -3.100001   0
## 5  2.300000   0
## 6  1.299999   0
  • $はデータフレームから特定の変数を抽出する演算子、データセット wage1 から特定の変数 wage を確認する
# 他の変数も試してみてください
# head()：先頭の6個の要素を表示する関数（head(x, 10)とすれば先頭10個を表示できる）
head(wage1$wage)
## [1] 3.10 3.24 3.00 6.00 5.30 8.75
# unique()：その列（またはベクトル）に含まれる重複しない値を取り出す関数。同じ値が複数あっても、1 回だけ表示される
head(unique(wage1$wage))
## [1] 3.10 3.24 3.00 6.00 5.30 8.75
  • データフレームから特定の行や列を抽出する方法
# 1行目を抽出
wage1[1, ]
##   wage educ exper tenure nonwhite female married numdep smsa northcen south
## 1  3.1   11     2      0        0      1       0      2    1        0     0
##   west construc ndurman trcommpu trade services profserv profocc clerocc
## 1    1        0       0        0     0        0        0       0       0
##   servocc    lwage expersq tenursq
## 1       0 1.131402       4       0
# 1列目を抽出
wage1[ , 1]
##   [1]  3.10  3.24  3.00  6.00  5.30  8.75 11.25  5.00  3.60 18.18  6.25  8.13
##  [13]  8.77  5.50 22.20 17.33  7.50 10.63  3.60  4.50  6.88  8.48  6.33  0.53
##  [25]  6.00  9.56  7.78 12.50 12.50  3.25 13.00  4.50  9.68  5.00  4.68  4.27
##  [37]  6.15  3.51  3.00  6.25  7.81 10.00  4.50  4.00  6.38 13.70  1.67  2.93
##  [49]  3.65  2.90  1.63  8.60  5.00  6.00  2.50  3.25  3.40 10.00 21.63  4.38
##  [61] 11.71 12.39  6.25  3.71  7.78 19.98  6.25 10.00  5.71  2.00  5.71 13.08
##  [73]  4.91  2.91  3.75 11.90  4.00  3.10  8.45  7.14  4.50  4.65  2.90  6.67
##  [85]  3.50  3.26  3.25  8.00  9.85  7.50  5.91 11.76  3.00  4.81  6.50  4.00
##  [97]  3.50 13.16  4.25  3.50  5.13  3.75  4.50  7.63 15.00  6.85 13.33  6.67
## [109]  2.53  9.80  3.37 24.98  5.40  6.11  4.20  3.75  3.50  3.64  3.80  3.00
## [121]  5.00  4.63  3.00  3.20  3.91  6.43  5.48  1.50  2.90  5.00  8.92  5.00
## [133]  3.52  2.90  4.50  2.25  5.00 10.00  3.75 10.00 10.95  7.90  4.72  5.84
## [145]  3.83  3.20  2.00  4.50 11.55  2.14  2.38  3.75  5.52  6.50  3.10 10.00
## [157]  6.63 10.00  2.31  6.88  2.83  3.13  8.00  4.50  8.65  2.00  4.75  6.25
## [169]  6.00 15.38 14.58 12.50  5.25  2.17  7.14  6.22  9.00 10.00  5.77  4.00
## [181]  8.75  6.53  7.60  5.00  5.00 21.86  8.64  3.30  4.44  4.55  3.50  6.25
## [193]  3.85  6.18  2.91  6.25  6.25  9.05 10.00 11.11  6.88  8.75 10.00  3.05
## [205]  3.00  5.80  4.10  8.00  6.15  2.70  2.75  3.00  3.00  7.36  7.50  3.50
## [217]  8.10  3.75  3.25  5.83  3.50  3.33  4.00  3.50  6.25  2.95  5.71  3.00
## [229] 22.86  9.00  8.33  3.00  5.75  6.76 10.00  3.00  3.50  3.25  4.00  2.92
## [241]  3.06  3.20  4.75  3.00 18.16  3.50  4.11  1.96  4.29  3.00  6.45  5.20
## [253]  4.50  3.88  3.45 10.91  4.10  3.00  5.90 18.00  4.00  3.00  3.55  3.00
## [265]  8.75  2.90  6.26  3.50  4.60  6.00  2.89  5.58  4.00  6.00  4.50  2.92
## [277]  4.33 18.89  4.28  4.57  6.25  2.95  8.75  8.50  3.75  3.15  5.00  6.46
## [289]  2.00  4.79  5.78  3.18  4.68  4.10  2.91  6.00  3.60  3.95  7.00  3.00
## [301]  6.08  8.63  3.00  3.75  2.90  3.00  6.25  3.50  3.00  3.24  8.02  3.33
## [313]  5.25  6.25  3.50  2.95  3.00  4.69  3.73  4.00  4.00  2.90  3.05  5.05
## [325] 13.95 18.16  6.25  5.25  4.79  3.35  3.00  8.43  5.70 11.98  3.50  4.24
## [337]  7.00  6.00 12.22  4.50  3.00  2.90 15.00  4.00  5.25  4.00  3.30  5.05
## [349]  3.58  5.00  4.57 12.50  3.45  4.63 10.00  2.92  4.51  6.50  7.50  3.54
## [361]  4.20  3.51  4.50  3.35  2.91  5.25  4.05  3.75  3.40  3.00  6.29  2.54
## [373]  4.50  3.13  6.36  4.68  6.80  8.53  4.17  3.75 11.10  3.26  9.13  4.50
## [385]  3.00  8.75  4.14  2.87  3.35  6.08  3.00  4.20  5.60 10.00 12.50  3.76
## [397]  3.10  4.29 10.92  7.50  4.05  4.65  5.00  2.90  8.00  8.43  2.92  6.25
## [409]  6.25  5.11  4.00  4.44  6.88  5.43  3.00  2.90  6.25  4.34  3.25  7.26
## [421]  6.35  5.63  8.75  3.20  3.00  3.00 12.50  2.88  3.35  6.50 10.38  4.50
## [433] 10.00  3.81  8.80  9.42  6.33  4.00  2.90 20.00 11.25  3.50  6.00 14.38
## [445]  6.36  3.55  3.00  4.50  6.63  9.30  3.00  3.25  1.50  5.90  8.00  2.90
## [457]  3.29  6.50  4.00  6.00  4.08  3.75  3.05  3.50  2.92  4.50  3.35  5.95
## [469]  8.00  3.00  5.00  5.50  2.65  3.00  4.50 17.50  8.18  9.09 11.82  3.25
## [481]  4.50  4.50  3.71  6.50  2.90  5.60  2.23  5.00  8.33  2.90  6.25  4.55
## [493]  3.28  2.30  3.30  3.15 12.50  5.15  3.13  7.25  2.90  1.75  2.89  2.90
## [505] 17.71  6.25  2.60  6.63  3.50  6.50  3.00  4.38 10.00  4.95  9.00  1.43
## [517]  3.08  9.33  7.50  4.75  5.65 15.00  2.27  4.67 11.56  3.50
# 1～3行目、1～3列目を抽出
wage1[1:3, 1:3]
##   wage educ exper
## 1 3.10   11     2
## 2 3.24   12    22
## 3 3.00   11     2
# 1～3行目、'wage'列だけを抽出
wage1[1:3, "wage"]
## [1] 3.10 3.24 3.00
# 1～3行目、'wage'と'educ'列を抽出、'c()'は列名をまとめて指定するための関数
wage1[1:3, c("wage", "educ")]
##   wage educ
## 1 3.10   11
## 2 3.24   12
## 3 3.00   11

2.5 練習課題

練習：wage1から1～10行目のwageとexper列を抽出してwage1_1に代入してみましょう。

2.6 データをインポート

RとRstudioがパソコンにインストール済みの場合

# 現在の作業ディレクトリを確認（今どこで作業しているかを確認）
getwd()

# 作業ディレクトリを設定（データファイルが保存されているフォルダのパスを指定）
# Windowsの例： setwd("C:/Users/YourName/Desktop/folder")
# Macの例：     setwd("/Users/YourName/Desktop/folder")
setwd("")

#手動でSession>Set Working Directory>Choose Directoryで作業ディレクトリを設定することも可能


# ---------------------------------------------
# データファイルを作業ディレクトリに入れる
# readrパッケージのread_excel関数を使って、xlsxやxlsファイルを読み込む
# データファイル名を正しく指定してください
# データファイル`令和2年国勢調査都道府県.xlsx`を読み込み、`population`という名前として保存する
# ---------------------------------------------
population <- read_excel("令和2年国勢調査都道府県.xlsx") 
## New names:
## • `15歳未満（人）` -> `15歳未満（人）...15`
## • `15～64歳（人）` -> `15～64歳（人）...16`
## • `65歳以上（人）` -> `65歳以上（人）...17`
## • `15歳未満（％）` -> `15歳未満（％）...18`
## • `15～64歳（％）` -> `15～64歳（％）...19`
## • `65歳以上（％）` -> `65歳以上（％）...20`
## • `15歳未満（人）` -> `15歳未満（人）...21`
## • `15～64歳（人）` -> `15～64歳（人）...22`
## • `65歳以上（人）` -> `65歳以上（人）...23`
## • `15歳未満（％）` -> `15歳未満（％）...24`
## • `15～64歳（％）` -> `15～64歳（％）...25`
## • `65歳以上（％）` -> `65歳以上（％）...26`
## • `15歳未満（人）` -> `15歳未満（人）...27`
## • `15～64歳（人）` -> `15～64歳（人）...28`
## • `65歳以上（人）` -> `65歳以上（人）...29`
## • `15歳未満（％）` -> `15歳未満（％）...30`
## • `15～64歳（％）` -> `15～64歳（％）...31`
## • `65歳以上（％）` -> `65歳以上（％）...32`
## • `一般世帯（世帯）` -> `一般世帯（世帯）...37`
## • `一般世帯（世帯）` -> `一般世帯（世帯）...40`
# 必要に応じてデータフレームに変換し、データの内容を確認する
population %>% 
  data.frame() %>% 
  head()
##   都道府県名 都道府県.市区町村名 都道府県.市区町村名_英語
## 1    00_全国          00000_全国                    Japan
## 2  01_北海道        01000_北海道                 Hokkaido
## 3  01_北海道        01100_札幌市              Sapporo-shi
## 4  01_北海道  01101_札幌市中央区      Sapporo-shi Chuo-ku
## 5  01_北海道    01102_札幌市北区      Sapporo-shi Kita-ku
## 6  01_北海道    01103_札幌市東区   Sapporo-shi Higashi-ku
##   市などの別.._地域識別コード      総数       男       女
## 1                           a 126146099 61349581 64796518
## 2                           a   5224614  2465088  2759526
## 3                           1   1973395   918682  1054713
## 4                           0    248680   112853   135827
## 5                           0    289323   136596   152727
## 6                           0    265379   126023   139356
##   X2015年.平成27年.の人口.組替..人. X5年間の人口増減数.人.
## 1                         127094745                -948646
## 2                           5381733                -157119
## 3                           1952356                  21039
## 4                            237627                  11053
## 5                            285321                   4002
## 6                            261912                   3467
##   X5年間の人口増減率... 面積.参考. 人口密度.人.km2.       平均年齢.歳.
## 1  -0.74641000000000002  377976.41            338.2 47.627789999999997
## 2   -2.9194900000000001   83424.44             66.6 49.782380000000003
## 3               1.07762    1121.26           1760.0 47.722880000000004
## 4    4.6514100000000003      46.42           5357.2           46.57555
## 5               1.40263      63.57           4551.3 46.706919999999997
## 6    1.3237300000000001      56.97           4658.2 46.764470000000003
##       年齢中位数.歳. X15歳未満.人....15 X15.64歳.人....16 X65歳以上.人....17
## 1           48.50685           15031602          75087865           36026632
## 2 51.342840000000002             556526           2988800            1679288
## 3 48.491399999999999             215386           1208858             549151
## 4 46.773389999999999              23869            165634              59177
## 5 47.418979999999998              33526            177388              78409
## 6 47.222630000000002              29912            165626              69841
##    X15歳未満......18   X15.64歳......19  X65歳以上......20 X15歳未満.人....21
## 1 11.916029999999999 59.524520000000003 28.559449999999998            7700036
## 2 10.651999999999999 57.206139999999998 32.141860000000001             284897
## 3 10.914490000000001 61.257779999999997 27.827729999999999             110206
## 4 9.5982800000000008 66.605279999999993           23.79645              12161
## 5           11.58774 61.311410000000002 27.100850000000001              17127
## 6 11.271430000000001 62.411119999999997 26.317460000000001              15190
##   X15.64歳.人....22 X65歳以上.人....23  X15歳未満......24   X15.64歳......25
## 1          38008577           15640968 12.551080000000001 61.954090000000001
## 2           1477750             702441           11.55728 59.947150000000001
## 3            580412             228064            11.9961           63.17877
## 4             77465              23227 10.775969999999999 68.642390000000006
## 5             86609              32860           12.53843           63.40522
## 6             81649              29184           12.05336 64.788970000000006
##    X65歳以上......26 X15歳未満.人....27 X15.64歳.人....28 X65歳以上.人....29
## 1 25.494820000000001            7331566          37079288           20385664
## 2           28.49558             271629           1511050             976847
## 3 24.825130000000001             105180            628446             321087
## 4           20.58164              11708             88169              35950
## 5 24.056339999999999              16399             90779              45549
## 6 23.157679999999999              14722             83977              40657
##    X15歳未満......30   X15.64歳......31  X65歳以上......32           人口性比
## 1           11.31475 57.224200000000003           31.46105 94.680369999999996
## 2 9.8433200000000003           54.75759 35.399090000000001 89.330119999999994
## 3 9.9723799999999994           59.58455 30.443069999999999 87.102559999999997
## 4 8.6197900000000001 64.912719999999993 26.467490000000002 83.085840000000005
## 5           10.73746 59.438740000000003 29.823799999999999 89.438019999999995
## 6 10.564310000000001 60.260770000000001           29.17492 90.432419999999993
##   日本人.人. 外国人.人. 総世帯.世帯. 一般世帯.世帯....37 施設等の世帯.世帯.
## 1  123398962    2747137     55830154            55704949             125205
## 2    5188441      36173      2476846             2469063               7783
## 3    1959523      13872       969161              967372               1789
## 4     246027       2653       141429              141223                206
## 5     286191       3132       139675              139449                226
## 6     263485       1894       131188              130904                284
##   X2015年.平成27年.の世帯数.組替..世帯. 一般世帯.世帯....40
## 1                              53448685            55704949
## 2                               2444810             2469063
## 3                                921837              967372
## 4                                132006              141223
## 5                                133662              139449
## 6                                124425              130904
##   うち..核家族世帯.世帯. 夫婦のみの世帯.世帯. 夫婦と子供..から成る世帯.世帯.
## 1               30110571             11158840                       13949190
## 2                1324406               584819                         511571
## 3                 496691               204909                         204131
## 4                  54941                24467                          20812
## 5                  71323                28025                          30992
## 6                  65521                25570                          27792
##   男親と子供..から成る世帯.世帯. 女親と子供..から成る世帯.世帯.
## 1                         738006                        4264535
## 2                          29921                         198095
## 3                          10175                          77476
## 4                            878                           8784
## 5                           1592                          10714
## 6                           1406                          10753
##   うち..単独世帯.世帯. うち..65歳以上の単独世帯.世帯.
## 1             21151042                        6716806
## 2               999825                         361735
## 3               422160                         121789
## 4                80425                          16447
## 5                60279                          16690
## 6                59006                          16049
##   X.再掲...夫65歳以上.妻60歳以上の夫婦のみの世帯.世帯. X.再掲...3世代世帯.世帯.
## 1                                              6533895                  2337703
## 2                                               345741                    59601
## 3                                               110890                    15727
## 4                                                11059                     1111
## 5                                                15765                     2573
## 6                                                13364                     2071

RStudio Cloudを使用する場合 (https://posit.cloud/)

# データファイル（例：rent_0521.csv）をプロジェクトにアップロードするには、
# RStudio Cloud の画面右下にある「Files」タブを開いてください。

# 手順：
# 1. 「Files」タブをクリック
# 2. 「Upload」ボタン（黄色い矢印）をクリック
# 3. 「Choose File」を選択し、パソコンから‘rent_0521.csv’を選ぶ
# 4. 「OK」を押すとアップロード完了
# 5.  アップされたcsvファイルを選択し、「Import」ボタンを押すと、Rデータオブジェクトとして読み込める

# アップされたcsvファイルを選択し、「Import」ボタンを押すと、Rデータオブジェクトとして読み込める
# Importした場合、オブジェクトデータ名はファイル名と同じ'rent_0521'になる
# データ名を「Name」項目での変更可能（'rent_0521'→'rent'）

2.7 データの内容を確認

<bdl>や<chr>という記号が出てくるが、これは変数の種類を示す。<bdl>だと数値、<int>だと整数、<chr>は文字情報であることを示す

データwage1

# Environmentにあるデータ名をクリックでデータ自体を確認できる
# `wage1'を確認
glimpse(wage1)
## Rows: 526
## Columns: 24
## $ wage     <dbl> 3.10, 3.24, 3.00, 6.00, 5.30, 8.75, 11.25, 5.00, 3.60, 18.18,…
## $ educ     <int> 11, 12, 11, 8, 12, 16, 18, 12, 12, 17, 16, 13, 12, 12, 12, 16…
## $ exper    <int> 2, 22, 2, 44, 7, 9, 15, 5, 26, 22, 8, 3, 15, 18, 31, 14, 10, …
## $ tenure   <int> 0, 2, 0, 28, 2, 8, 7, 3, 4, 21, 2, 0, 0, 3, 15, 0, 0, 10, 0, …
## $ nonwhite <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ female   <int> 1, 1, 0, 0, 0, 0, 0, 1, 1, 0, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1, 1…
## $ married  <int> 0, 1, 0, 1, 1, 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0…
## $ numdep   <int> 2, 3, 2, 0, 1, 0, 0, 0, 2, 0, 0, 0, 2, 0, 1, 1, 0, 0, 3, 0, 0…
## $ smsa     <int> 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ northcen <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ south    <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ west     <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ construc <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ ndurman  <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ trcommpu <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ trade    <int> 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ services <int> 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ profserv <int> 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 1…
## $ profocc  <int> 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 0, 1, 0, 0, 1…
## $ clerocc  <int> 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0…
## $ servocc  <int> 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0…
## $ lwage    <dbl> 1.1314021, 1.1755733, 1.0986123, 1.7917595, 1.6677068, 2.1690…
## $ expersq  <int> 4, 484, 4, 1936, 49, 81, 225, 25, 676, 484, 64, 9, 225, 324, …
## $ tenursq  <int> 0, 4, 0, 784, 4, 64, 49, 9, 16, 441, 4, 0, 0, 9, 225, 0, 0, 1…
glimpse(wage1$wage)
##  num [1:526] 3.1 3.24 3 6 5.3 ...

データpopulation

# `population'を確認
population <- read_excel("令和2年国勢調査都道府県.xlsx") 
## New names:
## • `15歳未満（人）` -> `15歳未満（人）...15`
## • `15～64歳（人）` -> `15～64歳（人）...16`
## • `65歳以上（人）` -> `65歳以上（人）...17`
## • `15歳未満（％）` -> `15歳未満（％）...18`
## • `15～64歳（％）` -> `15～64歳（％）...19`
## • `65歳以上（％）` -> `65歳以上（％）...20`
## • `15歳未満（人）` -> `15歳未満（人）...21`
## • `15～64歳（人）` -> `15～64歳（人）...22`
## • `65歳以上（人）` -> `65歳以上（人）...23`
## • `15歳未満（％）` -> `15歳未満（％）...24`
## • `15～64歳（％）` -> `15～64歳（％）...25`
## • `65歳以上（％）` -> `65歳以上（％）...26`
## • `15歳未満（人）` -> `15歳未満（人）...27`
## • `15～64歳（人）` -> `15～64歳（人）...28`
## • `65歳以上（人）` -> `65歳以上（人）...29`
## • `15歳未満（％）` -> `15歳未満（％）...30`
## • `15～64歳（％）` -> `15～64歳（％）...31`
## • `65歳以上（％）` -> `65歳以上（％）...32`
## • `一般世帯（世帯）` -> `一般世帯（世帯）...37`
## • `一般世帯（世帯）` -> `一般世帯（世帯）...40`
glimpse(population)
## Rows: 1,965
## Columns: 49
## $ 都道府県名                                                   <chr> "00_全国", …
## $ `都道府県・市区町村名`                                       <chr> "00000_全国…
## $ `都道府県・市区町村名_英語`                                  <chr> "Japan", …
## $ `市などの別\r\n_地域識別コード`                              <chr> "a", "a",…
## $ 総数                                                         <chr> "12614609…
## $ 男                                                           <chr> "61349581…
## $ 女                                                           <chr> "64796518…
## $ `2015年（平成27年）の人口（組替）（人）`                     <chr> "12709474…
## $ `5年間の人口増減数（人）`                                    <dbl> -948646, …
## $ `5年間の人口増減率（％）`                                    <chr> "-0.74641…
## $ `面積(参考)`                                                 <chr> "377976.4…
## $ `人口密度（人/km2）`                                         <dbl> 338.2, 66…
## $ `平均年齢（歳）`                                             <chr> "47.62778…
## $ `年齢中位数（歳）`                                           <chr> "48.50685…
## $ `15歳未満（人）...15`                                        <chr> "15031602…
## $ `15～64歳（人）...16`                                        <chr> "75087865…
## $ `65歳以上（人）...17`                                        <chr> "36026632…
## $ `15歳未満（％）...18`                                        <chr> "11.91602…
## $ `15～64歳（％）...19`                                        <chr> "59.52452…
## $ `65歳以上（％）...20`                                        <chr> "28.55944…
## $ `15歳未満（人）...21`                                        <chr> "7700036"…
## $ `15～64歳（人）...22`                                        <chr> "38008577…
## $ `65歳以上（人）...23`                                        <chr> "15640968…
## $ `15歳未満（％）...24`                                        <chr> "12.55108…
## $ `15～64歳（％）...25`                                        <chr> "61.95409…
## $ `65歳以上（％）...26`                                        <chr> "25.49482…
## $ `15歳未満（人）...27`                                        <chr> "7331566"…
## $ `15～64歳（人）...28`                                        <chr> "37079288…
## $ `65歳以上（人）...29`                                        <chr> "20385664…
## $ `15歳未満（％）...30`                                        <chr> "11.31475…
## $ `15～64歳（％）...31`                                        <chr> "57.22420…
## $ `65歳以上（％）...32`                                        <chr> "31.46105…
## $ 人口性比                                                     <chr> "94.68036…
## $ `日本人（人）`                                               <chr> "12339896…
## $ `外国人（人）`                                               <chr> "2747137"…
## $ `総世帯（世帯）`                                             <chr> "55830154…
## $ `一般世帯（世帯）...37`                                      <chr> "55704949…
## $ `施設等の世帯（世帯）`                                       <chr> "125205",…
## $ `2015年（平成27年）の世帯数（組替）（世帯）`                 <chr> "53448685…
## $ `一般世帯（世帯）...40`                                      <chr> "55704949…
## $ `うち\r\n核家族世帯（世帯）`                                 <chr> "30110571…
## $ `夫婦のみの世帯（世帯）`                                     <chr> "11158840…
## $ `夫婦と子供\r\nから成る世帯（世帯）`                         <chr> "13949190…
## $ `男親と子供\r\nから成る世帯（世帯）`                         <chr> "738006",…
## $ `女親と子供\r\nから成る世帯（世帯）`                         <chr> "4264535"…
## $ `うち\r\n単独世帯（世帯）`                                   <chr> "21151042…
## $ `うち\r\n65歳以上の単独世帯（世帯）`                         <chr> "6716806"…
## $ `（再掲）\r\n夫65歳以上，妻60歳以上の夫婦のみの世帯（世帯）` <chr> "6533895"…
## $ `（再掲）\r\n3世代世帯（世帯）`                              <chr> "2337703"…

2.8 データの基本操作

パイプ %>% は、tidyverseパッケージを呼び出すことで利用できるパイプと呼ばれる演算子で、「左のものを右の関数に入れる」という意味になる、WindowsではControl+Shift+M、MacではCommand+Shift+Mで入力できる

データwage1

# ---------------------------------------------
# まず`wage1`でやってみる
# ---------------------------------------------
# 変数名を１列で表示
ls(wage1)
##  [1] "clerocc"  "construc" "educ"     "exper"    "expersq"  "female"  
##  [7] "lwage"    "married"  "ndurman"  "nonwhite" "northcen" "numdep"  
## [13] "profocc"  "profserv" "services" "servocc"  "smsa"     "south"   
## [19] "tenure"   "tenursq"  "trade"    "trcommpu" "wage"     "west"
# ---------------------------------------------
# `select()`関数は列の選択に利用、この関数を用いて一部の変数を取り出す
# `選んだ変数をnew_wage`に代入し、新しいデータセットとして保存する
# ---------------------------------------------
new_wage <- wage1 %>% 
  select(exper, female, married, nonwhite, educ, wage)
# ---------------------------------------------
# `filter()`関数は行の選択に利用、この関数を用いて一部の変数を取り出す
# `new_wage_1`というデータセットとして保存する
# 条件を満たすデータだけを取り出し、`new_wage_1`というデータセットとして保存する
# ---------------------------------------------
 new_wage_1 <- new_wage %>% 
   filter(exper > 5)

# うまくいく場合下記のNOTEを無視してください
# ---------------------------------------------
# NOTE!
# tidyverseを読み込むとdplyr::filter() が使えるようになるが、
# R Markdownではchunkの実行順や読み込みエラーによって、
# 他の関数が優先されてしまうことがある。
# その場合、"object 'xxx' not found" のようなエラーが出るため、
# 下記のdplyr::filter()のようにパッケージ名を明示的に指定するとうまく実行できる
# ---------------------------------------------
 new_wage_1 <- new_wage %>% 
   dplyr::filter(exper > 5)
# ---------------------------------------------
# ダミー変数'female'でデータサンプルの男女数を確認する
# table() は R における度数分布表を作成する関数、指定したベクトル（数値、文字列、因子など）の中で、# 各要素が何回出現するかを数え上げる時に使用
# ---------------------------------------------
table(new_wage_1$female)
## 
##   0   1 
## 209 181
# 変数`exper`を昇順（A→Z）で表示
sort(new_wage_1$exper)
##   [1]  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6  7  7  7  7  7  7  7
##  [26]  7  7  7  7  7  7  7  7  7  7  7  7  8  8  8  8  8  8  8  8  8  8  8  8  8
##  [51]  8  8  9  9  9  9  9  9  9  9  9  9  9  9  9  9  9  9  9  9 10 10 10 10 10
##  [76] 10 10 10 10 10 10 10 10 10 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11
## [101] 11 11 11 11 11 12 12 12 12 12 12 13 13 13 13 13 13 13 13 13 13 13 13 13 13
## [126] 13 13 14 14 14 14 14 14 14 14 14 14 14 14 14 14 14 14 14 14 14 14 15 15 15
## [151] 15 15 15 15 15 15 15 15 15 15 15 15 16 16 16 16 16 16 16 16 17 17 17 17 17
## [176] 17 17 17 17 18 18 18 18 18 18 18 18 18 19 19 19 19 19 19 19 19 19 19 19 20
## [201] 20 20 20 20 20 21 21 21 21 21 21 21 21 22 22 22 22 22 22 22 22 22 22 22 22
## [226] 23 23 23 23 23 23 23 23 23 23 23 23 24 24 24 24 24 24 25 25 25 25 25 25 26
## [251] 26 26 26 26 26 26 26 26 26 26 27 27 27 28 28 28 28 28 28 29 29 29 29 29 29
## [276] 29 29 29 29 30 30 30 31 31 31 31 31 31 31 31 31 31 32 33 33 33 33 33 33 33
## [301] 33 34 34 34 34 34 34 34 34 35 35 35 35 35 35 35 35 35 36 36 36 36 36 36 36
## [326] 36 36 36 37 37 37 37 37 37 38 38 38 38 38 39 39 39 39 39 39 40 40 40 40 40
## [351] 41 41 41 41 41 42 42 42 42 42 43 43 43 43 44 44 44 44 45 45 45 45 45 45 46
## [376] 47 47 47 47 47 47 48 48 48 49 49 49 49 50 51
# 降順でソート
sort(new_wage_1$exper, decreasing = TRUE)
##   [1] 51 50 49 49 49 49 48 48 48 47 47 47 47 47 47 46 45 45 45 45 45 45 44 44 44
##  [26] 44 43 43 43 43 42 42 42 42 42 41 41 41 41 41 40 40 40 40 40 39 39 39 39 39
##  [51] 39 38 38 38 38 38 37 37 37 37 37 37 36 36 36 36 36 36 36 36 36 36 35 35 35
##  [76] 35 35 35 35 35 35 34 34 34 34 34 34 34 34 33 33 33 33 33 33 33 33 32 31 31
## [101] 31 31 31 31 31 31 31 31 30 30 30 29 29 29 29 29 29 29 29 29 29 28 28 28 28
## [126] 28 28 27 27 27 26 26 26 26 26 26 26 26 26 26 26 25 25 25 25 25 25 24 24 24
## [151] 24 24 24 23 23 23 23 23 23 23 23 23 23 23 23 22 22 22 22 22 22 22 22 22 22
## [176] 22 22 21 21 21 21 21 21 21 21 20 20 20 20 20 20 19 19 19 19 19 19 19 19 19
## [201] 19 19 18 18 18 18 18 18 18 18 18 17 17 17 17 17 17 17 17 17 16 16 16 16 16
## [226] 16 16 16 15 15 15 15 15 15 15 15 15 15 15 15 15 15 15 14 14 14 14 14 14 14
## [251] 14 14 14 14 14 14 14 14 14 14 14 14 14 13 13 13 13 13 13 13 13 13 13 13 13
## [276] 13 13 13 13 12 12 12 12 12 12 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11
## [301] 11 11 11 11 11 11 10 10 10 10 10 10 10 10 10 10 10 10 10 10  9  9  9  9  9
## [326]  9  9  9  9  9  9  9  9  9  9  9  9  9  8  8  8  8  8  8  8  8  8  8  8  8
## [351]  8  8  8  7  7  7  7  7  7  7  7  7  7  7  7  7  7  7  7  7  7  7  6  6  6
## [376]  6  6  6  6  6  6  6  6  6  6  6  6  6  6  6
# データを`exper`の降順に並び替えし、`new_wage_1`に上書き
# new_wage_1 %>% arrange(-exper) → new_wage_1 のデータを exper列の値で並べ替える。-exper としているので、降順に並べ替えられる。

new_wage_1 <- new_wage_1 %>% 
  arrange(-exper)

データpopulation

  • 新しい変数を作成する
# ---------------------------------------------
# 次に外部から読み込んだ`population`でやってみる
# ---------------------------------------------
# データセット`population`にある変数名を１列で表示
ls(population)
##  [1] "（再掲）\r\n3世代世帯（世帯）"                             
##  [2] "（再掲）\r\n夫65歳以上，妻60歳以上の夫婦のみの世帯（世帯）"
##  [3] "15～64歳（％）...19"                                       
##  [4] "15～64歳（％）...25"                                       
##  [5] "15～64歳（％）...31"                                       
##  [6] "15～64歳（人）...16"                                       
##  [7] "15～64歳（人）...22"                                       
##  [8] "15～64歳（人）...28"                                       
##  [9] "15歳未満（％）...18"                                       
## [10] "15歳未満（％）...24"                                       
## [11] "15歳未満（％）...30"                                       
## [12] "15歳未満（人）...15"                                       
## [13] "15歳未満（人）...21"                                       
## [14] "15歳未満（人）...27"                                       
## [15] "2015年（平成27年）の世帯数（組替）（世帯）"                
## [16] "2015年（平成27年）の人口（組替）（人）"                    
## [17] "5年間の人口増減数（人）"                                   
## [18] "5年間の人口増減率（％）"                                   
## [19] "65歳以上（％）...20"                                       
## [20] "65歳以上（％）...26"                                       
## [21] "65歳以上（％）...32"                                       
## [22] "65歳以上（人）...17"                                       
## [23] "65歳以上（人）...23"                                       
## [24] "65歳以上（人）...29"                                       
## [25] "うち\r\n65歳以上の単独世帯（世帯）"                        
## [26] "うち\r\n単独世帯（世帯）"                                  
## [27] "うち\r\n核家族世帯（世帯）"                                
## [28] "一般世帯（世帯）...37"                                     
## [29] "一般世帯（世帯）...40"                                     
## [30] "人口密度（人/km2）"                                        
## [31] "人口性比"                                                  
## [32] "外国人（人）"                                              
## [33] "夫婦と子供\r\nから成る世帯（世帯）"                        
## [34] "夫婦のみの世帯（世帯）"                                    
## [35] "女"                                                        
## [36] "女親と子供\r\nから成る世帯（世帯）"                        
## [37] "市などの別\r\n_地域識別コード"                             
## [38] "平均年齢（歳）"                                            
## [39] "年齢中位数（歳）"                                          
## [40] "施設等の世帯（世帯）"                                      
## [41] "日本人（人）"                                              
## [42] "男"                                                        
## [43] "男親と子供\r\nから成る世帯（世帯）"                        
## [44] "総世帯（世帯）"                                            
## [45] "総数"                                                      
## [46] "都道府県・市区町村名"                                      
## [47] "都道府県・市区町村名_英語"                                 
## [48] "都道府県名"                                                
## [49] "面積(参考)"
# ---------------------------------------------
#  <chr>（文字列型）変数を<bdl>（数値型）変数に変換し、新しい変数を作成、新しいデータセット`population_1`にああああ代入
#　`gsub()`は「文字列置換関数」
#  `[^0-9]`は正規表現で、「数字以外のすべての文字」という意味
#  ""は置換後の文字列、ここでは「空文字」に置き換え、数字以外の文字を削除という意味。例: "12,345人" → "12345"
#　`as.numeric(gsub())`はgsub()によって作られた「数字のみの文字列」を、数値型(numeric)に変換

# ---------------------------------------------
#　mutate: `dplyr`パッケージに含まれている関数、新しい変数の追加、既存変数の加工などに使う
#　変数`総数`の文字列から数字だけ取り出して、数値型に変換し、新しい変数`totalpop`に指定
# ---------------------------------------------
population_1 <- population %>%
    select(1:11) %>%
  mutate(
    totalpop = as.numeric(gsub("[^0-9]", "", `総数`)),
    male = as.numeric(gsub("[^0-9]", "", `男`)),
    area = as.numeric(gsub("[^0-9]", "", `面積(参考)`))
  )

# 新しい変数を複数作成
population_1 <- population_1 %>%
  mutate(
    pop_density = totalpop / area,
    male_share = male / totalpop
  )

# IDをつける
population_1 <- population_1 %>%
mutate(ID = row_number())

glimpse(population_1)
## Rows: 1,965
## Columns: 17
## $ 都道府県名                               <chr> "00_全国", "01_北海道", "01_北海道", …
## $ `都道府県・市区町村名`                   <chr> "00000_全国", "01000_北海道", "011…
## $ `都道府県・市区町村名_英語`              <chr> "Japan", "Hokkaido", "Sapporo…
## $ `市などの別\r\n_地域識別コード`          <chr> "a", "a", "1", "0", "0", "0",…
## $ 総数                                     <chr> "126146099", "5224614", "1973…
## $ 男                                       <chr> "61349581", "2465088", "91868…
## $ 女                                       <chr> "64796518", "2759526", "10547…
## $ `2015年（平成27年）の人口（組替）（人）` <chr> "127094745", "5381733", "1952…
## $ `5年間の人口増減数（人）`                <dbl> -948646, -157119, 21039, 1105…
## $ `5年間の人口増減率（％）`                <chr> "-0.74641000000000002", "-2.9…
## $ `面積(参考)`                             <chr> "377976.41", "83424.44", "112…
## $ totalpop                                 <dbl> 126146099, 5224614, 1973395, …
## $ male                                     <dbl> 61349581, 2465088, 918682, 11…
## $ area                                     <dbl> 3.779764e+07, 8.342444e+06, 1…
## $ pop_density                              <dbl> 3.337407e+00, 6.262690e-01, 1…
## $ male_share                               <dbl> 0.4863375, 0.4718220, 0.46553…
## $ ID                                       <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10…
# 変数を削除する。例："男"という変数を削除
population_1 <- population_1 %>%
  select(-c(男,女,総数))
  • 加工後のデータをCSV、XLSXファイルとして書き出し、現在の作業ディレクトリに保存される
# ---------------------------------------------
# データセット`population_1`をCSVファイルとして書き出し、現在の作業ディレクトリを確認してください。
# ---------------------------------------------
write_csv(population_1, "population_1.csv")

# ---------------------------------------------
# データセット`population_1`をXLSXファイルとして書き出、現在の作業ディレクトリを確認してください。
# ---------------------------------------------
write_xlsx(population_1, "population_1.xlsx")

2.9 基本統計量の算出と統計表の作成

複数の変数の最小値、最大値、中位数、平均値などを出力、見やすい表に整理

データwage1

# データセット`wage1`を使って統計量を確認する
summary(wage1)
##       wage             educ           exper           tenure      
##  Min.   : 0.530   Min.   : 0.00   Min.   : 1.00   Min.   : 0.000  
##  1st Qu.: 3.330   1st Qu.:12.00   1st Qu.: 5.00   1st Qu.: 0.000  
##  Median : 4.650   Median :12.00   Median :13.50   Median : 2.000  
##  Mean   : 5.896   Mean   :12.56   Mean   :17.02   Mean   : 5.105  
##  3rd Qu.: 6.880   3rd Qu.:14.00   3rd Qu.:26.00   3rd Qu.: 7.000  
##  Max.   :24.980   Max.   :18.00   Max.   :51.00   Max.   :44.000  
##     nonwhite          female          married           numdep     
##  Min.   :0.0000   Min.   :0.0000   Min.   :0.0000   Min.   :0.000  
##  1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:0.000  
##  Median :0.0000   Median :0.0000   Median :1.0000   Median :1.000  
##  Mean   :0.1027   Mean   :0.4791   Mean   :0.6084   Mean   :1.044  
##  3rd Qu.:0.0000   3rd Qu.:1.0000   3rd Qu.:1.0000   3rd Qu.:2.000  
##  Max.   :1.0000   Max.   :1.0000   Max.   :1.0000   Max.   :6.000  
##       smsa           northcen         south             west       
##  Min.   :0.0000   Min.   :0.000   Min.   :0.0000   Min.   :0.0000  
##  1st Qu.:0.0000   1st Qu.:0.000   1st Qu.:0.0000   1st Qu.:0.0000  
##  Median :1.0000   Median :0.000   Median :0.0000   Median :0.0000  
##  Mean   :0.7224   Mean   :0.251   Mean   :0.3555   Mean   :0.1692  
##  3rd Qu.:1.0000   3rd Qu.:0.750   3rd Qu.:1.0000   3rd Qu.:0.0000  
##  Max.   :1.0000   Max.   :1.000   Max.   :1.0000   Max.   :1.0000  
##     construc          ndurman          trcommpu           trade       
##  Min.   :0.00000   Min.   :0.0000   Min.   :0.00000   Min.   :0.0000  
##  1st Qu.:0.00000   1st Qu.:0.0000   1st Qu.:0.00000   1st Qu.:0.0000  
##  Median :0.00000   Median :0.0000   Median :0.00000   Median :0.0000  
##  Mean   :0.04563   Mean   :0.1141   Mean   :0.04373   Mean   :0.2871  
##  3rd Qu.:0.00000   3rd Qu.:0.0000   3rd Qu.:0.00000   3rd Qu.:1.0000  
##  Max.   :1.00000   Max.   :1.0000   Max.   :1.00000   Max.   :1.0000  
##     services         profserv         profocc          clerocc      
##  Min.   :0.0000   Min.   :0.0000   Min.   :0.0000   Min.   :0.0000  
##  1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:0.0000  
##  Median :0.0000   Median :0.0000   Median :0.0000   Median :0.0000  
##  Mean   :0.1008   Mean   :0.2586   Mean   :0.3669   Mean   :0.1673  
##  3rd Qu.:0.0000   3rd Qu.:1.0000   3rd Qu.:1.0000   3rd Qu.:0.0000  
##  Max.   :1.0000   Max.   :1.0000   Max.   :1.0000   Max.   :1.0000  
##     servocc           lwage            expersq          tenursq       
##  Min.   :0.0000   Min.   :-0.6349   Min.   :   1.0   Min.   :   0.00  
##  1st Qu.:0.0000   1st Qu.: 1.2030   1st Qu.:  25.0   1st Qu.:   0.00  
##  Median :0.0000   Median : 1.5369   Median : 182.5   Median :   4.00  
##  Mean   :0.1407   Mean   : 1.6233   Mean   : 473.4   Mean   :  78.15  
##  3rd Qu.:0.0000   3rd Qu.: 1.9286   3rd Qu.: 676.0   3rd Qu.:  49.00  
##  Max.   :1.0000   Max.   : 3.2181   Max.   :2601.0   Max.   :1936.00
# 一部の変数と基本統計量のみ出力する
# 行の部分 → 空白（= すべての行を意味する）
# 列の部分 → c('wage','educ','exper','tenure','nonwhite')（= その列だけ取る）

summary(wage1[, c('wage','educ', 'exper', 'tenure', 'nonwhite')], fast=TRUE)　
##       wage             educ           exper           tenure      
##  Min.   : 0.530   Min.   : 0.00   Min.   : 1.00   Min.   : 0.000  
##  1st Qu.: 3.330   1st Qu.:12.00   1st Qu.: 5.00   1st Qu.: 0.000  
##  Median : 4.650   Median :12.00   Median :13.50   Median : 2.000  
##  Mean   : 5.896   Mean   :12.56   Mean   :17.02   Mean   : 5.105  
##  3rd Qu.: 6.880   3rd Qu.:14.00   3rd Qu.:26.00   3rd Qu.: 7.000  
##  Max.   :24.980   Max.   :18.00   Max.   :51.00   Max.   :44.000  
##     nonwhite     
##  Min.   :0.0000  
##  1st Qu.:0.0000  
##  Median :0.0000  
##  Mean   :0.1027  
##  3rd Qu.:0.0000  
##  Max.   :1.0000
# 一部の変数を統計表の形式に整理する
describe(wage1[, c('wage','educ', 'exper', 'tenure', 'nonwhite')], fast=TRUE)　
##          vars   n  mean    sd median  min   max range  skew kurtosis   se
## wage        1 526  5.90  3.69   4.65 0.53 24.98 24.45  2.00     4.94 0.16
## educ        2 526 12.56  2.77  12.00 0.00 18.00 18.00 -0.62     1.87 0.12
## exper       3 526 17.02 13.57  13.50 1.00 51.00 50.00  0.70    -0.65 0.59
## tenure      4 526  5.10  7.22   2.00 0.00 44.00 44.00  2.10     4.63 0.32
## nonwhite    5 526  0.10  0.30   0.00 0.00  1.00  1.00  2.61     4.83 0.01
# 統計表のCSVファイルとしての外部出力
describe(wage1[, c('wage','educ', 'exper', 'tenure', 'nonwhite')]) %>% 
  write_csv("summary.csv")

# 統計表のExcelファイルとしての外部出力
describe(wage1[, c('wage','educ', 'exper', 'tenure', 'nonwhite')]) %>% 
  write_xlsx("summary.xlsx")

2.10 練習課題

練習：データpopulation_1で統計記述表を作成してみましょう。

3 Rを使った練習：可視化１

3.1 関連パッケージ

library(wooldridge)
library(tidyverse)
library(readxl)
library(writexl)
library(psych)
library(ggplot2)

3.2 データの可視化：分布図

1つの変数の値の分布を可視化する

  • ヒストグラムを作成
# ---------------------------------------------
# データ`wage1`の変数`wage` を使って、ヒストグラム（度数分布図）を作る
# geom_histogram() の設定はすべてデフォルト。「geom」は geometric object の略、グラフの形を指定する関数。
# 他のgeom_*()の例→ geom_point()：散布図；　geom_line()：折れ線グラフ；　geom_bar()：棒グラフ；　geom_boxplot()：箱ひげ図
# ---------------------------------------------
library(ggplot2)
ggplot(data=wage1, aes(x=wage)) +
  geom_histogram()
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

# 変数`wage`の要約統計量を表示
summary(wage1$wage)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   0.530   3.330   4.650   5.896   6.880  24.980
# ---------------------------------------------
# ビン幅(binwidth)を1に指定してヒストグラムを作る
# ---------------------------------------------
ggplot(wage1, aes(wage)) +
  geom_histogram(binwidth = 1)

# ---------------------------------------------
# ビン幅を1に、棒グラフの内部(fill)を「スチールブルー」に塗り、枠線(color)を白に設定
# タイトル(title)を「Histogram」、x軸ラベルを「Wage」、y軸ラベルを「Count」に設定
# ---------------------------------------------
ggplot(wage1, aes(wage)) +
  geom_histogram(binwidth = 1, fill = "steelblue", color = "white") +
  labs(title = "Histogram", x = "Wage", y = "Count")

  • Log変換してヒストグラムを作成
# ---------------------------------------------
# `wage`変数の自然対数(log)を計算して、`wage1`に新たに`logwage`という変数列を追加
# ---------------------------------------------
wage1 <- wage1 %>%
  mutate(
    logwage = log(wage)
  )

# ---------------------------------------------
# 変数`logwage`の要約統計量を表示
# ---------------------------------------------
summary(wage1$logwage)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -0.6349  1.2030  1.5369  1.6233  1.9286  3.2181
# ---------------------------------------------
# `logwage`を使ってヒストグラムを作る
# ---------------------------------------------
ggplot(wage1, aes(logwage)) +
  geom_histogram(binwidth = 0.1, fill = "steelblue", color = "white", alpha = 0.7) +
  labs(title = "Histogram", x = "Log(Wage)", y = "Count")

3.3 練習課題

  • ステップ１：パッケージ｀Wooldridge｀に含まれるデータセット｀wage1｀を読み込む
  • ステップ２：データ｀wage1｀に含まれる変数｀exper｀のヒストグラムをggplot2で描く
  • ステップ３：変数｀exper｀をログ変換して新しい変数logexperを作成し、logexperのヒストグラムをggplot2で描く

3.4 データの可視化：散布図

2つの変数の関係性（相関）を可視化する

  • 相関係数
# ---------------------------------------------
# データ`wage1`の中の変数`educ`と`wage`の相関係数を計算
# ---------------------------------------------
cor.test(wage1$educ,wage1$wage)
## 
##  Pearson's product-moment correlation
## 
## data:  wage1$educ and wage1$wage
## t = 10.167, df = 524, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.3319280 0.4749166
## sample estimates:
##       cor 
## 0.4059033
# ---------------------------------------------
# データ`wage1`中の変数`educ`と`logwage`の相関係数を計算
# ---------------------------------------------
cor.test(wage1$educ,wage1$logwage)
## 
##  Pearson's product-moment correlation
## 
## data:  wage1$educ and wage1$logwage
## t = 10.935, df = 524, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.3587807 0.4981874
## sample estimates:
##       cor 
## 0.4310528
  • educとwage相関関係を示す散布図を作成
# ---------------------------------------------
# データ`wage1`の中の変数`educ`をx軸、`wage`をy軸として、散布図を作る
# ---------------------------------------------
library(ggplot2)
ggplot(wage1, aes(educ, wage)) + 
geom_point() 

# ---------------------------------------------
# グラフを`corrplot.pdf`という名前で、横5インチ × 縦3インチのサイズで現在の作業ディレクトリに保存
# ---------------------------------------------
ggsave("corrplot.pdf", width = 5, height = 3) 


# ---------------------------------------------
# `educ`と`wage`の散布図に、線形回帰の回帰直線(method = "lm")を赤(col = "red")で重ねて表示
# ---------------------------------------------
ggplot(wage1, aes(x=educ, y=wage)) + 
  geom_point() + 
  geom_smooth(method = "lm", color = "red")
## `geom_smooth()` using formula = 'y ~ x'

# ---------------------------------------------
# グラフを`corrplot_1.pdf`という名前で、横5インチ × 縦3インチのサイズで現在の作業ディレクトリに保存
# ---------------------------------------------
ggsave("corrplot_1.pdf", width = 5, height = 3) 
## `geom_smooth()` using formula = 'y ~ x'
  • educとlogwage相関関係を示す散布図を作成
# ---------------------------------------------
# データ`wage1`の中の変数`educ`をx軸、`logwage`をy軸として、散布図を作る
# ---------------------------------------------
library(ggplot2)
ggplot(wage1, aes(educ, logwage)) + 
geom_point() 

# ---------------------------------------------
# `educ`と`logwage`の散布図に、線形回帰の回帰直線を赤で重ねて表示
#  geom_smooth(method = "lm", color = "red")で回帰直線を追加; method = "lm" → 線形回帰（linear model）をフィット; color = "red" → 線の色は赤
# ---------------------------------------------
ggplot(wage1, aes(x=educ, y=logwage)) + 
  geom_point() +
  geom_smooth(method = "lm", color = "red")
## `geom_smooth()` using formula = 'y ~ x'

3.5 練習課題

  • ステップ１：パッケージ｀Wooldridge｀に含まれるデータセット｀wage1｀を読み込む
  • ステップ２：データ｀wage1｀に含まれる変数｀exper｀と｀wage｀の相関係数を計算する
  • ステップ３：変数｀exper｀と｀wage｀の散布図をggplot2で描く

4 Rを使った練習：可視化２

4.1 関連パッケージ

library(wooldridge)
library(tidyverse)
library(readxl)
library(writexl)
library(psych)
library(ggplot2)
# インストールされていない場合はまずインストールしてください
library(stargazer)
## 
## Please cite as:
##  Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.
##  R package version 5.2.3. https://CRAN.R-project.org/package=stargazer
library(modelsummary)
## 
## Attaching package: 'modelsummary'
## The following object is masked from 'package:psych':
## 
##     SD

4.2 単回帰分析

# データを読み込む
# ---------------------------------------------
# csvファイルを読み込む場合
# library(readr) → readrパッケージのread_csv関数を使って、CSVファイルを読み込む
# rent <- read_csv("rent_0521.csv") → csvファイルを'rent'という名前のデータオブジェクトとして保存
# ---------------------------------------------

library(wooldridge)
data(wage1)

# 回帰: wageを目的変数とする回帰
# lm()：線形モデル(linear model)を作成する関数
result <- lm(wage ~ educ, data = wage1)

# resultの詳細な統計的要約を表示
summary(result)
## 
## Call:
## lm(formula = wage ~ educ, data = wage1)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.3396 -2.1501 -0.9674  1.1921 16.6085 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.90485    0.68497  -1.321    0.187    
## educ         0.54136    0.05325  10.167   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.378 on 524 degrees of freedom
## Multiple R-squared:  0.1648, Adjusted R-squared:  0.1632 
## F-statistic: 103.4 on 1 and 524 DF,  p-value: < 2.2e-16

4.3 重回帰分析

# データを読み込む
library(wooldridge)
data(wage1)

# 回帰1: wageを目的変数とする回帰
result1 <- lm(wage ~ educ + exper, data = wage1)

# lwage（対数変換）列を追加
wage1$lwage <- log(wage1$wage)

# 回帰2: log(wage) を目的変数とする回帰
result2 <- lm(lwage ~ educ + exper, data = wage1)

4.4 分析結果の可視化

参考資料

stargazerで出力結果をテーブルにまとめて書き出す

# パッケージ読み込み
library(stargazer)

# ---------------------------------------------
# 記述統計をHTML出力
# ---------------------------------------------

data(wage1)

stargazer(wage1, # データの基本統計量をテーブル形式に出力
  out = 'summary_wage1.html',      # 出力ファイル名
  type = "html" )            # 出力形式を明示的にHTMLにする（エラー防止）


# ---------------------------------------------
# 生成されたsummary_wage1.html は、Rの作業ディレクトリに保存される
# その場所を確認するには以下を実行してください：
# ---------------------------------------------
getwd()  # 現在の作業ディレクトリを確認

# ---------------------------------------------
# 回帰結果をHTML出力
# ---------------------------------------------
stargazer(result1, result2,
  digits = 3,               # 小数点以下3桁表示
  digits.extra = 0,         # 追加の桁数なし
  align = TRUE,             # 表の列を整列
  keep.stat = c('n', 'adj.rsq', 'f'), # 必要な統計量のみ表示
  df = FALSE,               # 自由度を表示しない
  dep.var.labels = c("Wage", "LogWage"), # 各モデルの目的変数名
  title = 'OLS Results',    # タイトル
  out = 'result.html',      # 出力ファイル名
  type = "html" )            # 出力形式を明示的にHTMLにする（エラー防止）


# ---------------------------------------------
# 生成された result.html は、Rの作業ディレクトリに保存される
# その場所を確認するには以下を実行してください：
# ---------------------------------------------
getwd()  # 現在の作業ディレクトリを確認

modelsummaryで出力結果をテーブルにまとめて書き出す

# データを読み込む
data(wage1)

# 回帰1: wageを目的変数とする単回帰分析
result1 <- lm(wage ~ educ + exper, data = wage1)


# 回帰2: wageを目的変数とする重回帰分析
result2 <- lm(wage ~ educ + exper, data = wage1)

# lwage（対数変換）列を追加
wage1$lwage <- log(wage1$wage)

# 回帰3: log(wage) を目的変数とする回帰
result3 <- lm(lwage ~ educ + exper, data = wage1)

# ---------------------------------------------
# 3つの回帰モデル (result1, result2, result3) を`models`という名前付きリストとして格納
# modelsummary() 関数などで複数モデルを比較するための準備
# ---------------------------------------------
models <- list("model 1"=result1,
               "model 2"=result2,
               "model 3"=result3)
# ---------------------------------------------
# 回帰結果をテーブルで表示
# ---------------------------------------------
# NOTE: `statistic = ""`この引数には、表示したい統計量のタイプを指定できる
#  "std.error": 標準誤差（Standard Error）
#  "statistic": 検定統計量（t値など）
#  "p.value": p値
#  "conf.int": 信頼区間
# ---------------------------------------------
modelsummary(models, # モデルリストの名前
             fmt = 3, # 小数点以下3桁まで表示
             statistic = "std.error", # 各係数の下に標準誤差を表示
             stars = TRUE, # 星印を付ける
             conf_level = 0.95, # 95%信頼区間を計算
             gof_omit='RMSE|AIC|BIC|Log.Lik.|F')  # 表示不要の指標を非表示
model 1 model 2 model 3
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001
(Intercept) -3.391*** -3.391*** 0.217*
(0.767) (0.767) (0.109)
educ 0.644*** 0.644*** 0.098***
(0.054) (0.054) (0.008)
exper 0.070*** 0.070*** 0.010***
(0.011) (0.011) (0.002)
Num.Obs. 526 526 526
R2 0.225 0.225 0.249
R2 Adj. 0.222 0.222 0.246 
# 回帰結果をテーブルで表示（星印なしのバージョン）
modelsummary(models, 
             fmt = 3, 
             statistic = "statistic", # 各係数の下に標準誤差ではなくt値を表示
             conf_level = 0.95, # 95%信頼区間を計算
             gof_omit='RMSE|AIC|BIC|Log.Lik.|F')
model 1 model 2 model 3
(Intercept) -3.391 -3.391 0.217
(-4.423) (-4.423) (1.997)
educ 0.644 0.644 0.098
(11.974) (11.974) (12.848)
exper 0.070 0.070 0.010
(6.385) (6.385) (6.653)
Num.Obs. 526 526 526
R2 0.225 0.225 0.249
R2 Adj. 0.222 0.222 0.246 
# ---------------------------------------------
# 回帰結果をWORD形式で出力
# ---------------------------------------------
modelsummary(models, gof_omit='RMSE|AIC|BIC|Log.Lik.|F', "results.docx")

modelsummaryで出力結果をグラフに描く

# ---------------------------------------------
# 回帰係数の可視化（モデル3）
# ---------------------------------------------
# `modelplot()`はモデル3の係数と信頼区間を可視化
# `coef_omit = "Interc"`は切片を表示しない
# `conf_level = .95`は95%信頼区間を表示
# `color = 'black'`は係数点と信頼区間バーの色を指定
# `coef_map`は変数名をわかりやすいラベルに変更
# `geom_vline(xintercept=0)`は参考線で係数のゼロラインを表示
# ---------------------------------------------
modelplot(result3, 
          coef_omit = "Interc", 
          conf_level = .95, 
          color = 'black',
          coef_map = c(
            "educ" = "Education",
            "exper" = "Experience")) + 
  geom_vline(xintercept = 0, color = "maroon") +
  labs(x = "Coefficients (95% CI)", 
       title = "A.モデル３の分析結果") +
       theme_bw(base_family = jpfont)

4.5 練習課題

Step 1~6 に従って演習を行い、散布図の作成と回帰分析をしてみましょう

  • Step 1：パッケージ｀AER｀をインストールする
  • Step 2：パッケージ｀AER｀に含まれるデータセット｀MASchools｀を読み込む
  • Step 3：データ｀MASchools｀に含まれる変数｀income｀と｀score8｀の散布図をggplot2で作成する
  • Step 4：できる範囲で、軸ラベルの設定、タイトルの追加、テーマの調整なども挑戦してみてください
  • Step 5：｀score8｀を被説明変数、｀income｀を説明変数として、単回帰分析を行い、結果をHTMLやWORD形式に出力してみてください
  • Step 6：｀income｀に加え、他の変数（参考資料にあるデータの説明を確認してください）も説明変数として取り入れて重回帰分析を試してみてください

参考資料

5 Rを使った練習：可視化３

参考資料

5.1 関連パッケージ

library(wooldridge)
library(tidyverse)
library(readxl)
library(writexl)
library(psych)
library(ggplot2)

5.2 データの可視化：表示

散布図の結合・比較

data(wage1)

wage1 <- wage1 %>% 
  mutate(logwage=log(wage1$wage))
install.packages("patchwork")
# 複数の`ggplot2`グラフを簡単に並べて表示・結合するためのパッケージ
library(patchwork)

# 2つのグラフを作成、それぞれp1、p2として保存
p1 <- ggplot(wage1, aes(x = educ, y = wage)) +
  geom_point() +
  geom_smooth(method = "lm", color = "red") +
  labs(title = "A. Wage")

p2 <- ggplot(wage1, aes(x = educ, y = logwage)) +
  geom_point(shape=1) +
  geom_smooth(method = "lm", color = "blue") +
  labs(title = "B. Log(Wage)")
# 横に並べて表示
p1 + p2
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'

# 縦に並べて表示
p1 / p2 
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'

# ---------------------------------------------
# グラフを`corrplot_combine.pdf`という名前で、横10インチ × 縦5インチのサイズで現在の作業ディレクトリに保存
# ---------------------------------------------
ggsave("corrplot_combine.pdf", p1+p2, width = 10, height = 5) 
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'

5.3 データの可視化：配色

ggplotのカラーパレットテーマの試し

# GitHubからパッケージ`scico`をインストール
devtools::install_github("thomasp85/scico")
# 利用可能なパレットを確認
library(scico)
scico_palette_show()

# ---------------------------------------------
# 男女色別で散布図を作り、変数femaleで分類する（femaleが数値型: 0 or 1）
# `aes(color = 変数名)`でポイントの色(color)を変数に応じて指定
# scale_color_scico()は連続値変数（数値型）の配色に使う
# ---------------------------------------------
ggplot(data = wage1, aes(x = educ, y = wage, color = female)) +
  geom_point() +
  scale_color_scico(palette = "berlin")  # カラーパレットを指定

# ---------------------------------------------
# グループごとの回帰線付き（femaleをカテゴリ変数として扱う）
# `scale_color_scico_d()`は離散変数（カテゴリ変数）の配色に使う
# `scale_color_scico_d()`でカテゴリ用の色を指定
# `labs()` を使って軸ラベルや凡例タイトルを変更する
# `theme_bw()`はシンプルな背景を使う
# `base_family = jpfont`で設定したフォントを使うように指定、そうしないと日本語が出てこない
# ---------------------------------------------
ggplot(data = wage1, aes(x = educ, y = wage, color = as.factor(female))) +
  geom_point() +
  geom_smooth(method = "lm", se = FALSE) + 
  scale_color_scico_d(palette = "berlin") +　# _d は「離散」色を使う指定
  labs( color = "性別", # 凡例タイトルを変更
        x= "教育",
        y= "賃金") +
  theme_bw(base_family = jpfont)　
## `geom_smooth()` using formula = 'y ~ x'

# ---------------------------------------------
# `theme_class()`や`theme_minimal()`テーマも指定してみてください
# ---------------------------------------------
# `theme_bw()`などのテーマを指定したくない場合、下記のように設定すれば日本語が表示できる
# `theme_bw(base_family = jpfont)`の代わりに使う
# ---------------------------------------------
theme(text = element_text(family = jpfont))　

5.4 各種グラフ

参考資料

library(wooldridge)
library(tidyverse)
library(readxl)
library(writexl)
library(psych)
library(ggplot2)
library(scico)　# カラーパレット拡張パッケージ
library(patchwork)　# 複数のggplotグラフを並べて表示するためのパッケージ

棒グラフの作成

# ---------------------------------------------
# `wooldridge`パッケージの付属時系列データ`intdef`を使う
# ---------------------------------------------
data("intdef")

# ---------------------------------------------
# `geom_bar()`関数で棒グラフを作成する
# `stat = "identity"`はyの値を棒の高さとして使用することを指定する
# `coord_flip()`関数で軸を入れ替えて横棒グラフにする
# ---------------------------------------------

# 縦棒グラフを作成し、p1に代入
p1 <- ggplot(intdef, aes(x = year, y = inf)) +
  geom_bar(stat = "identity", color="white", linewidth=0.1, fill="blue") +
  labs(y="Inflation rate")

# 横向きの棒グラフを作成し、p2に代入
p2 <- ggplot(intdef, aes(x = year, y = inf)) +
  geom_bar(stat = "identity") +
  coord_flip() +
  labs(y="Inflation rate")
# 並べて表示
p1+p2

折れ線グラフの作成

# ---------------------------------------------
# `geom_line()`関数で線グラフを作成する
# `geom_line()`関数で折れ線グラフの下を塗りつぶしたグラフになる
# 色(color)や線種(linetype)を指定することで、見やすくなる
# ---------------------------------------------
data("intdef")

# 線グラフを作成し、p1に代入
p1 <- ggplot(intdef, aes(x = year, y = inf)) +
  geom_line(color = "black", linewidth = 0.5, linetype = 1) +
  labs(y="Inflation rate")

# 折れ線＋塗りつぶしのグラフを作成し、p2に代入
p2 <- ggplot(intdef, aes(x = year, y = inf)) +
  geom_line(color = "black", linewidth = 0.5, linetype = 1) +
  geom_area(fill = "gray") +
  labs(y="Inflation rate") 
p1 + p2

# ---------------------------------------------
# 複数系列の折れ線グラフを作成する
# `geom_line()`を複数回使い、それぞれ異なる変数をy軸に指定して描く
# `aes(color = "ラベル名")`を使うと、凡例に表示する変数の名前を変更
# ---------------------------------------------
ggplot(intdef, aes(x = year)) +
  geom_line(aes(y = inf, color = "Inflation"), linewidth = 0.8, linetype = 1) +
  geom_line(aes(y = i3, color = "Treasury bill"), linewidth = 0.8, linetype = 2) +
  scale_color_scico_d(palette = "berlin") +
  labs(color = "Type", # 凡例タイトルを変更
       x = "Year",
       y = "Rate") +
  theme_minimal()

パネルデータを使う時折れ線グラフの作成

  • パネルデータは、観測単位(州など)×時系列(年)という構造になっているため、すべての観測値を一度にプロットすると系列が多すぎて図が読みにくくなる。そのため、「年ごとの平均値」をとって代表的な傾向を把握することがよく行われる
# ---------------------------------------------
# `AER`パッケージの付属パネルデータ`Guns`を使う
# ---------------------------------------------
data("Guns", package = "AER")

# ---------------------------------------------
# `group_by()`関数で年ごとにデータをグループ化し、violentとrobberyの平均値を計算
# `na.rm = TRUE`は、欠損値(NA)を除いて年ごと平均を計算
# ---------------------------------------------
guns <- Guns %>%
  group_by(year) %>%　
  mutate( year=as.numeric(as.character(year)),
          violent_mean = mean(violent, na.rm = TRUE),
          rob_mean = mean(robbery, na.rm = TRUE))

# ---------------------------------------------
# `aes(color = "ラベル名")`を使うと、凡例に表示する変数の名前を変更
# `scale_color_manual()`でggplot2の凡例で使われる色(color)を手動で指定できる
# `values=`でどのラベルにどの色を使うかを指定
# ---------------------------------------------
ggplot(guns, aes(x = year)) +
  geom_line(aes(y = violent_mean, color = "Violent"), linewidth = 1, linetype = 1) +
  geom_line(aes(y = rob_mean, color = "Robbery"), linewidth = 1, linetype = 2) +
  scale_color_manual(values = c("Violent" = "red", "Robbery" = "blue")) +
  labs(color = "Crime type", # 凡例タイトルを変更
       x = "Year",
       y = "Crime rate")  +
  theme_minimal()

マップの作成

# 関連パッケージをインストールする
install.packages(c("sf", "spData"))
# 関連パッケージを読み込む
library(sf) #空間データ(地図など)を扱うパッケージ
## Linking to GEOS 3.13.0, GDAL 3.8.5, PROJ 9.5.1; sf_use_s2() is TRUE
library(spData) #空間データのサンプルデータ集
## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`
# ---------------------------------------------
# `spData`パッケージのデータ`world`を使う
# ---------------------------------------------
data(world)

# ---------------------------------------------
# `aes(fill = )`で塗りつぶしの色(fill)を変数に応じて指定
# `geom_sf()`関数で地図の国境線・形状を描く
# `coord_sf()`で地図の座標系を指定（歪みを調整）
# `scale_fill_scico()`関数で塗りつぶしの色スケール(カラーパレット)を指定
# `labs(fill = "")`で凡例のタイトルを変更
# ---------------------------------------------
map <- ggplot(world, aes(fill = lifeExp)) + 
  geom_sf(data = world) + 
  coord_sf(crs = "+proj=robin") + 
  scale_fill_scico(palette = "berlin", direction = -1) +
  labs(fill = "Life expectancy") + # 凡例タイトルを変更
  theme_void()

map