TimeSeries_CVLAB
Takuma Mimura
4/27/2020
This file is going to cover the basics of cross validation method to time series prediction.
#install.packages("forecast")
#install.packages("ggfortify")
#install.packages("ggplot2")
#install.packages("gridExtra")
library(forecast)## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary(ggfortify)## Loading required package: ggplot2## Registered S3 methods overwritten by 'ggfortify':
## method from
## autoplot.Arima forecast
## autoplot.acf forecast
## autoplot.ar forecast
## autoplot.bats forecast
## autoplot.decomposed.ts forecast
## autoplot.ets forecast
## autoplot.forecast forecast
## autoplot.stl forecast
## autoplot.ts forecast
## fitted.ar forecast
## fortify.ts forecast
## residuals.ar forecastlibrary(ggplot2)
library(gridExtra)
# 対象データ
# Monthly Airline Passenger Numbers 1949-1960
# 対数変換したものを使う
log_air_pass <- log(AirPassengers)
class(log_air_pass)## [1] "ts"#[1] "ts"
ggtsdisplay(log_air_pass)
print(log_air_pass, digits = 4)## Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
## 1949 4.718 4.771 4.883 4.860 4.796 4.905 4.997 4.997 4.913 4.779 4.644 4.771
## 1950 4.745 4.836 4.949 4.905 4.828 5.004 5.136 5.136 5.063 4.890 4.736 4.942
## 1951 4.977 5.011 5.182 5.094 5.147 5.182 5.293 5.293 5.215 5.088 4.984 5.112
## 1952 5.142 5.193 5.263 5.198 5.209 5.384 5.438 5.489 5.342 5.252 5.147 5.268
## 1953 5.278 5.278 5.464 5.460 5.434 5.493 5.576 5.606 5.468 5.352 5.193 5.303
## 1954 5.318 5.236 5.460 5.425 5.455 5.576 5.710 5.680 5.557 5.434 5.313 5.434
## 1955 5.489 5.451 5.587 5.595 5.598 5.753 5.897 5.849 5.743 5.613 5.468 5.628
## 1956 5.649 5.624 5.759 5.746 5.762 5.924 6.023 6.004 5.872 5.724 5.602 5.724
## 1957 5.753 5.707 5.875 5.852 5.872 6.045 6.142 6.146 6.001 5.849 5.720 5.817
## 1958 5.829 5.762 5.892 5.852 5.894 6.075 6.196 6.225 6.001 5.883 5.737 5.820
## 1959 5.886 5.835 6.006 5.981 6.040 6.157 6.306 6.326 6.138 6.009 5.892 6.004
## 1960 6.033 5.969 6.038 6.133 6.157 6.282 6.433 6.407 6.230 6.133 5.966 6.068# サンプルサイズ
length(log_air_pass)## [1] 144# [1] 144
# 何年？
length(log_air_pass) / 12## [1] 12# [1] 12
# test and train
# first 10 years = train、last 2 year = test
train <- window(log_air_pass, end=c(1958, 12))
test <- window(log_air_pass, start=c(1959, 1))Model
mod_sarima <- auto.arima(train,
ic = "aic",
stepwise = F,
parallel = TRUE,
num.cores = 4)
mod_sarima## Series: train
## ARIMA(0,1,1)(0,1,1)[12]
##
## Coefficients:
## ma1 sma1
## -0.3424 -0.5405
## s.e. 0.1009 0.0877
##
## sigma^2 estimated as 0.001432: log likelihood=197.51
## AIC=-389.02 AICc=-388.78 BIC=-381#Series: train
#ARIMA(0,1,1)(0,1,1)[12]
#Coefficients:
# ma1 sma1
# -0.3424 -0.5405
#s.e. 0.1009 0.0877
#sigma^2 estimated as 0.001432: log likelihood=197.51
#AIC=-389.02 AICc=-388.78 BIC=-381
checkresiduals(mod_sarima)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,1)(0,1,1)[12]
## Q* = 20.34, df = 22, p-value = 0.5618
##
## Model df: 2. Total lags used: 24#Ljung-Box test
#data: Residuals from ARIMA(0,1,1)(0,1,1)[12]
#Q* = 20.34, df = 22, p-value = 0.5618
#Model df: 2. Total lags used: 24Forecast function
# 予測期間
h <- length(test)
# 予測
f_sarima <- forecast(mod_sarima, h = h)
# 予測結果の図示
autoplot(f_sarima, main = "Mimura prediction")
Test data
accuracy(f_sarima, test)## ME RMSE MAE MPE MAPE
## Training set -0.0001250111 0.03539809 0.02640383 -0.0003836216 0.4867206
## Test set 0.0895992058 0.09593236 0.08959921 1.4634770472 1.4634770
## MASE ACF1 Theil's U
## Training set 0.2145046 0.005834845 NA
## Test set 0.7279033 0.327959306 0.8735145# ME RMSE MAE MPE MAPE MASE
#Training set -0.0001250111 0.03539809 0.02640383 -0.0003836216 0.4867206 0.2145046
#Test set 0.0895992058 0.09593236 0.08959921 1.4634770472 1.4634770 0.7279033
# ACF1 Theil's U
#Training set 0.005834845 NA
#Test set 0.327959306 0.8735145
## 特別な技術を要しない単純な予測手法
# 訓練データの最終時点を予測値とする
f_rwf <- rwf(train, h = h)
# 全体のトレンドを組み込む
f_rwf_drift <- rwf(train, drift = TRUE, h = h)
# 1周期前(今回は12時点前)を予測値とする
f_snaive <- snaive(train, h = h)
p1 <- autoplot(f_rwf, main = "Prediction1")
p2 <- autoplot(f_rwf_drift, main = "Prediction(Naive)")
p3 <- autoplot(f_snaive, main = "Prediction(seasonal)")
grid.arrange(p1, p2, p3, nrow=3)
1．訓練データの最終時点を予測値とする 2．1にたいして、全体のトレンド(今回のデータでは右肩上がり)を組み込む 3．1周期前(今回は12時点前)を予測値とする
1．は平均に回帰しない、非定常な時系列データに対して用いられるナイーブ予測です。 2．は、さらに右肩上がりのトレンドを組みこみました。 3．は、周期性を持つデータに対してしばしば用いられる方法で、前年同期の値を予測値としてそのまま使っている
# 各々の予測結果のRMSEを比較する
accuracy(f_sarima, test)[, "RMSE"]## Training set Test set
## 0.03539809 0.09593236accuracy(f_rwf, test)[, "RMSE"]## Training set Test set
## 0.1058089 0.3230635accuracy(f_rwf_drift, test)[, "RMSE"]## Training set Test set
## 0.1054032 0.2204105accuracy(f_snaive, test)[, "RMSE"]## Training set Test set
## 0.1381999 0.1821230## 図示するために、整形をする
# 評価セットの凡例
types <- c("train", "test")
types## [1] "train" "test"#[1] "train" "test"
acc_df_sarima <- data.frame(forecast = rep("f_sarima", 2),
type = types,
RMSE = accuracy(f_sarima, test)[, "RMSE"],
row.names = NULL)
acc_df_sarima## forecast type RMSE
## 1 f_sarima train 0.03539809
## 2 f_sarima test 0.09593236# データを簡単に整形できる便利関数を自作する
make_acc_df <- function(forecast_name, types, rmse_vec){
# グラフにしやすい形式に、予測誤差を整形する関数
#
# Args:
# forecast_name: 予測の方法の名称
# types: 予測誤差のタイプ
# rmse_vec: RMSEが格納されたベクトル
#
# Returns:
# 整形された予測誤差のデータフレーム
vec_length <- length(rmse_vec)
acc_df <- data.frame(forecast = rep(forecast_name, vec_length),
type = types,
RMSE = rmse_vec,
row.names = NULL)
return(acc_df)
}
make_acc_df("f_sarima", types, accuracy(f_sarima, test)[, "RMSE"])## forecast type RMSE
## 1 f_sarima train 0.03539809
## 2 f_sarima test 0.09593236acc_df <- rbind(
make_acc_df("f_sarima", types, accuracy(f_sarima, test)[, "RMSE"]),
make_acc_df("f_rwf", types, accuracy(f_rwf, test)[, "RMSE"]),
make_acc_df("f_rwf_drift", types, accuracy(f_rwf_drift, test)[, "RMSE"]),
make_acc_df("f_snaive", types, accuracy(f_snaive, test)[, "RMSE"])
)
ggplot(acc_df, aes(x = forecast, y = RMSE, fill = type)) +
geom_bar(stat = "identity", position = "dodge") +
ggtitle("Comparing RMSE and model predictions")
Evaluating crossvalidation prediction methods.
mod_sarima## Series: train
## ARIMA(0,1,1)(0,1,1)[12]
##
## Coefficients:
## ma1 sma1
## -0.3424 -0.5405
## s.e. 0.1009 0.0877
##
## sigma^2 estimated as 0.001432: log likelihood=197.51
## AIC=-389.02 AICc=-388.78 BIC=-381f_sarimax_func <- function(ts_data, h){
# 訓練データと予測期間を指定すると、予測結果が出力される関数
#
# Args:
# ts_data: training data
# h:予測期間
#
# Returns:
# forecast関数による予測の結果
# 選ばれたモデルと同じ次数のARIMAモデルを推定
mod <- Arima(ts_data, order = c(0,1,1), seasonal = c(0,1,1))
# 予測結果を出力
return(forecast(mod, h = h))
}
model_conf <- f_sarimax_func(train, h = h)
head(model_conf$mean)## Jan Feb Mar Apr May Jun
## 1959 5.853880 5.803858 5.948464 5.920534 5.947908 6.116194head(f_sarima$mean)## Jan Feb Mar Apr May Jun
## 1959 5.853880 5.803858 5.948464 5.920534 5.947908 6.116194# クロスバリデーション法の実行
# h=3で実行
h <- 3
sarima_cv <- tsCV(log_air_pass, f_sarimax_func, h = h)
rwf_cv <- tsCV(log_air_pass, rwf, h = h)
rwf_drift_cv <- tsCV(log_air_pass, rwf, h = h, drift = TRUE)
snaive_cv <- tsCV(log_air_pass, snaive, h = h)
sarima_cv## h=1 h=2 h=3
## Jan 1949 NA NA NA
## Feb 1949 NA NA NA
## Mar 1949 NA NA NA
## Apr 1949 NA NA NA
## May 1949 NA NA NA
## Jun 1949 NA NA NA
## Jul 1949 NA NA NA
## Aug 1949 NA NA NA
## Sep 1949 NA NA NA
## Oct 1949 NA NA NA
## Nov 1949 NA NA NA
## Dec 1949 NA NA NA
## Jan 1950 NA NA NA
## Feb 1950 NA NA NA
## Mar 1950 -2.049961e-02 -0.0334477724 0.0326897899
## Apr 1950 -7.194391e-03 0.0589540035 0.0988687074
## May 1950 7.191144e-02 0.1118407031 0.1118574929
## Jun 1950 8.280028e-02 0.0828002756 0.0941542596
## Jul 1950 -3.139856e-03 0.0081735779 -0.0305860615
## Aug 1950 1.255202e-02 -0.0262045463 -0.0456688468
## Sep 1950 -4.528734e-02 -0.0648050747 0.0142346657
## Oct 1950 1.211607e-03 0.0803617959 0.1412139429
## Nov 1950 7.860547e-02 0.1394559241 0.0820226564
## Dec 1950 1.971271e-02 -0.0213144066 0.0374963931
## Jan 1951 -6.415114e-02 -0.0055931463 -0.0440372741
## Feb 1951 1.108546e-01 0.0660927861 0.1966632789
## Mar 1951 -5.000690e-02 0.0742281416 -0.0340960817
## Apr 1951 1.159393e-01 -0.0002373525 -0.0053838325
## May 1951 -8.529862e-02 -0.0944718471 -0.0944733337
## Jun 1951 -2.687042e-02 -0.0268752053 -0.0263726136
## Jul 1951 -1.311199e-02 -0.0126062597 0.0129367684
## Aug 1951 -5.771486e-03 0.0197715255 0.0602165488
## Sep 1951 2.277774e-02 0.0632226508 0.0257284893
## Oct 1951 5.153160e-02 0.0140375876 0.0355168110
## Nov 1951 -1.301198e-02 0.0086216327 0.0007684895
## Dec 1951 1.547489e-02 0.0076217335 -0.0545606535
## Jan 1952 -9.062326e-04 -0.0630886785 -0.0757806656
## Feb 1952 -6.261736e-02 -0.0753094469 -0.0352428182
## Mar 1952 -4.379977e-02 -0.0037330648 0.0648049178
## Apr 1952 2.170443e-02 0.0902419669 0.0320526178
## May 1952 7.796931e-02 0.0197803092 0.0706353172
## Jun 1952 -1.922473e-02 0.0316305841 -0.0362664612
## Jul 1952 4.159023e-02 -0.0263064833 0.0280013626
## Aug 1952 -4.657598e-02 0.0077316883 0.0339062667
## Sep 1952 2.963837e-02 0.0558135471 0.0228013716
## Oct 1952 4.266941e-02 0.0096573432 0.0055358655
## Nov 1952 -9.828902e-03 -0.0138699376 -0.0710528734
## Dec 1952 -9.317697e-03 -0.0665007373 0.0028465213
## Jan 1953 -6.212108e-02 0.0072261943 0.0576542816
## Feb 1953 3.536913e-02 0.0857973766 0.0789946670
## Mar 1953 6.994106e-02 0.0631388433 -0.0011266310
## Apr 1953 3.390825e-02 -0.0303572090 -0.0446947180
## May 1953 -4.672347e-02 -0.0610611678 -0.0439247325
## Jun 1953 -4.028690e-02 -0.0231501016 -0.0652104330
## Jul 1953 -6.816650e-03 -0.0488769729 -0.0342865805
## Aug 1953 -4.602346e-02 -0.0314330852 -0.0659235792
## Sep 1953 -1.146907e-02 -0.0459597799 -0.0807342375
## Oct 1953 -4.063591e-02 -0.0754104967 -0.0736802406
## Nov 1953 -5.614500e-02 -0.0542488131 -0.1816732836
## Dec 1953 -2.437729e-02 -0.1518018428 -0.0588978005
## Jan 1954 -1.375324e-01 -0.0446284251 -0.0346754871
## Feb 1954 3.806607e-02 0.0432214327 0.0877124815
## Mar 1954 1.814518e-02 0.0623957114 0.0803591495
## Apr 1954 5.171316e-02 0.0700169964 0.1172469628
## May 1954 4.152157e-02 0.0899594849 0.0349699412
## Jun 1954 6.423450e-02 0.0097983057 0.0047910062
## Jul 1954 -3.079117e-02 -0.0379175176 -0.0404207124
## Aug 1954 -1.928648e-02 -0.0212890557 -0.0085574795
## Sep 1954 -9.837372e-03 0.0030842947 -0.0058602384
## Oct 1954 8.956478e-03 -0.0003142929 0.0405420220
## Nov 1954 -6.050416e-03 0.0348685525 0.0007614486
## Dec 1954 3.850374e-02 0.0043455784 -0.0258973169
## Jan 1955 -1.989163e-02 -0.0493343318 -0.0033823646
## Feb 1955 -3.890192e-02 0.0062867980 0.0079697405
## Mar 1955 3.032337e-02 0.0332325757 0.0780532130
## Apr 1955 1.445654e-02 0.0594151604 0.0998501442
## May 1955 5.077295e-02 0.0910183326 0.0391241151
## Jun 1955 6.037290e-02 0.0088190015 0.0220282069
## Jul 1955 -2.932493e-02 -0.0152410067 -0.0243368348
## Aug 1955 3.411450e-03 -0.0057598576 -0.0220116454
## Sep 1955 -7.836480e-03 -0.0240804618 0.0105964239
## Oct 1955 -1.933558e-02 0.0152610915 0.0062600393
## Nov 1955 2.690699e-02 0.0180031877 0.0143556629
## Dec 1955 2.099094e-03 -0.0027886202 -0.0257959886
## Jan 1956 -4.269966e-03 -0.0272076557 -0.0191338180
## Feb 1956 -2.462559e-02 -0.0166561851 -0.0040262509
## Mar 1956 -1.878675e-03 0.0107387148 0.0471799220
## Apr 1956 1.206355e-02 0.0486290680 0.0284198476
## May 1956 4.113066e-02 0.0206944563 0.0167711436
## Jun 1956 -5.678157e-03 -0.0081621165 -0.0248952596
## Jul 1956 -4.921448e-03 -0.0217117590 -0.0459498298
## Aug 1956 -1.876803e-02 -0.0430057842 -0.0293432990
## Sep 1956 -3.174706e-02 -0.0180472451 -0.0349358096
## Oct 1956 1.190709e-03 -0.0157161027 -0.0137462507
## Nov 1956 -1.642951e-02 -0.0144299019 -0.0366109179
## Dec 1956 -4.419196e-03 -0.0267594797 -0.0082401237
## Jan 1957 -2.421573e-02 -0.0057060989 -0.0113548864
## Feb 1957 9.344043e-03 0.0032706388 0.0144819180
## Mar 1957 -2.329847e-03 0.0089762087 0.0413188304
## Apr 1957 1.046835e-02 0.0429179698 0.0276929657
## May 1957 3.642352e-02 0.0211751179 0.0434925285
## Jun 1957 -1.645139e-03 0.0214255283 -0.0015093003
## Jul 1957 2.239130e-02 -0.0005743203 -0.0184865190
## Aug 1957 -1.440420e-02 -0.0324334370 -0.0316863690
## Sep 1957 -2.357370e-02 -0.0228382912 -0.0576786990
## Oct 1957 -8.502672e-03 -0.0433061351 -0.0594772960
## Nov 1957 -3.812556e-02 -0.0542900364 -0.0880136324
## Dec 1957 -3.086907e-02 -0.0644696597 -0.0917010116
## Jan 1958 -4.511017e-02 -0.0723038609 -0.0929432809
## Feb 1958 -4.307181e-02 -0.0639172522 -0.0355979411
## Mar 1958 -3.505358e-02 -0.0067127857 0.0190821728
## Apr 1958 1.751391e-02 0.0436446557 0.0586519796
## May 1958 3.151985e-02 0.0465621083 0.0840264523
## Jun 1958 2.534698e-02 0.0628037977 -0.0281879948
## Jul 1958 4.568411e-02 -0.0452874344 -0.0210969956
## Aug 1958 -7.671706e-02 -0.0525252946 -0.0699282826
## Sep 1958 -1.177972e-03 -0.0186556009 -0.0508550588
## Oct 1958 -1.788584e-02 -0.0501035347 -0.0049848540
## Nov 1958 -3.841147e-02 0.0067093915 0.0050098251
## Dec 1958 3.222418e-02 0.0309522541 0.0578894909
## Jan 1959 9.695296e-03 0.0365552867 0.0394903515
## Feb 1959 3.006934e-02 0.0329409774 0.0647005209
## Mar 1959 1.302989e-02 0.0450930819 -0.0056166134
## Apr 1959 3.657371e-02 -0.0139689360 0.0227678876
## May 1959 -3.859335e-02 -0.0018985035 0.0105717980
## Jun 1959 2.383532e-02 0.0367652989 0.0201020950
## Jul 1959 2.107665e-02 0.0048423535 0.0077061715
## Aug 1959 -8.453780e-03 -0.0056209078 0.0142479160
## Sep 1959 -4.299025e-05 0.0198756286 0.0306671578
## Oct 1959 1.989131e-02 0.0306556329 0.0187929491
## Nov 1959 1.759474e-02 0.0058312380 -0.0084653939
## Dec 1959 -5.560612e-03 -0.0200554916 -0.1075567720
## Jan 1960 -1.645493e-02 -0.1039444384 0.0180127268
## Feb 1960 -9.325312e-02 0.0287140558 0.0115763680
## Mar 1960 9.032702e-02 0.0731918648 0.0534828936
## Apr 1960 1.614629e-02 -0.0034594058 0.0183069059
## May 1960 -1.334732e-02 0.0083055756 -0.0307554094
## Jun 1960 1.642469e-02 -0.0226284663 -0.0190726711
## Jul 1960 -3.268505e-02 -0.0289411469 0.0044228387
## Aug 1960 -9.610929e-03 0.0237984154 -0.0152183780
## Sep 1960 2.963461e-02 -0.0094006940 -0.0134263799
## Oct 1960 -2.720049e-02 -0.0312582943 NA
## Nov 1960 -1.500723e-02 NA NA
## Dec 1960 NA NA NAautoplot(sarima_cv, main = "Transitions of Sarima model errors.")
# 予測の評価からは訓練データが含まれる時点を除いたうえで、
# NULLの無い年を使う
sarima_cv_test <- na.omit(window(sarima_cv, start = c(1959, 1)))
rwf_cv_test <- na.omit(window(rwf_cv, start = c(1959, 1)))
rwf_drift_cv_test <- na.omit(window(rwf_drift_cv, start = c(1959, 1)))
snaive_cv_test <- na.omit(window(snaive_cv, start = c(1959, 1)))
# RMSEを計算する関数を作る
calc_rmse <- function(target_resid){
# 予測残差を入れるとRMSEが出力される関数
#
# Args:
# target_resid:予測残差(y - y_hat)
#
# Returns:
# RMSE
return(sqrt(mean(target_resid^2)))
}
apply(sarima_cv_test, 2, calc_rmse)## h=1 h=2 h=3
## 0.03497509 0.03579213 0.03478006# 凡例
types <- c("h=1", "h=2", "h=3")
# 整形されたデータ
acc_df_cv <- rbind(
make_acc_df("f_sarima", types, apply(sarima_cv_test, 2, calc_rmse)),
make_acc_df("rwf", types, apply(rwf_cv_test, 2, calc_rmse)),
make_acc_df("rwf_drift", types, apply(rwf_drift_cv_test, 2, calc_rmse)),
make_acc_df("snaive", types, apply(snaive_cv_test, 2, calc_rmse))
)
# 図示
ggplot(acc_df_cv, aes(x = forecast, y = RMSE, fill = type)) +
geom_bar(stat = "identity", position = "dodge") +
ggtitle("Comparing prediction errorsevaluated by Cross Validation")
sarima_cv_24 <- tsCV(log_air_pass, f_sarimax_func, h = 12)
sarima_cv_24 <- na.omit(window(sarima_cv_24, start = c(1959, 1)))
sarima_cv_24## h=1 h=2 h=3 h=4 h=5
## Jan 1959 9.695296e-03 0.036555287 0.039490352 0.071091263 0.019721211
## Feb 1959 3.006934e-02 0.032940977 0.064700521 0.013556406 0.050204492
## Mar 1959 1.302989e-02 0.045093082 -0.005616613 0.031095698 0.043646779
## Apr 1959 3.657371e-02 -0.013968936 0.022767888 0.035431202 0.019516340
## May 1959 -3.859335e-02 -0.001898504 0.010571798 -0.004811287 -0.002001832
## Jun 1959 2.383532e-02 0.036765299 0.020102095 0.022991303 0.042730875
## Jul 1959 2.107665e-02 0.004842354 0.007706171 0.027505575 0.037876836
## Aug 1959 -8.453780e-03 -0.005620908 0.014247916 0.024864722 0.013107243
## Sep 1959 -4.299025e-05 0.019875629 0.030667158 0.018785031 0.004594273
## Oct 1959 1.989131e-02 0.030655633 0.018792949 0.004585004 -0.083035399
## Nov 1959 1.759474e-02 0.005831238 -0.008465394 -0.096048096 0.025943442
## Dec 1959 -5.560612e-03 -0.020055492 -0.107556772 0.014404668 -0.002672862
## h=6 h=7 h=8 h=9 h=10
## Jan 1959 0.056335573 0.0684252585 0.054068670 0.056805356 0.076870581
## Feb 1959 0.062455295 0.0476676145 0.050435838 0.070439470 0.081526616
## Mar 1959 0.028042105 0.0308662271 0.050754055 0.061437683 0.049632708
## Apr 1959 0.022360381 0.0422045568 0.052734444 0.041038899 0.026680690
## May 1959 0.017917284 0.0287105591 0.016827504 0.002637848 -0.084990434
## Jun 1959 0.052888667 0.0414558505 0.026851054 -0.060606650 0.061339070
## Jul 1959 0.026293443 0.0118313334 -0.075683148 0.046283111 0.029154974
## Aug 1959 -0.001194799 -0.0887752281 0.033215459 0.015836775 -0.003794462
## Sep 1959 -0.083033533 0.0389750336 0.021418732 0.001847635 0.023536986
## Oct 1959 0.038970357 0.0214416788 0.001861036 0.023571438 -0.015441825
## Nov 1959 0.008556191 -0.0110722121 0.010745640 -0.028178377 -0.024906352
## Dec 1959 -0.022399438 -0.0003481843 -0.039078108 -0.036238004 -0.002793962
## h=11 h=12
## Jan 1959 0.088169803 0.07592215
## Feb 1959 0.069432339 0.05542497
## Mar 1959 0.035373482 -0.05222502
## Apr 1959 -0.060876313 0.06110563
## May 1959 0.037018313 0.01946024
## Jun 1959 0.044429652 0.02465342
## Jul 1959 0.009443843 0.03145714
## Aug 1959 0.018029882 -0.02088874
## Sep 1959 -0.015493751 -0.01198702
## Oct 1959 -0.011973389 0.02141582
## Nov 1959 0.008500399 -0.03052077
## Dec 1959 -0.041751766 -0.04591227types <- paste("h=", formatC(1:12, width=2, flag="0"), sep="")
types## [1] "h=01" "h=02" "h=03" "h=04" "h=05" "h=06" "h=07" "h=08" "h=09" "h=10"
## [11] "h=11" "h=12"acc_df_cv_24 <- make_acc_df("sarima_cv_24", types, apply(sarima_cv_24, 2, calc_rmse))
# 図示
ggplot(acc_df_cv_24, aes(x = type, y = RMSE)) +
geom_bar(stat = "identity", position = "dodge") +
ggtitle("Relationship between RMSE and duration of predictions")