From 4f28e1ff9a17f0c5b29e1009c64248af3263a182 Mon Sep 17 00:00:00 2001 From: Li Date: Tue, 12 Dec 2023 15:28:41 -0800 Subject: [PATCH 1/7] added ARIMA xreg --- .../R/univariate_forecast_methods.R | 41 ++++++++++ .../ext/r_forecast/_univariate_predictor.py | 25 ++++++- .../r_forecast/test_r_multi_seasonality.py | 75 +++++++++++++++++++ 3 files changed, 140 insertions(+), 1 deletion(-) diff --git a/src/gluonts/ext/r_forecast/R/univariate_forecast_methods.R b/src/gluonts/ext/r_forecast/R/univariate_forecast_methods.R index 39c98d5c41..a8a9efda09 100644 --- a/src/gluonts/ext/r_forecast/R/univariate_forecast_methods.R +++ b/src/gluonts/ext/r_forecast/R/univariate_forecast_methods.R @@ -46,6 +46,47 @@ arima <- function(ts, params) { handleModel(model, params) } + +fourier.arima.xreg <- function(ts, params, xreg_in, xreg_out){ + fourier.frequency.low.periods <- 4 + fourier.ratio.threshold.low.periods <- 18 + fourier.frequency.high.periods <- 52 + fourier.ratio.threshold.high.periods <- 2 + fourier.order <- 4 + + period <- frequency(ts) + len_ts <- length(ts) + fourier_ratio <- len_ts / period + + fourier <- FALSE + + if ((period > fourier.frequency.low.periods + && fourier_ratio > fourier.ratio.threshold.low.periods) + || (period >= fourier.frequency.high.periods + && fourier_ratio > fourier.ratio.threshold.high.periods)) { + # When the period is high, auto.arima becomes unstable + # per Rob's suggestion, we use Fourier series instead + # cf. https://robjhyndman.com/hyndsight/longseasonality/ + fourier <- TRUE + } + + if (fourier == TRUE) { + K <- min(fourier.order, floor(frequency(ts) / 2)) + seasonal <- FALSE + xreg <- forecast::fourier(ts, K=K) + xreg_in <- as.matrix(xreg_in, xreg) + model <- forecast::auto.arima(ts, seasonal = seasonal, xreg = xreg_in, trace=TRUE) + + xreg <- forecast::fourier(ts, K=K, h=params$prediction_length) + xreg_out <- as.matrix(xreg_out, xreg) + + handleModel(model, params, xreg_out) + } else { + model <- forecast::auto.arima(ts, xreg = xreg_in, trace=TRUE) + handleModel(model, params, xreg_out) + } +} + fourier.arima <- function(ts, params){ fourier.frequency.low.periods <- 4 fourier.ratio.threshold.low.periods <- 18 diff --git a/src/gluonts/ext/r_forecast/_univariate_predictor.py b/src/gluonts/ext/r_forecast/_univariate_predictor.py index c8d13d7478..8d4c97d61c 100644 --- a/src/gluonts/ext/r_forecast/_univariate_predictor.py +++ b/src/gluonts/ext/r_forecast/_univariate_predictor.py @@ -34,6 +34,7 @@ "thetaf", "stlar", "fourier.arima", + "fourier.arima.xreg", ] UNIVARIATE_POINT_FORECAST_METHODS = ["croston", "mlp"] SUPPORTED_UNIVARIATE_METHODS = ( @@ -136,7 +137,29 @@ def _get_r_forecast(self, data: Dict) -> Dict: r_params = self._robjects.vectors.ListVector(self.params) vec = self._robjects.FloatVector(data["target"]) ts = make_ts(vec, frequency=self.period) - forecast = self._r_method(ts, r_params) + + if ( + "feat_dynamic_real" in data + and self.method_name == "fourier.arima.xreg" + ): + import rpy2.robjects.numpy2ri + + rpy2.robjects.numpy2ri.activate() + data["feat_dynamic_real"] = np.transpose(data["feat_dynamic_real"]) + nrow, ncol = data["feat_dynamic_real"].shape + xreg_in = self._robjects.r.matrix( + data["feat_dynamic_real"][: -self.prediction_length], + nrow=nrow - self.prediction_length, + ncol=ncol, + ) + xreg_out = self._robjects.r.matrix( + data["feat_dynamic_real"][-self.prediction_length :], + nrow=self.prediction_length, + ncol=ncol, + ) + forecast = self._r_method(ts, r_params, xreg_in, xreg_out) + else: + forecast = self._r_method(ts, r_params) forecast_dict = dict(zip(forecast.names, map(unlist, list(forecast)))) diff --git a/test/ext/r_forecast/test_r_multi_seasonality.py b/test/ext/r_forecast/test_r_multi_seasonality.py index 2b88468af7..e611f86838 100644 --- a/test/ext/r_forecast/test_r_multi_seasonality.py +++ b/test/ext/r_forecast/test_r_multi_seasonality.py @@ -145,3 +145,78 @@ def test_compare_arimas(): ) assert fourier_arima_eval_metrics["MASE"] < arima_eval_metrics["MASE"] + + +dataset_xreg = [ + { + "start": pd.Period("1990-01-01 00", freq=freq), + "target": np.array( + [ + item + for i in range(21) + for item in np.sin( + 2 * np.pi / period * np.arange(1, period + 1, 1) + ) + ] + ) + + np.random.normal(0, 0.5, period * 21) + + np.array( + [ + item + for i in range(3) + for item in [0 for i in range(167)] + [8 for i in range(0, 1)] + ] + ), + "feat_dynamic_real": np.array( + [ + [ + item + for i in range(3) + for item in [0 for i in range(167)] + + [1 for i in range(0, 1)] + ] + ] + ), + } +] + + +def test_compare_arimas_xreg(): + evaluator = Evaluator(quantiles=[0.1, 0.5, 0.9]) + prediction_length = 24 * 7 + + fourier_arima = RForecastPredictor( + freq=freq, + prediction_length=prediction_length, + period=period, + params={ + "quantiles": [0.50, 0.10, 0.90], + "output_types": ["mean", "quantiles"], + }, + method_name="fourier.arima", + ) + fourier_arima_eval_metrics, _ = backtest_metrics( + test_dataset=dataset_xreg, predictor=fourier_arima, evaluator=evaluator + ) + + fourier_arima_xreg = RForecastPredictor( + freq=freq, + prediction_length=prediction_length, + period=period, + params={ + "quantiles": [0.50, 0.10, 0.90], + "output_types": ["mean", "quantiles"], + }, + method_name="fourier.arima.xreg", + ) + + fourier_arima_xreg_eval_metrics, _ = backtest_metrics( + test_dataset=dataset_xreg, + predictor=fourier_arima_xreg, + evaluator=evaluator, + ) + + assert ( + fourier_arima_xreg_eval_metrics["mean_wQuantileLoss"] + < fourier_arima_eval_metrics["mean_wQuantileLoss"] + ) From 81ef66e974fcceb0ea104c61fa543fd419471b15 Mon Sep 17 00:00:00 2001 From: Li Date: Wed, 13 Dec 2023 18:17:26 -0800 Subject: [PATCH 2/7] arima xreg test --- test/ext/r_forecast/test_r_multi_seasonality.py | 1 + 1 file changed, 1 insertion(+) diff --git a/test/ext/r_forecast/test_r_multi_seasonality.py b/test/ext/r_forecast/test_r_multi_seasonality.py index e611f86838..f6964265ae 100644 --- a/test/ext/r_forecast/test_r_multi_seasonality.py +++ b/test/ext/r_forecast/test_r_multi_seasonality.py @@ -146,6 +146,7 @@ def test_compare_arimas(): assert fourier_arima_eval_metrics["MASE"] < arima_eval_metrics["MASE"] +## Below shows improvement in metric when proper x_regressors are included dataset_xreg = [ { From 4af2b08dd14a28c7ba0a241cddb0d8abe07a59a1 Mon Sep 17 00:00:00 2001 From: Li Date: Wed, 13 Dec 2023 18:26:06 -0800 Subject: [PATCH 3/7] . --- test/ext/r_forecast/test_r_multi_seasonality.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/ext/r_forecast/test_r_multi_seasonality.py b/test/ext/r_forecast/test_r_multi_seasonality.py index f6964265ae..314eb86ef8 100644 --- a/test/ext/r_forecast/test_r_multi_seasonality.py +++ b/test/ext/r_forecast/test_r_multi_seasonality.py @@ -146,7 +146,7 @@ def test_compare_arimas(): assert fourier_arima_eval_metrics["MASE"] < arima_eval_metrics["MASE"] -## Below shows improvement in metric when proper x_regressors are included +## Below shows improvement in metric when proper x_regressors are included # dataset_xreg = [ { From 5ba3c5771275eb842000f0cf181be8eab8cbcf83 Mon Sep 17 00:00:00 2001 From: Li Date: Sat, 25 May 2024 21:52:53 -0700 Subject: [PATCH 4/7] fix arima xreg --- .../R/univariate_forecast_methods.R | 75 ++++++++++--------- 1 file changed, 40 insertions(+), 35 deletions(-) diff --git a/src/gluonts/ext/r_forecast/R/univariate_forecast_methods.R b/src/gluonts/ext/r_forecast/R/univariate_forecast_methods.R index a8a9efda09..ca48708154 100644 --- a/src/gluonts/ext/r_forecast/R/univariate_forecast_methods.R +++ b/src/gluonts/ext/r_forecast/R/univariate_forecast_methods.R @@ -48,42 +48,47 @@ arima <- function(ts, params) { fourier.arima.xreg <- function(ts, params, xreg_in, xreg_out){ - fourier.frequency.low.periods <- 4 - fourier.ratio.threshold.low.periods <- 18 - fourier.frequency.high.periods <- 52 - fourier.ratio.threshold.high.periods <- 2 - fourier.order <- 4 - - period <- frequency(ts) - len_ts <- length(ts) - fourier_ratio <- len_ts / period - - fourier <- FALSE - - if ((period > fourier.frequency.low.periods - && fourier_ratio > fourier.ratio.threshold.low.periods) - || (period >= fourier.frequency.high.periods - && fourier_ratio > fourier.ratio.threshold.high.periods)) { - # When the period is high, auto.arima becomes unstable - # per Rob's suggestion, we use Fourier series instead - # cf. https://robjhyndman.com/hyndsight/longseasonality/ - fourier <- TRUE - } - - if (fourier == TRUE) { - K <- min(fourier.order, floor(frequency(ts) / 2)) - seasonal <- FALSE - xreg <- forecast::fourier(ts, K=K) - xreg_in <- as.matrix(xreg_in, xreg) - model <- forecast::auto.arima(ts, seasonal = seasonal, xreg = xreg_in, trace=TRUE) - - xreg <- forecast::fourier(ts, K=K, h=params$prediction_length) - xreg_out <- as.matrix(xreg_out, xreg) - - handleModel(model, params, xreg_out) + if (missing(xreg_in)){ + fourier.arima(ts, params) } else { - model <- forecast::auto.arima(ts, xreg = xreg_in, trace=TRUE) - handleModel(model, params, xreg_out) + + fourier.frequency.low.periods <- 4 + fourier.ratio.threshold.low.periods <- 18 + fourier.frequency.high.periods <- 52 + fourier.ratio.threshold.high.periods <- 2 + fourier.order <- 4 + + period <- frequency(ts) + len_ts <- length(ts) + fourier_ratio <- len_ts / period + + fourier <- FALSE + + if ((period > fourier.frequency.low.periods + && fourier_ratio > fourier.ratio.threshold.low.periods) + || (period >= fourier.frequency.high.periods + && fourier_ratio > fourier.ratio.threshold.high.periods)) { + # When the period is high, auto.arima becomes unstable + # per Rob's suggestion, we use Fourier series instead + # cf. https://robjhyndman.com/hyndsight/longseasonality/ + fourier <- TRUE + } + + if (fourier == TRUE) { + K <- min(fourier.order, floor(frequency(ts) / 2)) + seasonal <- FALSE + xreg <- forecast::fourier(ts, K=K) + xreg_in <- as.matrix(xreg_in, xreg) + model <- forecast::auto.arima(ts, seasonal = seasonal, xreg = xreg_in, trace=TRUE) + + xreg <- forecast::fourier(ts, K=K, h=params$prediction_length) + xreg_out <- as.matrix(xreg_out, xreg) + + handleModel(model, params, xreg_out) + } else { + model <- forecast::auto.arima(ts, xreg = xreg_in, trace=TRUE) + handleModel(model, params, xreg_out) + } } } From fb52def81914b3e8cd9628999a058bf50b6a212f Mon Sep 17 00:00:00 2001 From: Lorenzo Stella Date: Fri, 14 Jun 2024 15:24:01 +0200 Subject: [PATCH 5/7] Update test/ext/r_forecast/test_r_multi_seasonality.py --- test/ext/r_forecast/test_r_multi_seasonality.py | 1 + 1 file changed, 1 insertion(+) diff --git a/test/ext/r_forecast/test_r_multi_seasonality.py b/test/ext/r_forecast/test_r_multi_seasonality.py index 314eb86ef8..c95ab395e0 100644 --- a/test/ext/r_forecast/test_r_multi_seasonality.py +++ b/test/ext/r_forecast/test_r_multi_seasonality.py @@ -146,6 +146,7 @@ def test_compare_arimas(): assert fourier_arima_eval_metrics["MASE"] < arima_eval_metrics["MASE"] + ## Below shows improvement in metric when proper x_regressors are included # dataset_xreg = [ From 31f726b354b9ee068aea434e0feb894c36150bf4 Mon Sep 17 00:00:00 2001 From: Lorenzo Stella Date: Fri, 14 Jun 2024 15:26:46 +0200 Subject: [PATCH 6/7] Update test/ext/r_forecast/test_r_multi_seasonality.py --- test/ext/r_forecast/test_r_multi_seasonality.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/ext/r_forecast/test_r_multi_seasonality.py b/test/ext/r_forecast/test_r_multi_seasonality.py index c95ab395e0..c130f03f74 100644 --- a/test/ext/r_forecast/test_r_multi_seasonality.py +++ b/test/ext/r_forecast/test_r_multi_seasonality.py @@ -147,7 +147,7 @@ def test_compare_arimas(): assert fourier_arima_eval_metrics["MASE"] < arima_eval_metrics["MASE"] -## Below shows improvement in metric when proper x_regressors are included # +## Below shows improvement in metric when proper x_regressors are included # dataset_xreg = [ { From f42156db6f3d9a9720adb181e62a0119ae8fc6ac Mon Sep 17 00:00:00 2001 From: Lorenzo Stella Date: Fri, 14 Jun 2024 15:29:00 +0200 Subject: [PATCH 7/7] Update test_r_univariate_predictor.py --- test/ext/r_forecast/test_r_univariate_predictor.py | 4 ---- 1 file changed, 4 deletions(-) diff --git a/test/ext/r_forecast/test_r_univariate_predictor.py b/test/ext/r_forecast/test_r_univariate_predictor.py index 4eefa1bbb0..eb5b59115c 100644 --- a/test/ext/r_forecast/test_r_univariate_predictor.py +++ b/test/ext/r_forecast/test_r_univariate_predictor.py @@ -48,10 +48,6 @@ def test_forecasts(method_name): "MLP currently does not work because " "the `neuralnet` package is not yet updated with a known bug fix in ` bips-hb/neuralnet`" ) - if method_name == "fourier.arima.xreg": - pytest.xfail( - "Method `fourier.arima.xreg` does not work because of a known issue." - ) dataset = datasets.get_dataset("constant")