prophetでガソリン価格の予測

prophetパッケージを利用してみようということで、 とりあえずvignetteにあるものをデータを変えて実行してみた。

ガソリン価格の予測を行う

資源エネルギー庁のHPから給油所小売価格調査(ガソリン、軽油、灯油)の週次データを取得し、 レギュラーガソリン価格の推移を利用して今後のガソリン価格の予測を行う。 当たり前だが、データ依存による影響がみられており、適切なデータ選定が必要なことがわかる。

library(tidyverse)
library(prophet)
library(magrittr)
file_url <- "https://github.com/ogwf/gasoline_value/blob/master/gasoline_prices.xls?raw=true"
tmp <- paste0(getwd(), "/gasoline_prices.xls")
download.file(url = file_url, destfile = tmp, mode = "wb")

gasoline <-
  readxl::read_xls(path = tmp, sheet = "レギュラー", col_types = "numeric") 

names(gasoline) %<>% 
  stringr::str_replace_all(pattern = "[:space:]", "")

gasoline <-
  gasoline %>% 
  filter(!is.na(調査日)) %>% 
  mutate(調査日 = 調査日 %>% as.Date(origin = "1900-01-01"))

これまでのガソリン価格の推移はこのような状況である。

ggplot(gasoline, aes(調査日, 全国)) + geom_path()

f:id:jerrarrdan:20171215000605p:plain

これをprophetを利用して予測してみると。。。

res <- 
gasoline %>% 
  transmute(ds = 調査日, y = 全国) %>% 
  prophet() 
Initial log joint probability = -9.47537
Optimization terminated normally: 
  Convergence detected: relative gradient magnitude is below tolerance
future <- make_future_dataframe(res, periods = 365)

plot(res, predict(res, future))

f:id:jerrarrdan:20171215000630p:plain

prophet_plot_components(res, predict(res))

f:id:jerrarrdan:20171215000653p:plain

2010年以降は、価格の変動が激しいこともありうまくトレースできておらず、近年は下落傾向となっている。

直近一年間のガソリン価格を使って予測してみる。

近年のガソリン価格の変動が激しくうまく予測できていないようなので、 期間を直近一年間に絞り同様に予測してみる。

直近のガソリン価格の推移は下記の通り。

gasoline %>% 
  arrange(desc(調査日)) %>% 
  slice(1:52) %>% 
  ggplot(aes(調査日, 全国)) + geom_path()

f:id:jerrarrdan:20171215000704p:plain

予測結果は

res <-
  gasoline %>% 
  arrange(desc(調査日)) %>% 
  slice(1:52) %>% 
  transmute(ds = 調査日, y = 全国) %>% 
  prophet() 
Initial log joint probability = -2.03264
Optimization terminated normally: 
  Convergence detected: relative gradient magnitude is below tolerance
future <- make_future_dataframe(res, periods = 60)

plot(res, predict(res, future))

f:id:jerrarrdan:20171215000726p:plain

prophet_plot_components(res, predict(res))

f:id:jerrarrdan:20171215000743p:plain

tmapパッケージはしばらく更新がないらしい

先日tmap(version 1.11)の更新通知が来ていた。 たまたまNewsを確認したところ下記のような記載があり、 sfパッケージがspに置き換わるまでしばらく更新はなさそうである。

NOTE: this will be the last version before the major update (in which sf fully replaces sp)

github.com

issuesには現時点で25項目上がっているが、これらは今後解決していくのだろうか。

tmapパッケージはしばらく更新がないらしい

先日tmap(version 1.11)の更新通知が来ていた。 たまたまNewsを確認したところ下記のような記載があり、 sfパッケージがspに置き換わるまでしばらく更新はなさそうである。

NOTE: this will be the last version before the major update (in which sf fully replaces sp)

github.com

issuesには現時点で25項目上がっているが、これらは今後解決していくのだろうか。

国立大学の職位別の人数と年齢構成

文科省では学校教員統計調査というのをやっているらしい。 最新だと平成28年の実施結果がe-statにて公開されているので、可視化してみた。

詳しいことはわからないが、 男性は職位と人数が反比例しているが、 近年の女性活躍させるという政策の下?、女性は助教の数が教授の人数を上回っている。

また、男性教授の年齢構成をみると、ある年齢の人から教授の数が減る(増える) ポイントがあることがわかる。 教授になる人が多い年齢ということなのか、昇進の抑制とかの影響なのかはよくわからない。

  • 職位別の年齢構成(縦棒は各職位の最頻値を表している) f:id:jerrarrdan:20171113002419p:plain

  • 職位別の人数の累積 f:id:jerrarrdan:20171113002429p:plain

高速にggplotする

ここ(Accelerating ggplot2: use a canvas to speed up plots creation)にあるように、 ベースマップを作成しておくことで高速な描画が期待できる。

library(tidyverse)
library(sf)

nc <- read_sf(system.file("shape/nc.shp", package="sf"))
nc_p <- st_sample(nc, size = 100) %>% st_sf()
base_map <- ggplot(data=nc) + geom_sf()

microbenchmark::microbenchmark(
  base_map + geom_sf(data = nc_p),
  ggplot(data=nc) + geom_sf() + geom_sf(data=nc_p)
)
Unit: milliseconds
                                                 expr      min       lq     mean   median       uq      max neval cld
                      base_map + geom_sf(data = nc_p) 2.378383 2.499209 2.937826 2.576802 2.827515 24.75987   100  a 
 ggplot(data = nc) + geom_sf() + geom_sf(data = nc_p) 4.982936 5.213071 5.761455 5.411302 6.064894 10.97817   100   b

行番号(名前)をデータフレームに追加する

これまではtibble::rowid_to_column()を利用してきたが、この処理をもう少し早くする必要が出てきてしらべていたら、dplyr::add_rownames()の存在を知り、速度比較したらあやかったので今後はこれを使うことにする。

library(tidyverse)
library(microbenchmark)
library(data.table)


iris_df <- iris %>% as.data.table()

microbenchmark(iris %>% add_rownames(),
               iris %>% rownames_to_column(),
               iris_df %>% add_rownames(),
               iris %>% mutate(id = seq.int(nrow(.))))
Unit: milliseconds
                                   expr      min       lq     mean   median       uq      max neval cld
                iris %>% add_rownames() 1.218453 1.262064 1.416237 1.310582 1.368730 4.750718   100 a  
          iris %>% rownames_to_column() 2.834875 2.919264 3.122361 2.996101 3.133352 5.571955   100   c
             iris_df %>% add_rownames() 1.384211 1.447833 1.556418 1.512777 1.599620 2.527902   100 a  
 iris %>% mutate(id = seq.int(nrow(.))) 2.348928 2.415761 2.593333 2.479005 2.568303 5.112818   100  b 

追記: dplyr::add_rownames()はduplicatedとなっており、代わりにtibble::rownames_to_column()を利用せよというwarningが出てくる。 少し遅いがtibbleの関数を利用するほうが無難という結論。

グラフをつくる前に読む本

グラフをつくる前に読む本を読んだ。 献本されたというツイッターを多く見たので、読んでみた。

下記に示すグラフの種類についてどのような場合に・どのようなことを表現したいときに利用すべきかが書いてある。

  • 棒グラフ
  • 折れ線グラフ
  • 円グラフ
  • レーダーチャート
  • ヒートマップ
  • 散布図
  • 積み上げグラフ

これまでにさまざまなグラフを見てきていると思うので、 自身の理解している(把握している)使い方と相違なければすらすらと読める。 ある意味当たり前のことが書かれている。ただし、自身のこれまで用いてきた使い方の再確認するという意味では一度読んでおいてよいと思う。

また、それぞれのグラフが用いられた起源について書かれており、どのような経緯でそれぞれのグラフで表現する必要性があり生まれたかが少しわかり面白いと思う。 また著者の方針ですべてのグラフにはオープンにされている実データを用いており、より分かりやすくなっていると思う。

ところどころに著者のデータに対する解釈が書かれており(その解釈が正しいかどうかは別として)、データの結果を単に示すことだけでなくそこからいかに考察しデータを読み解くことが重要であると考えさせられた。

利益が絡むと、よく見えるようにインチキ工夫しようと考えるが、 ぜひ本を読んでまじめなグラフづくりが行われるように願うばかりである。