K_DM Book
用語集
対象者 標準
対象者別の入口 →

対象者の設定はこのブラウザの中にだけ保存されます。

言語 JA

Prophetパラメーター

中級

04.15.3

Prophetパラメーター

最終更新 2020-02-26 読了時間 6 分
まとめ
Prophet のトレンド \(g(t)\) は、変化点 \(s_j\) ごとに傾きの変化量 \(\delta_j\) を持つ区分線形関数で、changepoint_prior_scale（既定 0.05）が \(\delta_j\) の大きさを抑える強さを決めます。このページでは、曜日の cos 波・月の偶奇の段差・対数トレンドに加えて 183 日目で水準が 0.5 跳ねる 365 日の合成データを作り、growth="linear" の Prophet で 90 日先を予測して、変化点がその跳ねをどう吸収するかを見ます。そのうえで Prophet の Stan コードを Python に写した linear_trend(k, m, delta, t, t_change) を定義し、変化点行列 \(A\) と \(\delta\) から折れ線が作られる仕組みを、ランダムな傾きで引いた線で可視化します。読み終えると、changepoint_prior_scale を上げ下げしたときに予測線がしなやかになったり硬くなったりする理由を式から説明できるようになります。

Prophetを実行してみる で fit → predict → plot の流れと、Prophetのモデルの中身 で growth="linear" と "flat" の違いを先に見ておくと、このページの linear_trend が何を再現しているのかがすぐ分かります。

Prophet にどのようなパラメーターがあるかを、トレンド \(g(t)\) を中心に整理します。

このページは書きかけです

直感 #

Prophet のトレンドは「何本かの直線をつないだ折れ線」です。学習期間の前 80% に 25 個の変化点候補 \(s_1, \dots, s_{25}\) を等間隔に置き、各候補で傾きがどれだけ変わるか（\(\delta_j\)）を推定します。\(\delta_j\) にはラプラス分布の事前分布がかかり、その幅が changepoint_prior_scale です。0.05 のような小さい値ならほとんどの \(\delta_j\) が 0 に押しつぶされて折れ線は直線に近づき、0.5 なら多くの候補で傾きが変わってデータの細かい揺れまで追います。このページのデータには 183 日目に 0.5 の跳ねを仕込んであるので、変化点がその付近で傾きを変えて吸収する様子を予測図で確認できます。

1
2
3
4
5

import numpy as np
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
from prophet import Prophet

実験に使用するデータ #

前のページと同じ構造（曜日の cos 波、偶数月の段差、対数トレンド、ノイズ）に、(i > 182) * 0.5 で 183 日目以降の水準を 0.5 上げる構造変化を足した 365 日のデータです。各成分の係数は前のページより小さく（cos は 1/3、段差は 1/5）してあるので、0.5 の跳ねが目立ちます。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26

date = pd.date_range("2020-01-01", periods=365, freq="D")

# 予測対象
y = [
    np.cos(di.weekday()) / 3
    + di.month % 2 / 5
    + np.log(i + 1) / 5.0
    + (i > 182) * 0.5
    + np.random.rand() / 10
    for i, di in enumerate(date)
]

# トレンド成分
x = [18627 + i - 364 for i in range(365)]
trend_y = [np.log(i + 1) / 3.0 for i, di in enumerate(date)]
weekly_y = [np.cos(di.weekday()) for i, di in enumerate(date)]
seasonal_y = [di.month % 2 / 2 for i, di in enumerate(date)]
noise_y = [np.random.rand() / 10 for i in range(365)]

df = pd.DataFrame({"ds": date, "y": y})
df.index = date

# 実験に使用するデータ
plt.title("サンプルデータ")
sns.lineplot(data=df)
plt.show()

実験に使用するデータの図

何が起きているか（ステップで見る） #

コードを読む前に、この手法が図の上で何をしているかを追ってみましょう。「再生」を押すか、スライダーを動かすと 1 手順ずつ進みます。

見るポイント：水準 → 傾き → 季節 と要素を足すほど、未来への伸ばし方が豊かになる。

growthパラメータ #

growth=“linear"とした場合の予測 #

週次季節性だけを有効にし、growth="linear" で 90 日先を予測します。changepoint_prior_scale は既定の 0.05 のままなので、変化点は控えめにしか使われません。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16

# モデルを訓練
m = Prophet(
    yearly_seasonality=False,
    weekly_seasonality=True,
    daily_seasonality=False,
    growth="linear",
)
m.fit(df)

# 将来を予測
future = m.make_future_dataframe(periods=90)
forecast = m.predict(future)
fig = m.plot(forecast)
plt.axvspan(18627, 18627 + 90, color="coral", alpha=0.4, label="予測期間")
plt.legend()
plt.show()

Initial log joint probability = -16.6231
    Iter      log prob        ||dx||      ||grad||       alpha      alpha0  # evals  Notes 
      81       792.152    0.00138864       97.2266   2.443e-05       0.001      145  LS failed, Hessian reset 
      99       792.228   7.40168e-07       55.8404     0.03797     0.03797      168   
    Iter      log prob        ||dx||      ||grad||       alpha      alpha0  # evals  Notes 
     138       792.388   3.51776e-05       59.7965    6.64e-07       0.001      263  LS failed, Hessian reset 
     199        792.62    3.0252e-05       72.8985      0.4243           1      340   
    Iter      log prob        ||dx||      ||grad||       alpha      alpha0  # evals  Notes 
     247       792.799   0.000165183       55.2724   3.174e-06       0.001      440  LS failed, Hessian reset 
     299       792.972   0.000162358       62.7075           1           1      508   
    Iter      log prob        ||dx||      ||grad||       alpha      alpha0  # evals  Notes 
     347       792.991   3.19534e-05       38.7196   4.536e-07       0.001      634  LS failed, Hessian reset 
     392       792.994   7.68689e-07       61.3331   1.015e-08       0.001      733  LS failed, Hessian reset 
     399       792.994     9.628e-08       55.3164      0.2134           1      742   
    Iter      log prob        ||dx||      ||grad||       alpha      alpha0  # evals  Notes 
     404       792.994   3.74603e-08       39.3307      0.2405      0.7727      749   
Optimization terminated normally: 
  Convergence detected: relative gradient magnitude is below tolerance

Convergence detected: relative gradient magnitude is below t…の図

実行結果の読み方 #

青い予測線は 7 月前後で傾きを急にして 0.5 の跳ねを追いかけ、その後また緩やかな傾きに戻る折れ線になります。瞬間的な 0.5 の段差を直線の組み合わせで表すため、跳ねの直前で予測が低すぎ、直後で高すぎる数週間が生じます。changepoint_prior_scale を 0.5 に上げると傾きの変化が大きく許されて段差に近い形になり、0.01 に下げるとほぼ 1 本の直線になって跳ねを無視します。coral の予測期間では、最後の変化点以降の傾きで直線的に伸び、帯は変化点が未来にも起こると仮定して広がります。

growth=“linear"で表現できるトレンド成分 #

growth="linear" でどのようなトレンドを表現できるかを、Prophet が growth="linear" のときにトレンドの計算に使っている Stan のコードを Python で再現して確かめます。

※参考文献のコードと記号を合わせました。ミスがある場合はお手数ですが issue に指摘いただけると大変助かります。

prophet/R/inst/stan/prophet.stan

行列 \(A\) の各行の次元数は変化点の数です。また、各次元の値は \( a_j(t) \left{ \begin{array}{ll} 1, & \text{if}　\ t \ge s_j \ 0, & \text{otherwise} \end{array} \right. \) となっています。時刻 \(t\) の \(a\) の \(j\) 次元目は、\(t\) が \(j\) 番目の変化点 \(s_j\) よりあとの時刻ならば 1 になる、ということです。つまり \(A\) の行は「時刻 \(t\) までに通過した変化点」を 0/1 で並べたもので、\(A \delta\) の和が「それまでに積み上がった傾きの変化の合計」になります。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19

def linear_trend(
    k: float, m: float, delta: np.array, t: np.array, t_change: list
) -> np.array:
    """線形トレンドを作成する

    Args:
        k (float): 係数
        m (float): 係数
        delta (np.array): delta
        t (np.array): タイムスタンプ
        t_change (list): トレンド変化点

    Returns:
        np.array: トレンド
    """
    A = np.vstack([np.where(t_change < ti, 1, 0) for ti in t])
    return (k + (A * delta).sum(axis=1)) * t + (
        m + (A * (-t_change * delta)).sum(axis=1)
    )

linear_trendで表現できるトレンドを可視化する #

変化点を 6 個等間隔に置き、各区間の傾きの変化量 \(\delta_j\) を −0.5〜0.5 の一様乱数で与えて、linear_trend がどんな線を引くかを描きます。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

# 適当にトレンドのような線を作成するため、変化点の数と各区間での線の傾きを指定します
change_point_num = 5
t = np.linspace(-10, 10, 100)
t_change = np.linspace(-10, 10, change_point_num + 1)
delta = np.array([np.random.rand() - 0.5 for _ in range(change_point_num + 1)])

# linear_trendで線を引く
trend_y = linear_trend(0.1, 0, delta, t, t_change)

# できた線をプロットしてみる
[plt.axvline(x=tci, color="red", alpha=0.5) for tci in t_change]  # 変化点をプロット
plt.plot(t, trend_y)
plt.show()

linear_trendで表現できるトレンドを可視化するの図

実行結果の読み方 #

折れ線は 6 本の赤い縦線（変化点）の位置でだけ傾きが変わり、それ以外の区間は直線です。これが Prophet の growth="linear" が表現できるトレンドのすべてで、滑らかな曲線や瞬間的な段差は、変化点を密に置いて直線で近似するしかありません。前の予測図で 183 日目の跳ねが「急な傾きの区間」として現れたのは、この制約のためです。changepoint_prior_scale は \(\delta_j\) の大きさの上限のように働くので、上げれば折れ方が激しくデータ追従的に、下げれば直線的で頑健になります。

使いどころと落とし穴 #

changepoint_prior_scale は既定の 0.05 から、トレンドが鈍すぎるなら 0.1〜0.5 に上げ、予測が暴れるなら 0.01〜0.03 に下げる、という調整が基本です。変化点候補は既定で学習期間の前 80%（changepoint_range=0.8）にしか置かれないので、この例のように跳ねが学習期間の後半 20% にあると吸収されず、changepoint_range=0.95 に広げるか、changepoints に日付を直接渡す必要があります。また予測期間では最後の傾きがそのまま延長されるため、学習末尾でたまたま急な傾きになっていると 90 日先で大きく外れます。長期の予測ほど changepoint_prior_scale を小さめにし、uncertainty_samples で帯の広がりを確認して、どこまで先を信じるかを決めてください。