Prophetでシンプルなモデルから始めるWeb系データの時系列予測

はじめまして、eurekaのBIチームでデータアナリストをしているクールガイ田中です。このtech blogで記事を書くのは初めてです。緊張します。

好きなアニメは「さくら荘のペットな彼女」です。よろしくお願いします。

この記事は eureka Engineering Advent Calendar 2017 — Qiitaの15日目の記事です。
14日目は コジスティックサーチことeurekaのESマスター小島さん の「【応用編】Elasticsearchの検索クエリを使いこなそう」でした。

はじめに

僕は、普段の業務ではSQLのクエリを書いてデータ抽出したり、GASをいじってスプレッドシートで色々やったりしていることが多いのですが、先日、時系列データの予測依頼が舞い込んできました。

Tokyo.Rでhoxo_mさんの発表を聞いてからずっと、Prophetを業務で使ってみたいと思っていた僕は、ここぞとばかりにProphetを使ってみることにしました。
 
 
今回の記事では、僕がそのとき実際にやった作業の流れをふまえつつ、Web系やアプリ系サービスの実データに近そうなものをProphetで予測するためのシンプルなモデルを作成してみたいと思います。これから初めてProphetを触ってみようという方の参考になれば幸いです。

Prophetとは何か

facebook/prophet
 
今さら説明する必要はないかもしれませんが、ProphetとはFacebookから今年リリースされた時系列予測ツールです。RとPythonから使えて、とても簡単に時系列データの予測ができます。9月にv0.2になって、予測に外部変数も反映させることができる「Extra Regressors」など新たな機能が追加されたりもしているようです。

Prophetの予測モデルは以下の式のようになっています。

それぞれ、g(t)はトレンドを表す関数、s(t)が周期性を表す関数、h(t)が祝日やイベントの効果を表す関数、そして誤差項εです。
これら時間tの関数の和としてモデルを表現できることがProphetの特徴であり、それは平滑化スプラインや一般化加法モデルのように、時系列予測を曲線フィッティングの問題として取り扱うということで…
みたいな話らしいのですが、しかしまあ、今回はそういう難しい話は置いておきましょう。Prophetの良いところは何といっても「使うのが簡単」なところです。この簡単さによって、時系列分析にそこまで詳しくない人でもドメイン知識などを駆使して、良い予測を作ることができるわけです。使っていくうちに色々といじりたくなってくると思うので、必要な知識はその時に勉強すれば良いと思います。
 
というわけで、本記事ではあくまでシンプルな予測モデルを上手いこと作るというところに主眼を置いて進めていきたいと思います。

Prophetでシンプルな予測

そもそもどんな依頼だったのか

本題に入っていきたいと思います。
元々、今回きた予測依頼の内容は、主に以下のようなものでした。

  • 自社サービスへの登録者数の時系列予測をしたい
  • 自社サービスの累計会員数の伸びを予測したい
  • 自社サービスのDAUの時系列予測をしたい

これらの目的に対して、Prophetで予測をしていったわけですが、今回はExtra Regressorsなどの新機能は使わず、予測したいデータと日付を基本データとして、そこにドメイン知識を組み込んでいくというシンプルなProphetの使い方で予測をおこないました。
社内でのProphetの利用自体が初めてだったのと、僕自身まだまだ時系列データの扱いや時系列分析に詳しくないため、まずはシンプルでわかりやすいものから始めることにしました。

どうやったのか

実際にやったProphetでの予測モデル作成の流れは、以下になります。シンプルです。

  1. デフォルトで突っ込んで予測してみる
  2. シーズナリティの有無を考えて調整する
  3. 新機能リリースなどの変化点をモデルに組み込んでみる
  4. (モデルの評価 ※今回はここの話はしません)

予測に使った実際のデータやその処理過程をお見せすることはできないので、本記事用にサンプルデータを用意してみました(このサンプルデータの準備がかなり辛かった…)。
このデータを、ある架空のWebサービスの日次での会員登録者数として、時系列の予測を実際にしてみようと思います。

こんな感じのデータありそうですよね。

しかし、簡単にプロットしてみただけでは、よくわかりませんね。
では、上記の手順を順番になぞってProphetでの予測をしてみましょう。何かわかるかもしれません。

1. デフォルトで突っ込んで予測してみる

Rからprophetパッケージを使っていきます。
データはすでにdaily_dataという名前で読み込み済みとします。

# インストールはこれで一発です
install.packages('prophet')
# パッケージの読み込み
library(prophet)
# カラム名はds,yのみ受け付けられる
colnames(daily_data) <- c("ds", "y")
# 予測モデルの作成
m <- prophet(daily_data)
# 予測する未来のデータフレームを用意
future_df <- make_future_dataframe(m, 365)
# 予測してみる
pred <- predict(m, future_df)
# 図示する
plot(m, pred) + ggtitle("Prophet Forecast Plot")
たったこれだけで、未来の予測値とプロットを得ることができました。特に何も考えずにデータを突っ込んだだけなので、とても簡単でした。
デフォルトでもなんとなく予測はできそうですが、もうちょっといじってみたいところです。
2. シーズナリティの有無を考えて調整する
次に、上で作ってきたデフォルトのモデルを少しいじってみましょう。
Prophetでは、時系列データの周期性として、主に
  • 年周期 (yearly.seasonality)
  • 週周期 (weekly.seasonality)
の2つを簡単に指定することができます。
このシーズナリティを表しているのは、冒頭のモデル式のs(t)の部分です。年周期は月や季節ごとに登録者数に違いがあるのか、週周期は曜日ごとに登録者数に違いがあるのかを考慮します。
早速、モデルを修正してみましょう。このWebサービスでは、年末年始などの特定の時期と土日に登録者数が多くなるということが経験的に知られているとします。このことをモデルに反映してみましょう。
  
Prophetでモデルをいじるときには、基本的に始めのprophet()関数で引数を指定していきます。
# 変わったのはここです。
m <- prophet(daily_data, yearly.seasonality = TRUE, weekly.seasonality = TRUE)
future_df <- make_future_dataframe(m,365)
pred <- predict(m,future_df)
plot(m,pred) + ggtitle("Prophet Forecast Plot")
# 変わったのはここです。
prohet_plot_components(m, pred)
先ほどと同じようなプロットが返ってきましたね。これはなぜかというと、Prophetはデフォルトで周期性を検知して、良い感じにモデルを作ってくれるからです。便利ですね。
さらに、先ほどよりも図が増えています。2〜4枚目のプロットは、1枚目にあるような予測を、冒頭で紹介したモデル式の各要素に分解して見られるようにしたものです。上から順にトレンド、週周期、年周期を表しています。prohet_plot_components()で簡単にプロットできます。たしかに土日や年末年始あたりで値が大きくなっている感じがします。
  
  
ちなみに、経験的に知られているシーズナリティがあれば、それをそのままモデルに組み込んでしまっても良いのですが、一応、実データをみながら本当に週周期や年周期がありそうなのかを確認してみるのも良いと思います。
その場合は、以下のようなRコードで曜日の影響をわかりやすくしたプロットを描くと良いかもしれません。
# 日本語文字化けが嫌なので変えておきます
Sys.setlocale("LC_TIME", "en_US")
# dfに曜日を追加
daily_data_weekday <- daily_data
daily_data_weekday$weekdays <- weekdays(as.POSIXct(daily_data$ds))
# 曜日で色分けしてプロット
p <- ggplot(daily_data_weekday, aes(as.POSIXct(ds),y, group=1))
p <- p + geom_point(aes(colour=weekdays))
p <- p + theme(axis.title.x = element_blank())
p
色分けの仕方を特に指定しなかったので少し見にくいですが、「緑っぽい色の土曜や日曜あたりに登録者数が多くなっているというパターンがありそうだ」と言えそうです。
また、季節や月ごとにも若干変わっていそうなので、今回は上記のモデルで問題なさそうですね(約2年分のデータしか用意していないので、なんとも言えませんが…)。
3. 新機能リリースや競合の動向などの変化点をモデルに組み込んでみる
さて、次は変化点を手動で設定して、トレンドを表すg(t)を変化させてみましょう。
Prophetのペーパーにもあるように、「成長率は常に一定である」と考えることはほとんどの場合、不自然なことです。たとえば、ある時点でサービスの新機能などがリリースされたり、競合他社の動向に変化があれば、その影響で成長率自体が変化することはあり得ることです。
  
Prophetでは、この成長率の変化も簡単にモデルに組み込むことができます。ですが、シーズナリティの組み込みなどとは違い、テキトーにやると結構エラーを吐かれてしまうので、注意しましょう。先ほどのseasonalityを調整したモデルなどにchangepoint=で変化点を手動で設定することができます。では、3/5に競合他社のサービスで新機能のリリースがあったとして、先ほどまでのモデルにこの変化点も追加してみます。
# この場合もここだけ変わっています。
m <- prophet(daily_data, yearly.seasonality = TRUE, weekly.seasonality = TRUE, changepoints = c('2017-03-05'))
future_df <- make_future_dataframe(m,365)
pred <- predict(m,future_df)
plot(m,pred) + ggtitle("Prophet Forecast Plot")
prohet_plot_components(m, pred)
雰囲気が変わりましたね。先ほどまでよりもトレンドを表す関数の勢いがすごい感じになり、かつマイナス傾向になってしまいました。この変化はどういったものなのか少し考えてみます。
  
Prophetでは、変化点に関しても手動で設定せずともデフォルトで自動検知してくれるのですが、この自動検知は全データの前半80%から候補点が選ばれるようになっているようです。
対して、ここでは明示的に2017年の3月5日のみを指定したので、prohet_plot_componentsの関数g(t)の変化も一時点でのみ発生していて、直線と直線を繋げたようなグラフになっていることがわかると思います。ここでマイナスの方向に変化したトレンドがひたすら続いている、という状況になっているわけですね(ちなみにデフォルトではこの少しあとに上昇のトレンド変化が検知されているようです)。
  
あらためて、変化点が複数選ばれているデフォルトでのg(t)と変化点を明示的に1つだけ指定した場合のg(t)を並べてみます。
変化点の数によって、かなり雰囲気が違いますね。
実際の時系列データにおいては、変化点が1つのみと考えるのも微妙なのですが、ここではわかりやすさのために、あえて変化点として2017/3/5という一時点のみを指定してみました。
実務においては、この辺りは実際のデータを見つつ、ドメイン知識を駆使しながら調整していけば良さそうです。あるいはm$changepointsm$paramsなどで、Prophetが自動検知した変化点やその変化の量を見ることができるので、それらを見ながら調整していくというのも良いと思います。
  
ところで、新機能のリリースや競合の動向とはいえ、「それほど大きな効果は期待されないだろうな」という場合や「もっと大きな変化がありそうだ」という場合には、この変化点の影響度を調整したくなると思います。
Prophetでは、changepoint.prior.scaleで変化の度合いを調整できます(デフォルトは0.05)。
m <- prophet(daily_data, yearly.seasonality = TRUE, weekly.seasonality = TRUE, changepoint.prior.scale = 0.1)
試しに、デフォルトの予測モデルに対して、変化点の影響を少し大きくしてみました。先ほど見たようにデフォルトで検知されたトレンドの変化点は2017年3月後半頃を境に上昇していますから、この変化の影響が大きくなり、デフォルトのモデルよりも未来の会員登録者数を多く予測するモデルになりました。
  
サンプルで生成したデータなので、どうなったら正解かというのはなんとも言い難いのですが、今回はこの辺りでお終いにしておこうと思います。
(プロットをみる限り、2016年の年末あたりに外れ値っぽい値が発生しているので、これをイベント効果として扱うのか、単に異常値として扱うのかでさらにモデルが変わってきたりはしますが、今回はここで力尽きました…)
まとめ
以上、簡単にではありますが、難しい理論などは一旦抜きにして、ProphetでWeb系サービスやアプリの実データっぽいものをシンプルに時系列予測してみました。
  
あらためて流れをおさらいしつつ、まとめると
  • まず、実際のデータを確認する(この記事中では周期性をみるところで曜日などを可視化したプロットをしましたが、基本的に分析のスタート時にこのプロットで確認するくらいでも良いかもしれません)
  • 次に、デフォルトのprophet()でとりあえずモデルに突っ込んでみる
  • 返ってきた結果や始めに作成したプロットを眺めつつ、周期性をモデルに取り入れていく(yearly.seasonality, weekly.seasonality)
  • トレンドを変化させるような新機能のリリースなどがあった日付があれば、モデルに組み込む(changepoints)
  • 変化量を確認したり、changepoint.prior.scaleで変化の度合いを調整する
こんな感じでしょうか。
ちなみに今回は触れられませんでしたが、Prophetで時系列予測をやっていく上では、そのうちぶつかることになりそうな壁はいくつかあります。一部だけ簡単に紹介しておくと、
  • 実務で予測をやる場合は、絶対にモデルの評価をした方が良い
絶対とか言うならここにも書けよって話なんですが、結構大変になりそうだったので、それはまた機会があればいずれ、ということで…
  • 変化点の手動指定とかでエラーが出たときに対処するために裏側でのstanの挙動を理解する
これはまあ、そんなに気にしなくても良いのかもしれません。
  
などなどです。
ですがまあ、そういう難しいことは使いながら後から考えれば良いと思うので、まずはシンプルなモデルをたくさん作りながら徐々にProphetを理解して、使いこなしていきましょう!
正直なところ、僕自身もそんな感じでまだまだ勉強中です…(笑)

なので、間違っているところやおかしなところがあったら、むしろご指摘いただければ幸いです。
  
最後は少し話が逸れてしまいましたが、ここまで長々とお読みいただき、ありがとうございました!
明日は、eurekaのスゴい人集団SREチームのdatchこと原田くんです!乞うご期待!
One clap, two clap, three clap, forty?

By clapping more or less, you can signal to us which stories really stand out.