R で可視化〜時系列データ〜

2021年11月21日2022年1月12日R,Rで可視化,データサイエンス

ggplot2 を使って時系列データの可視化をやってみます。

データは R のベースのデータセットとして用意されている Seatbelts を使います。このデータは1969年から1984年の間で、運転手が重傷または死亡者数やその時のガソリンの値段などが月次で記録された多変量時系列データです。

可視化

データのクラスやどんな変数がいくつあるのかをみてみます。

library(dplyr); library(tidyverse); library(ggplot2) # スタメンライブラリ

ts = datasets::Seatbelts
class(ts); dim(ts); head(ts)

# 実行結果
[1] "mts" "ts" 
[1] 192   8
     DriversKilled drivers front rear   kms PetrolPrice VanKilled law
[1,]           107    1687   867  269  9059   0.1029718        12   0
[2,]            97    1508   825  265  7685   0.1023630         6   0
[3,]           102    1507   806  319  9963   0.1020625        12   0
[4,]            87    1385   814  407 10955   0.1008733         8   0
[5,]           119    1632   991  454 11823   0.1010197        10   0
[6,]           106    1511   945  427 12391   0.1005812        13   0

ここには表示されていませんが、実はこの ts クラスは日付の情報も含んでる場合が多いです。ggplot2 で可視化する際は日付のカラムが欲しいので、別で作成します。この際に使う lubridate パッケージは、日付に関する操作を大変簡単にしてくれる優れものです。

library(lubridate)
date = seq(ymd('1969/01/01'), by = 'month', length.out = nrow(ts))
data = data.frame(date = date, ts) %>% as_tibble()
data

# A tibble: 192 x 9
   date       DriversKilled drivers front  rear   kms PetrolPrice VanKilled   law
   <date>             <dbl>   <dbl> <dbl> <dbl> <dbl>       <dbl>     <dbl> <dbl>
 1 1969-01-01           107    1687   867   269  9059       0.103        12     0
 2 1969-02-01            97    1508   825   265  7685       0.102         6     0
 3 1969-03-01           102    1507   806   319  9963       0.102        12     0
 4 1969-04-01            87    1385   814   407 10955       0.101         8     0
 5 1969-05-01           119    1632   991   454 11823       0.101        10     0
 6 1969-06-01           106    1511   945   427 12391       0.101        13     0
 7 1969-07-01           110    1559  1004   522 13460       0.104        11     0
 8 1969-08-01           106    1630  1091   536 14055       0.104         6     0
 9 1969-09-01           107    1579   958   405 12106       0.104        10     0
10 1969-10-01           134    1653   850   437 11372       0.103        16     0
# … with 182 more rows

1変数の可視化

とりあえず死亡者数をデフォルトでプロットしてみます。

ggplot(data = data, aes(x = date, y = DriversKilled)) + 
  geom_line()

これを基本として、色々といじっていきます。まず、日付の表示形式を”年/月”にして、表示する日付の数を増やし、文字が被らないように斜めにしてみます。

ggplot(data = data, aes(x = date, y = DriversKilled)) + 
  geom_line() + 
  scale_x_date(date_labels = '%Y/%m', date_breaks = '2 year') + 
  theme(axis.text.x = element_text(angle = 50, hjust = 1))

scale_x_date() の中の date_labels で日付の表示形式を指定し、date_breaks で日付を表示する間隔を指定してます（この場合2年ごと）。このとき引数に limit = c(as.Date(1969/01/01), as.Date(1975/01/01)) と書けば、表示させる x軸の範囲を変更して局所的に系列をみることができます。theme() の中の element_text() の angle で傾ける角度、hjust で y 軸方向の位置を指定してます。

DriversKilled が大きいところ（>160）だけ図中に日付をプロットし、シートベルト着用が義務付けられた日（1983/01/31）がわかりやすいように縦線を入れてみます。

data_high = data %>% filter(DriversKilled > 160)
ggplot(data = data, aes(x = date, y = DriversKilled)) + 
  geom_line() + 
  scale_x_date(date_labels = '%Y/%m/%d', date_breaks = '1 year') + 
  theme(axis.text.x = element_text(angle = 50, hjust = 1)) +
  geom_vline(aes(xintercept = as.Date('1983/02/01')), alpha = 0.7, col = 2, lwd = 1.5) + 
  geom_text(data = data_high, aes(x = date, y = DriversKilled, label = date))

geom_vline() でどこに縦線を引くかを指定し、geom_text() で図中のどこにラベルを記入するかを指定してます。x 軸がちゃんと対応しているカラムを用意すれば、このように元データからいくつかデータを切り取って、それを図中のラベル表示用に使うことができます。

データに関する印象を述べると、11~12月に死亡事故が多いようです。シートベルト着用義務化の数ヶ月後までは大幅に減少していますが、後の増加具合を見ると「義務化に伴い、それを取り締まるためにパトロールが増えた」といった原因があるのかもしれません。一応義務化以降に過去と比べると水準が下がったように見えますが、義務化後の観測期間が比較的短いのでまだ何とも言い難いです。

月に関する季節成分を丸めてみるために、DriversKilled を年ごとに合計した値をプロットしてみます。点の数が減ったので、geom_point() も足してみやすくしました。

data_year = data %>% 
  mutate(year = floor_date(`date`, unit = 'year')) %>% 
  group_by(year) %>% 
  summarise(DriversKilled_year = sum(DriversKilled))
ggplot(data = data_year, aes(x = year, y = DriversKilled_year)) + 
  geom_line() + 
  geom_point() + 
  scale_x_date(date_labels = '%Y', date_breaks = '1 year')

年でまとめてみると、シートベルト義務化以降の1983年、1984年は劇的に減少してますね。

次に、月毎の特徴を見るため、月別の死者数を棒グラフでみてみます。

data_month = data %>% 
  mutate(month = month(date, label = TRUE)) %>% 
  group_by(month) %>% 
  summarise(DriversKilled_byMonth = sum(DriversKilled))
ggplot(data = data_month, aes(x = month, y = DriversKilled_byMonth)) + 
  geom_bar(stat = 'identity')

時系列で見た印象通り、12月に近づくにつれて死亡者数が多くなっています。ところが1月になると急に下がっていることがわかります。

理由は分かりませんが、「路面凍結によるスリップ等で事故の数が増える」「日没が早くなるので、視界が悪い時間帯が長くなるため」など、色々と仮説を立てるのには役立ちそうです。

多変量を一気に表示

あまり多くてもごちゃごちゃしてしまうので、今回は front と rear を同時にプロットしてみます。ちなみに、株価のシミュレーションなど、ある系列（確率過程など）がどのような振る舞いをするかなどをみたい場合は多数の系列をいっぺんにプロットしますが、同様の方法でできます。

library(reshape2)
data_melt = data[,c('date', 'front', 'rear')] %>% melt(id = 'date')
ggplot(data = data_melt, aes(x = date, y = value)) + 
  geom_line(aes(col = variable))