Rメモ
田中 鮎夢
2024-02-11
This version: 2025-04-19
1 Rのグラフでの日本語の表示
- 日本語を表示するためには、フォントを指定する必要がある。
- しかし、Rマークダウンをknitしたとき、HTMLでは日本語が表示されても、PDFに変換するときにエラーが発生することがある。
- その場合、以下の対処策を実施。
#install.packages("showtext")
library(showtext)
# 日本語フォントの追加と有効化
font_add(family = "MSMincho", regular = "MSMincho.ttf")
showtext_auto()
# プロット時にフォント指定
par(family = "MSMincho")2 変数ラベル
- 変数ラベルを付ける方法
attr関数を使う。
library(wooldridge)
data("wage1")
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 64wage1データセットのwage変数にラベルを付ける。
# ラベルを付与する
attr(wage1$wage, "label") <- "賃金"
# ラベルを表示
attr(wage1$wage, "label")## [1] "賃金"- ラベルを消す
attr(wage1$wage, "label") <- NULL3 Rでデータを接続する方法
rbind関数を使うと、データを縦に結合できる。cbind関数を使うと、データを横に結合できる。merge関数を使うと、key変数を用いてデータを結合できる。- 同じく
left_join関数を使うと、key変数を用いてデータを結合できる。
library(wooldridge)
data(murder)
library(datasets)
data(state)
# cbindを使って、アメリカの州の略号と名前を接続
us_state <- cbind(state.abb, state.name)
# us_state をデータフレームにしておく。
us_state <- as.data.frame(us_state)
# mergeを使ってmurderデータとアメリカの州の略号と名前を接続
muder_merged <- merge(murder, us_state,
by.x = "state", by.y = "state.abb")
library(dplyr)
# left_join関数を用いて、murderとus_stateを接続する。
muder_left_join <- left_join(murder, us_state,
by = c("state" = "state.abb"))4 Rでデータをcollapseする方法
4.1 方法1
aggregate関数
aggregate関数を使えば、集計量を計算できる。
state.x77 <- state.x77
state.x77.collapsed <- aggregate(state.x77, list(Region = state.region), mean)
warpbreaks <- warpbreaks
warpbreaks_collapsed <- aggregate(breaks ~ wool + tension, data = warpbreaks, mean)4.2 方法2
dplyrパッケージ
dplyrパッケージを使うと、データの集計が簡単にできる。group_by関数でグループ化し、summarise関数で集計する。- データが大きい場合は、
dplyrパッケージでは処理が不可能になることがある。その場合は、data.tableパッケージを使うと良い。
library(wooldridge)
data(wage1)
library(dplyr)
# 男女・人種別賃金の記述統計を計算
collapse_dplyr <- wage1 %>%
group_by(female, nonwhite) %>%
summarise(mean_wage = mean(wage, na.rm = TRUE),
sum_wage = sum(wage, na.rm = TRUE),
first_wage = first(wage),
n_worker = n())4.3 方法3
data.tableパッケージ
data.tableパッケージは、大規模データの読み込みや処理に適している。data.tableパッケージを使うと、データの集計がdplyrパッケージより高速にできる。dplyrパッケージでは、group_byとsummarise関数を使っていたが、data.tableパッケージでは、by引数を使って集計する。- ()の前に.をつけることで集計することができる。
data.tableパッケージの使い方は、以下のリンクを参照。
library(wooldridge)
data(wage1)
library(data.table)
# データ・テーブルに変換する
wage_df <- setDT(wage1)
# `data.table`を使って男女・人種別賃金の記述統計を計算
library(data.table)
# Group by female and nonwhoite and calculate mean and sum of wages
collapse_data.table <- wage_df[, .(mean_wage = mean(wage, na.rm = TRUE),
sum_wage = sum(wage, na.rm = TRUE),
n_worker = .N,
first_wage = first(wage)),
by = .(female, nonwhite)]4.4 方法4
doBy
library(doBy)
hsb2 <- read.table("https://stats.idre.ucla.edu/wp-content/uploads/2016/02/hsb2-1.csv", header=T, sep=",")
collapse_doBy <- summaryBy(socst + math ~ prog + ses + female, FUN=c(mean,sd), data=hsb2)
collapse_doBy## prog ses female socst.mean math.mean socst.sd math.sd
## 1 1 1 0 47.57143 46.71429 6.502747 8.118175
## 2 1 1 1 49.00000 48.22222 5.894913 7.546154
## 3 1 2 0 50.50000 51.10000 8.959787 9.267026
## 4 1 2 1 52.10000 51.10000 9.938142 6.026792
## 5 1 3 0 57.25000 54.00000 17.500000 4.320494
## 6 1 3 1 49.60000 50.40000 10.714476 7.700649
## 7 2 1 0 46.75000 50.00000 11.926860 6.683313
## 8 2 1 1 55.13333 55.00000 10.098562 10.281745
## 9 2 2 0 55.54545 58.13636 8.985318 10.086745
## 10 2 2 1 55.40909 54.50000 8.511388 7.255540
## 11 2 3 0 59.33333 57.42857 7.799573 7.743200
## 12 2 3 1 59.61905 59.42857 9.058014 8.152125
## 13 3 1 0 32.50000 46.75000 4.725816 5.251984
## 14 3 1 1 38.25000 42.25000 8.189715 3.918819
## 15 3 2 0 44.00000 48.20000 11.464230 9.540889
## 16 3 2 1 51.00000 46.06250 8.755950 7.540723
## 17 3 3 0 50.25000 42.25000 6.946222 3.947573
## 18 3 3 1 46.00000 55.66667 10.000000 10.503968注)
FUN = meanを省略すると、デフォルトで平均値が計算される。
library(wooldridge)
data(wage1)
library(doBy)
# 男女別賃金の平均値を計算
collapse_FUNmean <- summaryBy(wage ~ female, data = wage1, FUN = mean)
# `FUN = mean`は省略可能。 (デフォルトは平均値)
collapse_noFUNmean <- summaryBy(wage ~ female, data = wage1)5 欠損値の除去
- 変数のない行を削除
dat2 <- subset(dat, !(is.na(dat$変数)))- NAがある行はすべて削除
dat2 <- na.omit(dat)6
$と[[]]の使い方
$と[[]]のいずれを使っても、データの中の特定の変数を指定できる。- どちらもBase Rのオペレーター。
- 既存変数の対数を作成する例。
[[]]を使う場合、変数名を引用符で囲む必要があることに注意。例）“wage”$を使う場合、変数名を引用符で囲む必要はない。例）wage
# AERパッケージのデータを使う
library(AER)
data("CPS1985")
# `$`を用いた変数の指定
CPS1985$lnwage1 = log(CPS1985$wage)
head(CPS1985$lnwage1)## [1] 1.629241 1.599388 1.897620 1.386294 2.014903 2.570320# [[]]を用いた変数の指定
CPS1985[["lnwage2"]] = log(CPS1985[["wage"]])
head(CPS1985$lnwage2)## [1] 1.629241 1.599388 1.897620 1.386294 2.014903 2.5703207 複数の変数の対数をとる。
- Stataの
foreachのように、複数の変数の対数をとる方法。
# 既存データを一旦削除
rm(list =ls())
# AERパッケージのデータを使う
library(AER)
data("CPS1985")
# 対数変換したい変数のリスト
variables <- c("age", "wage")
for (var in variables) {
# 新しい変数名
new_var_name <- paste0("ln", var)
# 対数変換
CPS1985[[new_var_name]] <- log(CPS1985[[var]])
}
head(CPS1985$lnwage)## [1] 1.629241 1.599388 1.897620 1.386294 2.014903 2.570320head(CPS1985$lnage)## [1] 3.555348 4.043051 2.944439 3.091042 3.555348 3.3322058 Rマークダウンからマークダウンを生成
YAMLの例
---
title: Rメモ
output: github_document
---
9 RマークダウンからLaTexを生成
- 例
---
title: "Rの基礎03 with Penn World Table"
author: "Ayumu Tanaka"
date: '2022-06-14'
output:
pdf_document:
latex_engine: xelatex
documentclass: bxjsarticle
classoption: xelatex,ja=standard
geometry: no
---- Quatroの場合
---
title: "Carolina Caetano and Brantly Callaway"
subtitle: "Difference-in-Differences when Parallel Trends Holds Conditional on Covariates"
author: "資料作成: 田中鮎夢"
format:
pdf:
toc: true
number-sections: false
colorlinks: true
documentclass: bxjsarticle
classoption: xelatex,ja=standard
---10 グローバルチャンクの例
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>",
fig.path = "figures/fig-",
out.width = "100%",
dpi = 300
)
knitr::opts_chunk$set(echo = TRUE, #コードを表示
cache = FALSE, #キャッシュを残さない
fig.path = "figures/fig-",
out.width = "100%",
dpi = 300,
message=FALSE, warning=FALSE)11 パイプ・オペレータの使い方
1+1## [1] 2a <- 1:3 |> sum()
a## [1] 6Base R では |> がパイプ演算子。
1:2 |> mean()## [1] 1.5モダンなR（tidyverse）では %>%
がパイプ演算子。Control + Shift + m がショートカット。
library(dplyr)
1:2 %>% mean()## [1] 1.5ともに左辺が右辺の第一引数となる。
12 グラフ作成・保存の基本
plot(1:5, 6:10)
library(dplyr)
state <- state.x77
state <- data.frame(state)
library(ggplot2)
ggplot(data = state) +
geom_point(mapping = aes(x = Income, y = Population) )
13 回帰分析
install.packages("wooldridge")
install.packages("AER")
library(wooldridge)
data("wage1")
wageModel <- lm(lwage ~ educ + exper + tenure, data = wage1)
summary(wageModel)##
## Call:
## lm(formula = lwage ~ educ + exper + tenure, data = wage1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.05802 -0.29645 -0.03265 0.28788 1.42809
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.284360 0.104190 2.729 0.00656 **
## educ 0.092029 0.007330 12.555 < 2e-16 ***
## exper 0.004121 0.001723 2.391 0.01714 *
## tenure 0.022067 0.003094 7.133 3.29e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4409 on 522 degrees of freedom
## Multiple R-squared: 0.316, Adjusted R-squared: 0.3121
## F-statistic: 80.39 on 3 and 522 DF, p-value: < 2.2e-1614 データ名称を変更する:
assign
assign("data_new", data)
rm(data)15
Stataのforvaluesと同じようにループで連番データを読み込む
# 2019年から2023年の役員データを読み込む
for (i in 2019:2023) {
# Construct file names
file_name <- paste0("yakuin_data_", i, ".csv")
file_name_new <- paste0("yakuin_data_", i, ".RData")
# Read the data
data <- read.csv(file_name, header = T, fileEncoding = "Shift-JIS")
# Assign data to a variable with its original name
assign(paste0("data", i), data)
}16
tapplyの使い方
tapplyを使うと、属性ごとの平均値を計算できる。
library(dplyr)
starwars <- starwars
height_mean <- tapply(starwars$height, starwars$sex, mean, na.rm=TRUE)
head(height_mean)## female hermaphroditic male none
## 171.5714 175.0000 179.1228 131.200017 交差項
?formulaで説明がある。
- a*bは、a + b + a:bと解釈される。
- in%演算子は、 a + b %in% a は a + a:b という式に展開される。
- 演算子 / は省略記法を提供し、a / b は a + b %in% a と等価である。つまり、a + a:bと解釈できる。
- 演算子 - は指定された項を削除するので、(a+b+c)^2 - a:b は a + b + c + b:c + a:c と同じである。
18 変数名の変更
データの読み込み
library(wooldridge)
data("wage1")18.1 方法1
names()関数
names(wage1)[names(wage1)=="educ"] <- "Education"18.2 方法2
colnames(wage1)[colnames(wage1)=="wage"] <- "Wage"18.3 方法3
rename()関数
library(tidyverse)
# Rename the variable "tenure" to "Tenure" and "exper" to "Exper" with `rename`
wage1 <- wage1 %>%
rename(Tenure = tenure, Exper = exper)18.4 方法4
rename_with()関数
# Rename the variables ending with "sq" to "2" with `rename_with()`
wage1 <- wage1 %>%
rename_with(~str_replace(., "sq", "2"), ends_with("sq"))19
mutateの使い方
mutateは既存の変数に操作を加えて、新しい変数を作成する関数。
library(dplyr)
starwars %>%
select(name, mass) %>%
mutate(
mass2 = mass * 2,
mass2_squared = mass2 * mass2
)## # A tibble: 87 × 4
## name mass mass2 mass2_squared
## <chr> <dbl> <dbl> <dbl>
## 1 Luke Skywalker 77 154 23716
## 2 C-3PO 75 150 22500
## 3 R2-D2 32 64 4096
## 4 Darth Vader 136 272 73984
## 5 Leia Organa 49 98 9604
## 6 Owen Lars 120 240 57600
## 7 Beru Whitesun Lars 75 150 22500
## 8 R5-D4 32 64 4096
## 9 Biggs Darklighter 84 168 28224
## 10 Obi-Wan Kenobi 77 154 23716
## # ℹ 77 more rowslibrary(dplyr)
starwars <- starwars
head(starwars) ## # A tibble: 6 × 14
## name height mass hair_color skin_color eye_color birth_year sex gender
## <chr> <int> <dbl> <chr> <chr> <chr> <dbl> <chr> <chr>
## 1 Luke Sky… 172 77 blond fair blue 19 male mascu…
## 2 C-3PO 167 75 <NA> gold yellow 112 none mascu…
## 3 R2-D2 96 32 <NA> white, bl… red 33 none mascu…
## 4 Darth Va… 202 136 none white yellow 41.9 male mascu…
## 5 Leia Org… 150 49 brown light brown 19 fema… femin…
## 6 Owen Lars 178 120 brown, gr… light blue 52 male mascu…
## # ℹ 5 more variables: homeworld <chr>, species <chr>, films <list>,
## # vehicles <list>, starships <list>starwars2 <- starwars %>%
select(name, mass, species) %>%
group_by(species) %>%
mutate(mass_norm = mass / mean(mass, na.rm = TRUE))
head(starwars2)## # A tibble: 6 × 4
## # Groups: species [2]
## name mass species mass_norm
## <chr> <dbl> <chr> <dbl>
## 1 Luke Skywalker 77 Human 0.947
## 2 C-3PO 75 Droid 1.08
## 3 R2-D2 32 Droid 0.459
## 4 Darth Vader 136 Human 1.67
## 5 Leia Organa 49 Human 0.603
## 6 Owen Lars 120 Human 1.48ここで、group_byは、指定した変数ごとに後のmutateで平均を計算することを指示するもので、それ自体ではデータに変化を及ぼさない。
20 x軸のラベルの不要な小数点の除去
20.1 plot()関数
Rのplot関数でグラフ作成すると、x軸のラベルに不要な小数点がついてしまうことがある。例）2005.0, 2005.5, 2006.0
year <- c(2005,2006,2007)
y <- c(1,2,3)
plot(year,y)
x軸の変数の値を増やすことで、小数点を避けられる場合もあるようだが、以下のようにうまくいかない場合もある。
year2 <- c(2005,2006,2007,2008)
y2 <- c(1,2,3,4)
plot(year2,y2)
その他の方法として、x軸のラベルをplot関数では指定せずに、後から、axis関数で追加する方法がある。
year <- c(2005,2006,2007)
y <- c(1,2,3)
plot(year,y, xaxt = "n")
# axisの最初の引数1は、x軸を意味する。
axis(1, 2005:2007)
yearを時間変数とする方法もある。
20.2 ggplot()関数
- コード例
20.2.1 データを作成
x <- 2010:2020
y <- 1:11
data1 <- data.frame(x, y)20.2.2 まず普通のグラフを作成
library(ggplot2)
g1 <- ggplot(data = data1) +
geom_point(mapping = aes(x = x, y = y))
g1 
20.2.3 x軸の小数点を消す
g2 <- g1 + scale_x_continuous(breaks = ~ axisTicks(., log = FALSE))
g2
21 Stata to R
- Stata to R :: CHEAT SHEET チートシート
- Getting Started in R<->Stata 2010年版。昔ながらのRコード。
- R FOR STATA USERS tidyverseを使った現代的なRコード。
- Translating Stata to R data.tableによるデータ処理、fixestを使った推定を中心とした解説。
- Translating Stata to R アルファベット順のコード解説。
22 RMarkdownとQuatroでの目次自動生成
- 「toc: true」を記載。オプションは種々ある。参考
- 「html_document:」のように「html_document」の後ろに「:」を入れる。「:」を入れ忘れるとエラーになるので注意。
例) RMarkdownの場合
---
title: "伝統的な重力方程式の推定"
author: "田中 鮎夢"
date: "2025-04-19"
output:
html_document:
toc: yes
toc_float: yes
number_sections: yes
bibliography: ref.bib
link-citations: yes
---例) Quatroの場合
---
title: "ユニクロの縫製工場の地図 (ggplot2)"
author: "Ayumu Tanaka"
format: html
toc: true
toc_float: true
number-sections: true
editor: visual
---23 多数の固定効果含むモデルの推定
- Laurent Berge and Grant McDermott. fixest – Fast Fixed-Effects Estimation: Short Introduction
24 Excelデータの読み込み、書き出し
library(readxl)
ABCD <- read_excel("ABCD.xlsx")
library("openxlsx")
write.xlsx(ABCD, "ABCD.xlsx")26 特定の文字の前後の文字列の切り出し
例) 「stringr」パッケージを用いて、「//」の後の文字列を切り出す。
library(stringr)
JPN <- str_split_i("The Antique: Secret of the Old Books // ビブリア古書堂の事件手帖", "//", i = -1)27 R MarkdownでのBibTexの使用
- 冒頭に、「bibliography: ref.bib」と書き入れておく。
- 「ref.bib」ファイルに文献情報を記載しておく。
- 本文で文献を引用するときは、「@joseph2022effect」のように「@key」とする。
- 自動的に文献リストが生成される。
- 参考) 4.5 参考文献と引用
例）
---
title: "Bibtexの練習"
author: "田中鮎夢"
date: "2023-12-12"
output:
html_document: default
bibliography: ref.bib
---
@joseph2022effect はハイチの研究
# 参考文献28 here
パッケージhereは、相対パスを使うためのもの。
インストール
install.packages("here")- here パッケージは、トップレベル・ディレクトリからの相対パスを作成する。このパッケージは、ロード時や here() を呼び出すたびに、現在のプロジェクトのトップ・レベルを表示する：
- Posit Cloudでは、hereはトップディレクトリーを”/cloud/project”と認識。
> library(here)
here() starts at /cloud/project- 例えば、“/cloud/project/qss/INTRO/UNpop.csv”を読み込む時は、以下のようにする。
#データの読み込み例
library(here)
a <- read.csv(here("qss", "INTRO", "UNpop.csv"))- 作業ディレクトリ自体をサブフォルダ”/cloud/project/qss/INTRO”に変更するには、以下のようにする。
#作業ディレクトリの変更(1)
library(here)
setwd(here("qss", "INTRO"))
getwd()#作業ディレクトリの変更(2)
library(here)
here::i_am("qss/INTRO/chap01.Rmd") #Rマークダウンがある場所をトップディレクトリからの相対パスとして表現
setwd(here("qss", "INTRO"))
getwd()- Posit Cloudでは、“/cloud/project”がトップディレクトリーなので、hereを使う利点はあまりない。