1. Rによる二次的データの分析

1.1 はじめに

Rの使い方をまとめてみようと思う。Rの使い方についてはおよそ6カ月間毎日取り組んだので、自分が使う可能性がある手法についてはだいたいマスターすることができたはずである。理解したことをまとめてみて検討してみようと思う。Rと統計学に関する知識を整理してみる。

また、このファイル（Rmd形式）を編集することで、RStudioが原稿などを書くのに使えるかどうかを確認してみたい。グラフや表のある原稿はRStudioが便利のようである。

（html_notebookを使って見たのだが使い方がややこしいので、元のhtml_documentに戻した。）
（次の行は空けておかないと##が有効にならない。空白行は改行と同じ効果。改行は半角スペース２個でも可能。）

1.2 long formatの活用

long formatは、縦持ち形式とか縦長形式とか訳される。これについて最初に検討してみよう。 いままで使ってきて、視覚的にとらえやすいのは変数列が横長に広がるwide format1であるが、long formatの方が便利な面があるということなので試している。wide formatをlong formatに変換する関数について確認してみる（gatherとpivot_longer）。（注は、ドキュメントの最後に表示される。）

ここで、 RStudioで文章を書いたり数式を示したりする場合の約束事について確認してみる。 まず、数式で使う記号の確認。**が^として使える。Rでの実行結果で確認できるように、２の3乗は次のどちらでもよい。（この部分は、実際には次のように入力している。）

\*\*が^として使える。**２の3乗**は次のどちらでもよい。

> 2**3
[1] 8
> 2^3
[1] 8

字体について。 2つのアスタリスク(*)で囲むとボールドになる。ひとつのアスタリスクで囲むとイタリック。

記号**を使うとボールド体の指定のはじまりと解釈されてしまうばあいがある（後から**が出てくる場合）。それを避けるためにはその前に記号\を前に1つ入れるとよい。（この文の実際の入力は以下のとおり。）

記号\*\*を使うとボールド体の指定のはじまりと解釈されてしまうばあいがある（後から\*\*が出てくる場合）。それを避けるためにはその前に記号\\を前に*1つ*入れるとよい。

「`」の使い方。この記号で囲むと囲まれた部分をコマンドなどとして表記できる。複数行にわたるものの場合は、最初と最後に、この記号を３個並べた「```」を使う。

1.3 gather関数の使い方

wide formatをlong formatに変換するこの関数は、次のような使い方をする。

  1. gather(df,"key","value",-x)
  2. gather(df,"key","value",y1,y2,y3)
gather(
  data,
  key = "key",
  value = "value",
  ...,
  na.rm = FALSE,
  convert = FALSE,
  factor_key = FALSE
)

gather関数を使うとwide formatをlong formatに変換することができる。これはtydyrパッケージに含まれている。 (1)の”key”と”value”のところは、別の名前をつけてよい。新変数であり、日本語を使ってもよい。新変数は引用符で囲むことが推奨されている。

(1)の場合、yがこの形式では、設定されていない。マイナス記号を付けた変数は、インデックス変数として扱われ、それ以外が新変数の”key”と”value”を使って集められる。

(2)の場合、y1,y2,y3以外は、インデックス変数として扱われる。 y1,y2,y3をy1:y3と表現できる。どちらで表記してもよい。ここに出てこない変数(ここではx)はインデックス変数として扱われる。

想定されているのは以下のようなデータフレーム（オブジェクト名:df）である。これはx,y1,y2,y3の3変数からなるデータフレームである。

library(tidyr)
library(ggplot2)

x <- 1:10
y1 <- x**1.5
y2 <- x**2
y3 <- x**2.5
df <- data.frame(x,y1,y2,y3)
knitr::kable(df)
x y1 y2 y3
1 1.000000 1 1.000000
2 2.828427 4 5.656854
3 5.196152 9 15.588457
4 8.000000 16 32.000000
5 11.180340 25 55.901699
6 14.696938 36 88.181631
7 18.520259 49 129.641814
8 22.627417 64 181.019336
9 27.000000 81 243.000000
10 31.622777 100 316.227766

gather関数を使い横長のdfを縦長のdf_longに変えるには以下のようにする。変換されたデータフレームはdf_longという名前にしている。長くなるので、13行目まで表示する。

df_long <- gather(df,"key","value",-x)
knitr::kable(df_long[1:13,])
x key value
1 y1 1.000000
2 y1 2.828427
3 y1 5.196152
4 y1 8.000000
5 y1 11.180340
6 y1 14.696938
7 y1 18.520259
8 y1 22.627417
9 y1 27.000000
10 y1 31.622777
1 y2 1.000000
2 y2 4.000000
3 y2 9.000000

あるいは、2番目の形式で以下のようにしても同じ結果である。

df_long <- gather(df,"key","value",y1,y2,y3)
#knitr::kable(df_long[1:13,])

調べてみてわかったことだが、 gather関数はすでに開発終了(retired lifecyle)で、以下のように、pivot-longerに移行することが推奨されている。

Gather columns into key-value pairs
Description
Retired lifecycle
Development on gather() is complete, and for new code we recommend switching to pivot_longer(), which is easier to use, more featureful, and still under active development. df %>% gather("key", "value", x, y, z) is equivalent to df %>% pivot_longer(c(x, y, z), names_to = "key", values_to = "value")

See more details in vignette("pivot").

Usage
gather(
  data,
  key = "key",
  value = "value",
  ...,
  na.rm = FALSE,
  convert = FALSE,
  factor_key = FALSE
)

次のような式が推奨されている。これに従って書き換えてみると同じようなが得られる（並べ方には違いがある）。gather関数の場合と同様に、“key”と”value”のところは変更可能であり、日本語で”凡例”,“値”などとしてもよい。(keyのところは凡例にすれば凡例が凡例と表示されるので便利であるが、複数の変数をまとめる新変数に「凡例」という名前が付くのは感覚的には受け入れにくいことである。その意味ではあとで、labsのなかで指定する方がよいかもしれない。)

df %>% pivot_longer(c(x, y, z), names_to = "key", values_to = "value")
library(tidyr)

df_long <-df %>%
  pivot_longer(c(y1, y2, y3), names_to = "key", values_to = "value")

df_long
## # A tibble: 30 × 3
##        x key   value
##    <int> <chr> <dbl>
##  1     1 y1     1   
##  2     1 y2     1   
##  3     1 y3     1   
##  4     2 y1     2.83
##  5     2 y2     4   
##  6     2 y3     5.66
##  7     3 y1     5.20
##  8     3 y2     9   
##  9     3 y3    15.6 
## 10     4 y1     8   
## # ℹ 20 more rows

pivot_longerを使った場合には、xを軸にしてy1,y2,y3の値が示されている。

データを用いてグラフを描く場合、 wide formatの新変数名に対応してggplot2でのaesの設定を変える。 図１では、ggplot2のaesでx=x,y=value,col=keyと設定している。 x軸、y軸、凡例のタイトルについてlabsの中で設定してはいない。 図２では、凡例及び値という変数をpivot_longer関数の中で設定し、それらを使ってggplot2のaesの中でx=x,y=値,col=凡例と指定している。y軸のラベル、凡例のタイトルにそれらが自動的に反映している。

library(patchwork)

p1 <- ggplot(df_long,aes(x=x,y=value,col=key))+
  geom_line()+
labs(title="図1",subtitle = "aes(x=x,y=value,col=key)")+
theme_minimal(base_family="HiraKakuProN-W3")


df_long <-df %>%
  pivot_longer(c(y1, y2, y3), names_to = "凡例", values_to = "値")

p2 <- ggplot(df_long,aes(x=x,y=値,col=凡例))+
  geom_line()+
labs(title="図2",subtitle = "aes(x=x,y=値,col=凡例)")+
theme_minimal(base_family="HiraKakuProN-W3")
df_long <- gather(df,"凡例","y",y1:y3)

p3 <- ggplot(df_long,aes(x,y,linetype=凡例))+
  geom_line()+
labs(title="図3",subtitle = "aes(x,y,linetype=凡例)")+
theme_minimal(base_family="HiraKakuProN-W3")

p4 <- ggplot(df_long,aes(x,y,linetype=凡例))+
  geom_line(col="blue")+
labs(title="図4",subtitle = "aes(x,y,linetype=凡例)")+
theme_minimal(base_family="HiraKakuProN-W3")

図３では、新変数を凡例、yと名付け、gather(df,"凡例","y",y1:y3)を実行し、さらにggplot2のaesでaes(x,y,linetype=凡例)と指定している。これは、aes(x=x,y=y,linetype=凡例)ということと同じである。この図3は、モノクロで線の種類を使い分けている。図４では、さらに、aesという形式ではなく、geom_lineについて、geom_line(col="blue")として青色を指定している。１色だけであるが、様々な色を線に指定することができる。

1.4 ggplot2でグラフを描く

ggplot2を使うと、さまざまなグラフを自由に描くことできる。方程式を使って曲線と直線を同じグラフに描いてみたり、軸を反転させたりしてみよう。(図2は図1の軸を反転させたものである。)

図１では、新変数としてbucketとyを使っている。 そしてggplot2のaesでx,y,col=bucketと指定している。 しかし、labsの中でcol="凡例"と設定しているので凡例は凡例と表示される。 labsでタイトル、軸、凡例などの名前を指定できる。 y軸のラベルはここでは設定されていないが、その上のaesの設定でy=yということなのでy軸にはｙが表示される。

図2は、coord_flip()を使ってｘ軸とｙ軸を交換している。また、モノクロで線を描いている。新変数は凡例とyであり、linetype=凡例とaesの中で設定している。なお、図2では、captionを使い、グラフの右下に説明を付け加えている。キャプションは以下のように設定している。2行に分けて表示されるように「\n」で改行をおこなっている。

caption = "図1の軸を交換したもの。linetypeに凡例\nをわりあて、点ではなく線を描く。"
x <- ppoints(100)
y1 <- x**(-0.5)
y2 <- -8*x**2+14.23
y3 <- -14*x + 14.23
df <- data.frame(x,y1,y2,y3)

df_long <- gather(df,"bucket","y",-x)

library(patchwork)

p1 <- ggplot(df_long,aes(x,y,col=bucket))+
  geom_point()+
labs(title="図1",subtitle="y1 = x**(-0.5)\ny2 = -8*x**2+14.23\ny3 = -14*x + 14.23",col="凡例")+
theme_minimal(base_family="HiraKakuProN-W3")


df_long <- gather(df,"凡例","y",-x) # 新変数を"凡例","ｙ"と設定。
p2 <- ggplot(df_long,aes(x,y,linetype=凡例))+
  geom_line()+
labs(title="図2",subtitle="y1 = x**(-0.5)\ny2 = -8*x**2+14.23\ny3 = -14*x + 14.23",caption = "図1の軸を交換したもの。linetypeに凡例\nをわりあて、点ではなく線を描く。")+
theme_minimal(base_family="HiraKakuProN-W3")+
  coord_flip()

1.5 SSDSE

colnamesで元のデータフレームの変数名を変更すればグラフの凡例にも反映する(図1、図2)。

df <- data.frame(pop,d22,OutOverIn)
colnames(df) <- c("pop","財政力指数","転出超過率")
df_long <- gather(df,"bucket","y",-pop)

以下のcolour(colでもよい)の指定に注意。

labs(title = "人口と財政力指数、転出超過率",subtitle = "2つのデータ系列の同時表示",x="人口",y="値",colour="凡例")+

labsの中のこの場合は、凡例のタイトルが指定される。ここでは「凡例」としたのでそのように表示されている。指定しなければ新変数bucketとなる。

図３では、labsの中でcaptionを指定している。その内容はグラフの右下に表示される。

library(tidyr)
library(ggplot2)

dat <- read.csv("./data/dat_ssdse.csv")
pop <- dat$A1101
d22 <- dat$D2201
OutOverIn <- dat$OutOverIn

df <- data.frame(pop,d22,OutOverIn)
colnames(df) <- c("pop","財政力指数","転出超過率")

df_long <- gather(df,"bucket","y",-pop)

library(patchwork)

p1 <- ggplot(df_long,aes(pop,y,col=bucket))+
  geom_point()+
labs(title = "図１. 人口と財政力指数、転出超過率",subtitle = "2つのデータ系列の同時表示",x="人口",y="値",colour="凡例")+
theme_minimal(base_family="HiraKakuProN-W3")

p2 <- ggplot(df_long,aes(pop,y,col=bucket))+
  geom_point()+
xlim(0,1e+06)+
labs(title = "図２. 人口と財政力指数、転出超過率",subtitle = "2つのデータ系列の同時表示",x="人口",y="値",col="凡例")+
theme_minimal(base_family="HiraKakuProN-W3")

df <- data.frame(pop,d22)
df_long <- gather(df,"key","value",-pop)
p3 <- ggplot(df_long,aes(x=pop,y=value,colour="財政力指数"))+
  geom_point(colour="blue")+
xlim(0,1e+06)+
labs(title = "図３. 人口と財政力指数",subtitle = "図２の赤い部分だけを表示",x="人口",y="財政力指数",colour="凡例",caption = "2015年度のデータ")+
theme_minimal(base_family="HiraKakuProN-W3")

library(tidyr)
library(ggplot2)

dat <- read.csv("./data/dat_ssdse.csv")

x <- dat$Dnsty # 可住地人口密度
y1 <- dat$D2201 # 財政力指数
y2 <- dat$r_OnePerson # 単独世帯率
y3 <- dat$OutOverIn # 転出超過率

df <- data.frame(x,y1,y2,y3)
df_long <- gather(df,"凡例","y",-x)

ggplot(df_long,aes(x,y,col=凡例))+
  geom_point()+
labs(title = "可住地人口密度と財政力指数等",subtitle = "y1: 財政力指数\ny2: 単独世帯率\ny3: 転出超過率",x="可住地人口密度",y="y")+
theme_minimal(base_family="HiraKakuProN-W3")

2. 散布図を描く

2.1 散布図(全国)

点の色は、geom_point(colour="red")などとして指定できる。

library(tidyr)
library(ggplot2)
library(ggrepel)

dat <- read.csv("./data/dat_ssdse.csv")

x <- dat$Dnsty # 可住地人口密度
y1 <- dat$D2201 # 財政力指数

df <- data.frame(x,y1)
df_long <- gather(df,"凡例","y",-x)

ggplot(df_long,aes(x,y))+ 
  geom_point(colour="seagreen")+
labs(title = "全国",subtitle = "y: 財政力指数",x="x: 可住地人口密度",y="y")+
theme_minimal(base_family="HiraKakuProN-W3")

library(tidyr)
library(ggplot2)

dat <- read.csv("./data/dat_ssdse.csv")

x <- dat$Dnsty # 可住地人口密度
y1 <- dat$D2201 # 財政力指数

df <- data.frame(x,y1)
df_long <- gather(df,"凡例","y",-x)

ggplot(df_long,aes(x,y,linetype="凡例"))+ 
  geom_point()+
xlim(0,10)+
labs(title = "全国（xの範囲を限定）",subtitle = "y: 財政力指数",x="x: 可住地人口密度",y="y")+
theme_minimal(base_family="HiraKakuProN-W3")

2.2 雛形（全国用）

以下は、colnames(dat)を実行した結果であり、変数名のリストである。これによって各変数がデータフレームの何列目にあるかを調べる。

##  [1] "X"                   "Code"                "Prefecture"         
##  [4] "Municipality"        "A1101"               "A6108"              
##  [7] "B1103"               "A5101"               "A5102"              
## [10] "D2201"               "F1102"               "F2701"              
## [13] "A710201"             "A710202"             "A1303"              
## [16] "A4101"               "A4200"               "A710101"            
## [19] "A810101"             "A810102"             "A810103"            
## [22] "A810104"             "A810105"             "D2205"              
## [25] "D2206"               "Dnsty"               "In"                 
## [28] "Out"                 "OutOverIn"           "sc"                 
## [31] "r_instHm"            "r_nuclear"           "r_nuclear_relativeH"
## [34] "r_OnePerson"         "aged"                "naturalInc"

csvファイルを読み込み、datというオブジェクトとして保存する。Codeとともに、2変数を取り出す。gather式の中で、Codeとx(横軸となる変数)をインデックス変数として指定する。以下の例では、新変数「凡例」に 「sc」が入れられる。

library(tidyr)
library(ggplot2)

dat <- read.csv("./data/dat_ssdse.csv")
colnames(dat)[c(2,26,30)]
## [1] "Code"  "Dnsty" "sc"
df_long <- gather(dat[,c(2,26,30)],"凡例","y",-Code,-Dnsty)
head(df_long)
##   Code     Dnsty 凡例        y
## 1 1100 44.474828   sc 44.14595
## 2 1202 21.566448   sc 87.57994
## 3 1203 15.223374   sc 80.90107
## 4 1204  9.649514   sc 91.29442
## 5 1205 20.467761   sc 85.60854
## 6 1206  5.653251   sc 87.12453
tail(df_long)
##       Code      Dnsty 凡例        y
## 1736 47360  1.3077586   sc 99.32341
## 1737 47361  1.9431220   sc 99.46401
## 1738 47362 11.6450321   sc 34.29054
## 1739 47375  0.7514160   sc 98.22866
## 1740 47381  0.7072351   sc 94.78186
## 1741 47382  1.5671769   sc 99.62035
library(patchwork)
p1 <- ggplot(df_long,aes(Dnsty,y,col=凡例))+
  geom_point(colour="blue")+
labs(title = "散布図(1):",subtitle = "可住地人口密度と自地域就業率(全国)",x="可住地人口密度",y="自地域就業率",caption = "2015年")+
theme_minimal(base_family="HiraKakuProN-W3")
p2 <- ggplot(df_long,aes(Dnsty,y,col=凡例))+
  geom_point(colour="seagreen")+
labs(title = "散布図(2)",subtitle = "x: 可住地人口密度\ny: 自地域就業率",x="x",caption = "2015年")+
theme_minimal(base_family="HiraKakuProN-W3")
p3 <- ggplot(df_long,aes(Dnsty,y,col=凡例))+
  geom_point(colour="skyblue")+
labs(title = "散布図(3)",subtitle = paste("x=",colnames(dat)[26],"\ny=",colnames(dat)[30]),x="x",caption = "2015年")+
theme_minimal(base_family="HiraKakuProN-W3")

散布図(1)は、できるだけわかりやすく説明を日本語で入れるもの。散布図(3)は、できるだけ自動的に図を作成する場合の雛形。サブタイトルの部分が自動的に挿入される。データフレームの変数名が反映する。

次の散布図(4)のコードと実行結果はつぎの通り。

df_long <- gather(dat[,c(2,31,33)],"凡例","y",-Code,-r_nuclear_relativeH)

ggplot(df_long,aes(r_nuclear_relativeH,y,col=凡例))+
  geom_point(colour="skyblue")+
ylim(0,20)+
labs(title = "散布図(4)",subtitle = paste("x=",colnames(dat)[33],"\ny=",colnames(dat)[31]),x="x",caption = "2015年")+
theme_minimal(base_family="HiraKakuProN-W3")

下のボールドの部分を書き換えると同じようなグラフが描ける。取りあげる2変数の列番号、x軸に指定する変数の名前、サブタイトル用の変数番号の部分である。

df_long <- gather(dat[,c(2,31,33)],“凡例”,“y”,-Code,-r_nuclear_relativeH)
ggplot(df_long,aes(r_nuclear_relativeH,y,col=凡例))+
geom_point(colour=“skyblue”)+
ylim(0,20)+ labs(title = “散布図(4)”,subtitle =
paste(“x=”,colnames(dat)[33],“=”,colnames(dat)[31]),x=“x”,caption = “2015年”)+
theme_minimal(base_family=“HiraKakuProN-W3”)

2.4 福岡県

## Warning: ggrepel: 20 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps

## Warning: ggrepel: 34 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps

## Warning: ggrepel: 25 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps

## Warning: ggrepel: 14 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps

## Warning: ggrepel: 13 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps

2.5 群馬県

## Warning: ggrepel: 1 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps

## Warning: ggrepel: 1 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps