今週も特にありません

進捗どうですか?

lubridate::ymdでAll formats failed to parseに遭遇

珍しいケースに遭遇したので、ちょっと調べたことをメモします。

日付列にほんの一部だけが日付のような文字列が入っていて、その他ほとんどがNAが入っている数十万レコードのデータを渡されました。文字列型ではなく、 日付型として扱おうと、lubridate::ymdを用いたら、All formats failed to parseに遭遇しました。

調べてみると、parse_date_time(その派生のymdymd_hmsなど)は最大501個のデータから疎推定で日付フォーマットを判定しているようでした。 lubridate.tidyverse.org

501個が判定できないフォーマットであると、全てパース失敗となるようです。

日付フォーマットが確実に分かっている場合には、フォーマットを指定した上でparse_date_timeの引数exact = TRUEとするか、parse_date_time2を用いる必要があります。

素数を使って選び出すthe best irregular guesser I could come up withな実装はここに書かれています。

lubridate/constants.r at master · tidyverse/lubridate · GitHub github.com

「Rによるディープラーニング」のための情報まとめ

オライリーの「RとKerasによるディープラーニング」も出版され、Rユーザでもディープラーニングに手を出してみようという人が増えている気がします。これは、Rユーザにとってのディープラーニング関連情報のありかをまとめたものになります。

あらためて、RStudioさんが大変素晴らしい仕事をしていただいているため、その最新情報をチェックし続けるだけでも良いのかもしれません。他に良い情報を発信しているところがあれば、追記していきます。

TensorFlow

TensorFlowをRから利用するRStudio謹製パッケージになります。 tensorflow.rstudio.com

GitHub

github.com

Blog

RStudio、Business Scienceなどに所属している方によるためになる記事。 blogs.rstudio.com

Keras

こちらもRStudioによる、KerasをRから利用するものになります。 keras.rstudio.com

GitHub

github.com

Cheat Sheet

github.com

チートシートは日本語版も用意されていています。 github.com

Book

Keras開発者のフランソワによる「PythonとKerasによるディープラーニング」のJ. J. アレールによるR版書籍になります。 www.oreilly.co.jp

www.manning.com

Article

Torch

Rからtorchライブラリを利用するものになります。開発段階で今後に期待できます。 dfalbel.github.io

GitHub

github.com

CNTK

MicrosoftによるCNTKをRから利用するものになります。 microsoft.github.io

GitHub

2019年5月現在積極的には開発されていない模様です。CNTKの人気がイマイチだからなのかもしれませんが、MicrosoftはRユーザに対しての利用拡大を推進していって欲しいところです。 github.com

Article

MXNet

mxnet.incubator.apache.org

Github

github.com

Article

H2O

docs.h2o.ai

Cheat Sheet

github.com

Article

その他

その他、適度に手入れがされているパッケージ。

nnet

CRAN - Package nnet

neuralnet

CRAN - Package neuralnet

github.com

他にも関連パッケージがありますが、開発が止まっているものが多いようです。

tidy時系列データに対する差分計算

以下の記事の通りで、差分計算することが多い方はすでにdplyr::lagを使っていると思います。ここでは、差分計算と変化率、対数収益率を計算する場合についてと、最近少し調べていたtsibbleの中に含まれる関数に関するメモになります。 notchained.hatenablog.com

変化率

tibbletimeパッケージに含まれているFacebookの株価データを利用します。

> library(tidyverse)
> library(tibbletime)
> 
> data(FB)
> FB
# A tibble: 1,008 x 8
   symbol date        open  high   low close    volume adjusted
   <chr>  <date>     <dbl> <dbl> <dbl> <dbl>     <dbl>    <dbl>
 1 FB     2013-01-02  27.4  28.2  27.4  28    69846400     28  
 2 FB     2013-01-03  27.9  28.5  27.6  27.8  63140600     27.8
 3 FB     2013-01-04  28.0  28.9  27.8  28.8  72715400     28.8
 4 FB     2013-01-07  28.7  29.8  28.6  29.4  83781800     29.4
 5 FB     2013-01-08  29.5  29.6  28.9  29.1  45871300     29.1
 6 FB     2013-01-09  29.7  30.6  29.5  30.6 104787700     30.6
 7 FB     2013-01-10  30.6  31.5  30.3  31.3  95316400     31.3
 8 FB     2013-01-11  31.3  32.0  31.1  31.7  89598000     31.7
 9 FB     2013-01-14  32.1  32.2  30.6  31.0  98892800     31.0
10 FB     2013-01-15  30.6  31.7  29.9  30.1 173242600     30.1
# … with 998 more rows

diffは差分計算できない最初の値を除いた値を返すため、mutateで追加しようとしてもエラーが出ます。

> FB %>%
+   select(symbol, date, adjusted) %>%
+   mutate(d = diff(adjusted))
 エラー: Column `d` must be length 1008 (the number of rows) or one, not 1007

dplyr::lagを用いることで差分と変化率は以下のように計算できます。

> FB %>%
+   select(symbol, date, adjusted) %>%
+   mutate(d = adjusted - lag(adjusted), 
+          cr = d / lag(adjusted) * 100)
# A tibble: 1,008 x 5
   symbol date       adjusted       d      cr
   <chr>  <date>        <dbl>   <dbl>   <dbl>
 1 FB     2013-01-02     28    NA      NA    
 2 FB     2013-01-03     27.8  -0.23   -0.821
 3 FB     2013-01-04     28.8   0.99    3.56 
 4 FB     2013-01-07     29.4   0.66    2.29 
 5 FB     2013-01-08     29.1  -0.360  -1.22 
 6 FB     2013-01-09     30.6   1.53    5.26 
 7 FB     2013-01-10     31.3   0.710   2.32 
 8 FB     2013-01-11     31.7   0.42    1.34 
 9 FB     2013-01-14     31.0  -0.770  -2.43 
10 FB     2013-01-15     30.1  -0.850  -2.75 
# … with 998 more rows

対数収益率

対数収益率については、diff(log(FB$adjusted))という感じで計算できるということが巷でよく書かれていますが、こちらも上記と同じ理由でmutateした場合にはエラーになります。

> FB %>%
+   select(symbol, date, adjusted) %>%
+   mutate(ld = diff(log(adjusted)))
 エラー: Column `ld` must be length 1008 (the number of rows) or one, not 1007

dplyr::lagを用いた場合には、多少冗長の感じもしますが、それぞれ対数をとった値の差分を計算することになります。

> FB %>%
+   select(symbol, date, adjusted) %>%
+   mutate(ld = log(adjusted) - log(lag(adjusted)))
# A tibble: 1,008 x 4
   symbol date       adjusted        ld
   <chr>  <date>        <dbl>     <dbl>
 1 FB     2013-01-02     28    NA      
 2 FB     2013-01-03     27.8  -0.00825
 3 FB     2013-01-04     28.8   0.0350 
 4 FB     2013-01-07     29.4   0.0227 
 5 FB     2013-01-08     29.1  -0.0123 
 6 FB     2013-01-09     30.6   0.0513 
 7 FB     2013-01-10     31.3   0.0229 
 8 FB     2013-01-11     31.7   0.0133 
 9 FB     2013-01-14     31.0  -0.0246 
10 FB     2013-01-15     30.1  -0.0278 
# … with 998 more rows

並び替え付きの差分計算

tsibbleパッケージにはdplyr::with_orderを利用したdifferenceという関数が定義されており、時系列で順序が並び替えられていないデータでも並び替えた上で差分計算ができます。

> library(tsibble)
>
> set.seed(12345)
> FB %<>%
+   select(symbol, date, adjusted) %>%
+   slice(sample(nrow(FB)))
> FB %>% 
+   mutate(d = difference(adjusted, order_by = date))
# A tibble: 1,008 x 4
   symbol date       adjusted       d
   <chr>  <date>        <dbl>   <dbl>
 1 FB     2015-11-18    108.   2.64  
 2 FB     2016-07-01    114.  -0.0900
 3 FB     2016-01-15     95.0 -3.40  
 4 FB     2016-07-15    117.  -0.43  
 5 FB     2014-10-27     80.3 -0.390 
 6 FB     2013-08-29     41.3  0.730 
 7 FB     2014-04-17     58.9 -0.780 
 8 FB     2015-01-09     77.7 -0.440 
 9 FB     2015-11-19    106.  -1.51  
10 FB     2016-12-02    115.   0.300 
# … with 998 more rows

結局、この後arrangeを使って並び替えて確認したりすることが多そうなので、あまり使う場面は限られるかもしれません。

tidy時系列データにおける相関計算

時系列データに対して相関を出す場面で、毎回どういう変換するんだっけを調べている気がするためメモします。

corrrパッケージでどのようなことができるのかについては、kazutanさんのページが大変参考になります。 kazutan.github.io

データ

店舗ごとの何らかの商品の売り上げのようなデータがあり、各店舗間の売り上げが相関しているのか分析したいということを想定しています。

library(tidyverse)
library(lubridate)
library(corrr)

set.seed(12345)
store_sales_tbl <- tibble(
  date = rep(seq(ymd("2019-01-01"), ymd("2019-12-31"), by = "1 week"), 5),
  store_name = c(rep("A", 53), rep("B", 53), rep("C", 53), rep("D", 53), rep("E", 53)),
  sales = rpois(53 * 5, 10)
)

データの形式としては、SQLでそのままの抽出してきたような状態が想定されます。

> store_sales_tbl
# A tibble: 265 x 3
   date       store_name sales
   <date>     <chr>      <int>
 1 2019-01-01 A             11
 2 2019-01-08 A             12
 3 2019-01-15 A              9
 4 2019-01-22 A              8
 5 2019-01-29 A             11
 6 2019-02-05 A              4
 7 2019-02-12 A             11
 8 2019-02-19 A              7
 9 2019-02-26 A              8
10 2019-03-05 A             11
# … with 255 more rows

相関計算

このようなデータに対して、tidyverse的にstats::corを適用するには意外と面倒です。corrr::correlateが適用できる形に持っていくのが早いです。

ここでは、時系列データの相関について、細かい考慮することはしませんが、変化率を計算する関数などを用意しておきます(ファイナンス関連だと対数変化率などの関数を利用することを想定します)。

cr <- function(x) {
  (x - lag(x)) / lag(x) * 100
}

あとは、このデータを一度、横持ちのデータに変換してから、変化率の系列に置き換え、corrr::correlateに繋げれば相関行列を算出できます。

store_sales_cor_df <- store_sales_tbl %>%
  spread(store_name, sales) %>%
  transmute_if(is.integer, cr) %>%
  correlate()

算出された相関行列を確認します。

> store_sales_cor_df
# A tibble: 5 x 6
  rowname        A        B        C        D       E
  <chr>      <dbl>    <dbl>    <dbl>    <dbl>   <dbl>
1 A        NA        0.118   -0.0512  -0.102    0.241
2 B         0.118   NA       -0.0853   0.0214   0.128
3 C        -0.0512  -0.0853  NA       -0.304   -0.212
4 D        -0.102    0.0214  -0.304   NA        0.158
5 E         0.241    0.128   -0.212    0.158   NA    

可視化

ここまでくると、ggplot2で可視化することもできますが、corrr::rplotを使えば、簡単に相関行列を可視化することができます。

store_sales_cor_df %>%
  rplot(print_cor = TRUE)

f:id:masaqol:20190414215924p:plain

corrplotの方が使い慣れているという場合には、corrr::as_matrixを利用すれば可視化まで持っていけます。

library(corrplot)

store_sales_cor_df %>%
  as_matrix() %>%
  corrplot.mixed(lower = "number", upper = "square")

f:id:masaqol:20190414215956p:plain

github.com

lubridate::floor_dateで任意の曜日スタートで切り下げ

タイトル通りのメモ。

library(tidyverse)
library(lubridate)

set.seed(12345)
order_tbl <- tibble(
  order_date = seq(ymd("2019-03-01"), ymd("2019-03-31"), by = "1 day"),
  order_num = rpois(31, 5)
)

日曜日と月曜日スタートする場合を追加します。

> order_tbl %>%
+   mutate(wd = weekdays(order_date),
+          sunday_start = floor_date(order_date, unit = "week"),
+          monday_start = floor_date(order_date, unit = "week", week_start = 1))
# A tibble: 31 x 5
   order_date order_num wd     sunday_start monday_start
   <date>         <int> <chr>  <date>       <date>      
 1 2019-03-01         6 金曜日 2019-02-24   2019-02-25  
 2 2019-03-02         8 土曜日 2019-02-24   2019-02-25  
 3 2019-03-03         6 日曜日 2019-03-03   2019-02-25  
 4 2019-03-04         8 月曜日 2019-03-03   2019-03-04  
 5 2019-03-05         5 火曜日 2019-03-03   2019-03-04  
 6 2019-03-06         3 水曜日 2019-03-03   2019-03-04  
 7 2019-03-07         4 木曜日 2019-03-03   2019-03-04  
 8 2019-03-08         5 金曜日 2019-03-03   2019-03-04  
 9 2019-03-09         6 土曜日 2019-03-03   2019-03-04  
10 2019-03-10        11 日曜日 2019-03-10   2019-03-04  
# … with 21 more rows

つまり、week_start = 1を指定すれば、月曜日スタートで切り下げられます。デフォルトでは、week_start = getOption("lubridate.week.start", 7)が設定されていますので、日曜日スタートになります。

round_datecelling_dateを含むfloor_dateの実装場所は以下になります。 github.com

lubridate::rollbackの使い所

日次データが存在して、その日次データとそれを月ごとにまとめ上げた月次データの可視化を行いたい場合があります。

library(tidyverse)
library(lubridate)

set.seed(12345)
order_tbl <- tibble(
  order_date = seq(ymd("2018-01-01"), ymd("2018-12-31"), by = "1 day"),
  order_num = rpois(365, 5)
)

日次データとしてプロットしたい場合は、以下のようにできます。

order_tbl %>%
  ggplot(aes(x = order_date, y = order_num)) +
  geom_point(size = 2, colour = "blue") +
  geom_line(size = 1, colour = "blue", alpha = 0.5) +
  scale_x_date(date_breaks = "3 months", date_minor_breaks = "1 months", date_labels = "%Y-%m") +
  scale_y_continuous(breaks = seq(0, 20, by = 2), minor_breaks = 1:20) +
  theme(axis.text = element_text(size = 15),
        axis.title = element_text(size = 20))

f:id:masaqol:20190305010115p:plain

次に、月ごとにまとめ上げたデータは以下のように計算できます (年と月のまとめ上げでも可能)。

order_tbl %>%
  mutate(y = year(order_date),
         m = month(order_date),
         ym = str_c(y, "-", sprintf("%02d", m))) %>%
  group_by(ym) %>%
  summarise(monthly_order_sum = sum(order_num))
#> # A tibble: 12 x 2
#>    ym      monthly_order_sum
#>    <chr>               <int>
#>  1 2018-01               161
#>  2 2018-02               133
#>  3 2018-03               161
#>  4 2018-04               163
#>  5 2018-05               176
#>  6 2018-06               171
#>  7 2018-07               156
#>  8 2018-08               169
#>  9 2018-09               161
#> 10 2018-10               144
#> 11 2018-11               146
#> 12 2018-12               170

しかし、これだと年月を表すカラムがcharactor型のため、上記の日次データを可視化した際に利用したscale_x_dateの部分を流用して同じスケールで可視化できません。

そこで、lubridate::rollbackを利用すれば、以下のようにできます。

order_tbl %>%
  mutate(first_date = rollback(order_date, roll_to_first = TRUE)) %>%
  group_by(first_date) %>%
  summarise(monthly_order_sum = sum(order_num)) %>%
  ggplot(aes(x = first_date, y = monthly_order_sum)) +
  geom_point(size = 2, colour = "blue") +
  geom_line(size = 1, colour = "blue", alpha = 0.5) +
  scale_x_date(date_breaks = "3 months", date_minor_breaks = "1 months", date_labels = "%Y-%m") +
  scale_y_continuous(breaks = seq(0, 200, by = 10), minor_breaks = seq(0, 200, by = 5)) +
  theme(axis.text = element_text(size = 15),
        axis.title = element_text(size = 20))

f:id:masaqol:20190305010155p:plain

ここでは、その月の最初の日のカラムを追加して、その日で月次にまとめ上げて可視化しています。こうすることで、日次データで利用したscale_x_dateの部分も流用できます。

rollback関数の中身自体は、day関数を使って、日付部分に1を入れるというシンプルな操作を行う便利関数になっています。 github.com

Rでモデルの評価指標を計算する

機械学習モデルの分類精度を評価したり、時系列モデルのn期先の予測精度を評価したりする際に、地味に評価指標の計算式をネットで調べて自前でRの関数として実装したりしていました。

しかし、(実際にはよく調べていなかっただけですが、)すでに基本的な指標はパッケージとしてまとめられていました。パッケージ名は、至ってそのままでMetricsパッケージになります。

cran.r-project.org

github.com

F1、Rrecision、Recallをはじめ、MAPEやRMSEなど基本的な評価指標は用意されています。

実装はシンプルで、関数に入力する実際の値と予測した値の個数が異なる場合でも、警告メッセージは出ますが値は出力される感じなので、その点は少し注意が必要です。