読者です 読者をやめる 読者になる 読者になる

hsuetsugu’s diary

ITの技術的なことに関して主に書きます。Rとpythonとd3.jsとAWSとRaspberryPiあたりを不自由なく使いこなせるようになりたいです。

(35連休)9日目:ARIMAモデルとdlm(動的線形モデル)による販売量予測モデルの検証②

前回(http://hsuetsugu.hatenablog.com/entry/2014/08/22/164140)の続きです。

前回のおさらいと予測誤差(RMSE)の算出

前回のおさらい

前回、RのBJsalesという150日間の日次販売量データを用いて、ARIMAとdlmにて、直近98日間(14週間)のデータでの翌週7日間の日次販売量予測を7週間ローリングして予測値を算出しました。その結果が下記です。
f:id:hsuetsugu:20140822163144p:plain
上記の予測誤差を計算していきます。

予測誤差指標について

一般的に使われる予測誤差の指標はいろいろありますが、在庫設定に反映しやすい、という点で、私は二乗平均平方根誤差(RMSE:Root Mean Squared Error)を用いることが多いです。これはこの名前の通り、誤差系列の2乗の平均をとったもので、「標準偏差」の計算式と似ています。
標準偏差は実績値の平均値に対する差分の二乗平均であるので、予測がどれだけ意味のあるものかを計る際に、実績の標準偏差とと予測値のRMSEを比較することで評価することができます。

RMSEの計算結果

下記の通りで、ARIMAもdlmも単純な平均値よりはもちろん優れており(=予測の意味がある)、今回のケースではARIMA < dlm という結果でした。

  • 実績のボラティリティ : 4.48
  • ARIMAでの予測結果のRMSE : 3.41
  • dlmでの予測結果のRMSE : 2.82

ちなみに、RMSEは誤差時系列の標準偏差とだいたい同じになる(誤差時系列の平均が0であれば完全に同じになる)ので、簡易的に誤差系列の標準偏差を求めることも多いです。

誤差を求めるコードは以下の通りです(前回のコードの続きと思ってください)。

> # 誤差評価(RMSEの算出)
> actual.sd <- sd(data[99:147]) #実績のボラティリリィ
> print(actual.sd)
[1] 4.478719
> 
> er.arima <- data[99:147] - frcst.arima[99:147] # 誤差時系列(arima)
> er.dlm <- data[99:147] - frcst.dlm[99:147]  # 誤差時系列(dlm)
> rmse.arima <- mean((er.arima)^2)^0.5 # RMSE(arima)
> rmse.dlm <- mean((er.dlm)^2)^0.5  # RMSE(dlm)
> print(rmse.arima)
[1] 3.417305
> print(rmse.dlm)
[1] 2.823517
> 
> er.arima.sd <- sd(er.arima)
> er.dlm.sd <- sd(er.dlm)
> print(er.arima.sd)
[1] 3.363156
> print(er.dlm.sd)
[1] 2.640505
> 
> er.arima.sdovermyu <- sd(er.arima) / mean(data[99:147])
> er.dlm.sdovermyu <- sd(er.dlm) / mean(data[99:147])
> print(er.arima.sdovermyu)
[1] 0.01312739
> print(er.dlm.sdovermyu)
[1] 0.01030667
> 
> # セッションインフォ  
> sessionInfo()
R version 3.1.1 (2014-07-10)
Platform: x86_64-apple-darwin10.8.0 (64-bit)

locale:
[1] ja_JP.UTF-8/ja_JP.UTF-8/ja_JP.UTF-8/C/ja_JP.UTF-8/ja_JP.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] dlm_1.1-3        forecast_5.5     timeDate_3010.98 zoo_1.7-11      

loaded via a namespace (and not attached):
 [1] colorspace_1.2-4 fracdiff_1.4-2   grid_3.1.1       lattice_0.20-29  nnet_7.3-8      
 [6] parallel_3.1.1   quadprog_1.5-5   Rcpp_0.11.2      tools_3.1.1      tseries_0.10-32

欠損値のある場合の予測結果

では、需要側の問題ではなく、物流に問題があってたまたま店舗欠品を過去にしていた場合、予測値にどう影響がでるでしょうか。
意図的に、50日目と51日目の売上を0として、全く同じスクリプトをまわした結果が下記です。
f:id:hsuetsugu:20140823124711p:plain

・・・ぶれますね、ここまでぶれますかってくらいぶれてしまいました。なんかdlmについては予測モデルの書き方が違うような気がします。元データを補正するしかないのか、書き方が悪いのか・・・。
わたしはこれまでは元データを補正して、異常値をのぞいてから予測していたのですが、dlmならうまいこと「状態」を推定してやってくれるのかなとか思ってましたがそんな甘くはないのかもしれません。
自分の理解が足りていないだけかもしれないので、引き続き勉強します。

※2014.8.25訂正
上記、計算方法に不適切な部分がありました。下記が訂正記事です。
http://hsuetsugu.hatenablog.com/entry/2014/08/25/152818