--- title: "utsf" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{utsf} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r, include = FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) ``` ```{r setup} library(utsf) ``` In this document, the **utsf** package is described. This package offers a meta-engine for applying different regression models to univariate time series forecasting using an autoregressive approach. One of the main features of the package is its extensibility, which allows you to use machine learning models not directly supported by the package (such as neural networks) or your own custom models. Currently, the package is mainly focused on regression tree models. An explanation of how to forecast univariate time series using regression trees can be found here: . # Univariate time series forecasting (autoregressive models) and recursive forecasts A univariate time series forecasting method is one in which the future values of a series are predicted using only information from the series itself—for example, using its historical mean value as a forecast. An advantage of this type of prediction is that, apart from the series being forecasted, there is no need to collect any further information to train the forecasting model. An autoregressive model is a type of univariate time series forecasting model in which a value of a time series is expressed as a function of some of its past values. That is, an autoregressive model is a regression model in which the independent variables are lagged values (previous values) of the response variable. For example, given a time series with the following historical values: $t = \{1, 3, 6, 7, 9, 11, 16\}$, suppose that we want to develop an autoregressive model in which a target "is explained" by its first, second, and fourth past values (in this context, a previous value is also called a *lag*, so lag 1 is the value immediately preceding a given value in the series). Given the series $t$ and lags (1, 2, and 4), the training set would be: | Lag 4 | Lag 2 | Lag 1 | Target | |-------|-------|-------|--------| | 1 | 6 | 7 | 9 | | 3 | 7 | 9 | 11 | | 6 | 9 | 11 | 16 | ## Recursive forecasts Given a model trained with the previous dataset, the next future value of the series is predicted as $f(\text{Lag}_4, \text{Lag}_2, \text{Lag}_1)$, where $f$ is the regression function and $\text{Lag}_4$, $\text{Lag}_2$, and $\text{Lag}_1$ are the fourth, second, and first lagged values of the next future value. Therefore, the next future value of the series $t$ is predicted as $f(7, 11, 16)$, producing a value that will be denoted as $F_1$. Suppose that the *forecast horizon* (the number of future values to be forecasted into the future) is greater than 1. If the regression function only predicts the next future value of the series, a recursive approach can be applied to forecast all the values within the horizon. Using this approach, the regression function is applied recursively until all multi-step-ahead values are forecasted. For instance, following the previous example, suppose that the forecast horizon is 3. As explained, to forecast the next value of the series (one-step-ahead) the regression function is fed with the vector $[7, 11, 16]$, producing $F_1$. To forecast the two-step-ahead value the regression function is fed with the vector $[9, 16, F_1]$. The forecast for the one-step-ahead value, $F_1$, is used as the first lag for the two-step-ahead forecast, because the actual value is unknown. Finally, to predict the three-step-ahead value, the regression function is fed with the vector [$11, F_1, F_2]$. This example of recursive forecast is summarized in the table below: |Steps ahead | Autoregressive values | Forecast | |------------|-----------------------|----------| | 1 | 7, 11, 16 | $F_1$ | | 2 | 9, 16, $F_1$ | $F_2$ | | 3 | 11, $F_1$, $F_2$ | $F_3$ | The recursive approach to multi-step-ahead forecasting is employed in classical statistical models, such as ARIMA and exponential smoothing. # The utsf package The **utsf** package simplifies the use of classical regression models for univariate time series forecasting, employing the autoregressive approach and the recursive prediction strategy explained in the previous section. All the supported models are applied using a uniform interface: * the `create_model()` function to build the forecasting model, and * the `forecast()` function for using the trained model to predict future values. Let us see an example in which a regression tree model is used to forecast the future values of a time series: ```{r} m <- create_model(AirPassengers, lags = 1:12, method = "rt") f <- forecast(m, h = 12) ``` In this example, an autoregressive tree model (`method = "rt"`, Regression Tree) is trained using the historical values of the `AirPassengers` time series and a forecast for its next 12 values (`h = 12`) is generated. The `create_model()` function returns an S3 object of class `utsf` with information about the trained model. The forecast information is included in the object returned by the `forecast()` function, as a component named `pred` of class `ts` (a time series): ```{r} f$pred library(ggplot2) autoplot(f) ``` The training set used to fit the model is built from the historical values of the time series using the autoregressive approach explained in the previous section. The `lags` parameter of the `create_model()` function specifies the autoregressive lags. In this example: `lags = 1:12`, meaning each target is a function of its 12 previous values. Next, we check the first targets (and their associated features) used to train the regression model: ```{r} head(m$targets) # first targets head(m$features) # and its associated features ``` The curious reader may have noticed that the features and targets are not on the same scale as the original time series. This is because, by default, a transformation is applied to the examples of the training set. This transformation will be explained later. Let us inspect the training set associated with the example from the previous section: ```{r} t <- ts(c(1, 3, 6, 7, 9, 11, 16)) out <- create_model(t, lags = c(1, 2, 4), trend = "none") cbind(out$features, Target = out$targets) ``` No transformation has been applied this time (`trend = "none"`). # Prediction intervals Prediction intervals for a forecast can be computed: ```{r} m <- create_model(USAccDeaths, method = "mt") f <- forecast(m, h = 12, PI = TRUE, level = 90) f library(ggplot2) autoplot(f) ``` Prediction intervals are calculated following the guidelines in [Forecasting: Principles and Practice](https://otexts.com/fpp3/nnetar.html#prediction-intervals-5). Random errors are assumed to follow a normal distribution. In this example, a 90% prediction interval (`level = 90`) was computed. # Supported models The `create_model()` and `forecast()` functions provide a unified interface for applying an autoregressive approach to time series forecasting with various regression models. These models are implemented across several R packages. Currently, our project focuses primarily on regression tree models, supporting the following approaches: * k-nearest neighbors: In this case no model is trained and the function `FNN::knn.reg()` is used, as regression function, to recursively predict the future values of the time series. * Linear models: The model is trained using the function `stats::lm()` and its associated method `stats::predict.lm()` is applied recursively for the forecasts, i.e., as regression function. * Regression trees: The model is trained using the function `rpart::rpart()` and its associated method `rpart::predict.rpart()` is used for the forecasts. * Model trees: The model is trained with the function `Cubist::cubist()` and its associated method `Cubist::predict.cubist()` is used for predictions. * Bagging: The model is trained with the function `ipred::bagging()` and its associated method `ipred::predict.regbagg()` is used for forecasting. * Random forest: The model is trained with the function `ranger::ranger()` and its associated method `ranger::predict.ranger()` is used for predictions. * Extreme gradient boosting: The model is trained with the function `xgboost::xgboost()` and its associated method `xgboost::predict.xgboost()` is used for predictions. The S3 object of class `utsf` returned by the `create_model()` function contains a component with the trained autoregressive model: ```{r} m <- create_model(fdeaths, lags = 1:12, method = "rt") m$model ``` In this case, the model is obtained by training a regression tree using the `rpart::rpart()` function on a training set composed of the features `m$features` and targets `m$targets`. Once trained, the `rpart::predict.rpart()` function can be used recursively within `forecast()` to predict future values of the time series. # Using your own models An appealing feature of the **utsf** package is its extensibility. Beyond the natively supported models, you can leverage the package's infrastructure to perform autoregressive time series forecasting with custom regression models. This allows user-defined models to benefit from built-in functionality, such as training set construction, recursive forecasting, prediction intervals, preprocessing, accuracy evaluation, and hyperparameter tuning. To use a custom regression model, pass a function capable of training your model to the `method` parameter of the `create_model()` function. We refer to this function as the *modeling function*, and it must return an object containing the trained regression model. The modeling function must accept three input parameters: * `X`: a data frame containing the features of the training examples. This data frame is constructed from the time series according to the autoregressive lags, as explained in a previous section. It is identical to the `features` component of the object returned by `create_model()`. * `y`: a vector containing the targets of the training examples, constructed as described in a previous section. It is identical to the `targets` component of the object returned by `create_model()`. * `...`: additional arguments passed to fine-tune model behavior. This parameter is not used in this section and will be detailed in the next. Furthermore, if the modeling function returns an object of class `model_class`, an S3 method with the signature `predict.model_class(object, new_value)` must be implemented. This method uses the model to predict a single new value (acting as the prediction function associated with the model). Its parameters are described below: * `object`: the object of class `model_class` returned by your modeling function. * `new_value`: a one-row data frame with the same structure as `X`, containing the features of the instance to be predicted. In summary, to use a custom regression model, you must provide: * A modeling function that trains the model given a training set. * A prediction function that predicts a single target value given a trained model and a new feature vector. The following examples demonstrate how to use this functionality to apply custom models to autoregressive time series forecasting, leveraging the capabilities of the **utsf** package. ## Example 1: k-nearest neighbors The k-nearest neighbors (k-NN) algorithm is a special case because it is a lazy learner and builds no explicit model. Instead, it simply stores the training set, leaving all the computation to the prediction function. Below is an initial example using the k-NN implementation from the **FNN** package: ```{r} # Modeling function: in this case (k-NN) just stores the training set my_knn_model <- function(X, y, ...) { structure(list(X = X, y = y), class = "my_knn") } # Regression function predict.my_knn <- function(object, new_value) { FNN::knn.reg(train = object$X, test = new_value, y = object$y)$pred } m <- create_model(AirPassengers, lags = 1:12, method = my_knn_model) f <- forecast(m, h = 12) f$pred autoplot(f) ``` The modeling function (`my_knn_model()`) creates a "model" that simply stores the training set. The regression function (`predict.my_knn()`) uses `FNN::knn.reg()` to find the $k$-nearest neighbors and compute their mean response value. The default number of neighbors in `FNN::knn.reg()` (3) is used. Later, we will modify this example to allow users to specify the value of $k$. The $k$-nearest neighbors algorithm is simple enough to be easily implemented from scratch, without relying on external R packages: ```{r} # Modeling function my_knn_model2 <- function(X, y, ...) { structure(list(X = X, y = y), class = "my_knn2") } # Regression function predict.my_knn2 <- function(object, new_value) { k <- 3 # number of nearest neighbors distances <- sapply(1:nrow(object$X), function(i) sum((object$X[i, ] - new_value)^2)) k_nearest <- order(distances)[1:k] mean(object$y[k_nearest]) } m2 <- create_model(AirPassengers, lags = 1:12, method = my_knn_model2) forecast(m2, h = 12)$pred ``` ## Example 2: random forest and neural networks The random forest algorithm is natively supported in the **utsf** package via the `ranger` package implementation. Here, we build an autoregressive model using the random forest implementation provided by the `randomForest` package: ```{r} library(randomForest) my_model <- function(X, y, ...) { randomForest(x = X, y = y) } m <- create_model(USAccDeaths, lags = 1:12, method = my_model) f <- forecast(m, h = 12) library(ggplot2) autoplot(f) print(m$model) ``` The modeling function (`my_model()`) simply calls `randomForest::randomForest()`, which returns an S3 object of class `randomForest` containing the trained model. In this case, we do not need to implement a regression function, as the existing `predict.randomForest()` S3 method already performs the required task. As another example, we implement an autoregressive neural network model using the **nnet** package: ```{r} library(nnet) my_model <- function(X, y, ...) { nnet(x = X, y = y, size = 5, linout = TRUE, trace = FALSE) } m <- create_model(USAccDeaths, lags = 1:12, method = my_model) f <- forecast(m, h = 12) library(ggplot2) autoplot(f) ``` In this case, `nnet()` returns an S3 object of class `nnet`. Once again, the class's default prediction method, `predict.nnet()`, provides the required functionality. Given that neural networks are sensitive to feature scales, preprocessing the time series would likely improve forecast accuracy. # Setting the parameters of the regression models Typically, a regression model can be tuned using hyperparameters. By default, the models supported by the package use specific default values defined by their underlying training functions (listed in a previous section). However, users can customize these hyperparameters to adjust model behavior. Below is an example: ```{r} # A bagging model set with default parameters m <- create_model(USAccDeaths, lags = 1:12, method = "bagging") length(m$model$mtrees) # number of regression trees (25 by default) # A bagging model set with 3 regression trees m2 <- create_model(USAccDeaths, lags = 1:12, method = "bagging", nbagg = 3 ) length(m2$model$mtrees) # number of regression trees ``` In the previous example, two bagging models (using regression trees) are trained with `create_model()`. The first model uses 25 trees (the default for `ipred::ipredbagg()`) while the second sets the number of trees to 3. To customize hyperparameters, consult the parameter documentation of the training function used internally by **utsf**. In this case, `ipred::ipredbagg()` uses the `nbagg` parameter to specify the number of trees in the ensemble. ## Setting the paramaters of a regression model supplied by the user When using a custom regression model, its parameters can also be tuned. As explained previously, custom models require a modeling function with three parameters: * `X`: a data frame containing the training features. * `y`: a vector containing the target values. * `...`: additional arguments passed to adjust model behavior. The `...` parameter allows the modeling function to accept an arbitrary number of arguments to tune model behavior. Below are several examples of passing arguments to the custom regression models implemented in the previous section: ### Example 1: k-nearest neighbors In this case, the model can be tuned by specifying the value of $k$. The modeling function accepts a parameter named `k` to set the number of nearest neighbors: ```{r} # Function to train the model my_knn_model <- function(X, y, ...) { param <- list(...) k <- if ("k" %in% names(param)) param$k else 3 structure(list(X = X, y = y, k = k), class = "my_knn") } # Regression function for object of class my_knn predict.my_knn <- function(object, new_value) { FNN::knn.reg(train = object$X, test = new_value, y = object$y, k = object$k)$pred } # The model is trained with default parameters (k = 3) m <- create_model(AirPassengers, lags = 1:12, method = my_knn_model) print(m$model$k) # The model is trained with k = 5 m2 <- create_model(AirPassengers, method = my_knn_model, k = 5) print(m2$model$k) ``` Additional arguments specified in `create_model()` (such as `k = 5`) are passed directly to the custom modeling function (`my_knn_model()`). ### Example 2: random forest As another example, suppose we create a random forest model using the **randomForest** package. By default, the model relies on the default parameters of `randomForest::randomForest()`, except for the number of trees, which is set to 200 (`ntree = 200`). However, users can override this and train the model with custom parameter values: ```{r} library(randomForest) my_model <- function(X, y, ...) { param <- list(...) args <- list(x = X, y = y, ntree = 200) # default parameters args <- args[!(names(args) %in% names(param))] args <- c(args, param) do.call(randomForest::randomForest, args = args) } # The random forest is built with our default parameters m <- create_model(USAccDeaths, method = my_model) print(m$model$ntree) print(m$model$mtry) # The random forest is built with ntree = 400 and mtry = 6 m <- create_model(USAccDeaths, method = my_model, ntree = 400, mtry = 6) print(m$model$ntree) print(m$model$mtry) ``` # Estimating forecast accuracy This section demonstrates how to estimate the forecast accuracy of a regression model when predicting a time series with **utsf**. Below is an example: ```{r} m <- create_model(UKgas, lags = 1:4, method = "knn") r <- efa(m, h = 4, type = "normal", size = 8) r$per_horizon r$global ``` To assess model forecast accuracy, use the `efa()` function on the fitted model. In this example, the forecasting method is a $k$-nearest neighbors algorithm (`method = "knn"`) with autoregressive lags 1 to 4, applied to the `UKgas` time series. The function estimates forecast accuracy across a 4-step-ahead horizon (`h = 4`) using multiple evaluation metrics. The `efa()` function returns a list with two components: * `per_horizon`: Contains accuracy estimates for each individual forecasting step. * `global`: Contains the overall accuracy estimates, calculated as the row means of the `per_horizon` matrix. Currently, the following forecast accuracy measures are computed: * MAE: mean absolute error * MAPE: mean absolute percentage error * sMAPE: symmetric MAPE * RMSE: root mean squared error Next, we describe how the forecasting accuracy measures are computed for a forecasting horizon $h$ ($y_t$ and $\hat{y}_t$ are the actual future value and its forecast for horizon $t$ respectively): $$ MAE = \frac{1}{h}\sum_{t=1}^{h} |y_t-\hat{y}_t| $$ $$ MAPE = \frac{1}{h}\sum_{t=1}^{h} 100\frac{|y_t-\hat{y}_t|}{y_t} $$ $$ sMAPE = \frac{1}{h}\sum_{t=1}^{h} 200\frac{\left|y_{t}-\hat{y}_{t}\right|}{|y_t|+|\hat{y}_t|} $$ $$ RMSE = \sqrt{\frac{1}{h}\sum_{t=1}^{h} (y_t-\hat{y}_t)^2} $$ Forecast accuracy is estimated using the well-known rolling origin evaluation strategy. This technique generates multiple training and test sets from the time series, as illustrated below for a 4-step-ahead forecasting horizon (`h = 4`): ```{r out.width = '85%', echo = FALSE} knitr::include_graphics("ro_normal.png") ``` In this setup, the last nine observations of the time series are reserved to construct the test sets. Six models are fitted using the expanding training sets, and their forecasts are evaluated against the corresponding test values. Consequently, the forecast error for each step of the horizon is averaged across six evaluation folds. The number of trailing observations used to form the test sets can be controlled via the `size` or `prop` parameters: ```{r, eval = FALSE} m <- create_model(UKgas, lags = 1:4, method = "knn") # Use the last 9 observations of the series to build test sets r1 <- efa(m, h = 4, type = "normal", size = 9) # Use the last 20% observations of the series to build test sets r2 <- efa(m, h = 4, type = "normal", prop = 0.2) ``` If the series is short, the training sets may be too small to fit a meaningful model. In such cases, an alternative rolling origin scheme can be applied by setting the `type` parameter to `"minimum"`: ```{r} m <- create_model(UKgas, lags = 1:4, method = "knn") r <- efa(m, h = 4, type = "minimum") r$per_horizon r$global ``` When this option is selected, the last `h` observations (4 in this example) are used to construct the test sets, as illustrated below for `h = 4`: ```{r out.width = '85%', echo = FALSE} knitr::include_graphics("ro_minimum.png") ``` In this case, the estimated forecast error for horizon 1 is based on 4 errors, for horizon 2 on 3 errors, for horizon 3 on 2 errors, and for horizon 4 on 1 error. The `efa()` function also returns the test sets (and their associated predictions) used during accuracy assessment: ```{r} timeS <- ts(1:25) m <- create_model(timeS, lags = 1:3, method = "mt") r <- efa(m, h = 5, size = 7) r$test_sets r$predictions ``` Each row in these tables represents a test set or a set of predictions. Below are the test sets generated when using the `"minimum"` evaluation scheme: ```{r} timeS <- ts(1:25) m <- create_model(timeS, lags = 1:3, method = "mt") r <- efa(m, h = 3, type = "minimum") r$test_sets r$predictions ``` # Parameter tuning Another key feature of **utsf** is parameter tuning. The `tune_grid()` function estimates the forecast accuracy of a model across various parameter combinations. The optimal configuration (the one achieving the best forecast accuracy) is then used to train a final model on the full historical series and generate future forecasts. Below is an example: ```{r} m <- create_model(UKgas, lags = 1:4, method = "knn") r <- tune_grid(m, h = 4, tuneGrid = expand.grid(k = 1:7), type = "normal", size = 8) # To check the estimated forecast accuracy with the different configurations r$tuneGrid ``` In this example, the `tuneGrid` parameter specifies (as a data frame) the set of hyperparameters to evaluate. The forecast accuracy for each parameter combination is estimated as described in the previous section, using the most recent observations of the time series as a validation set. The `size` parameter controls the length of this test set. The `tuneGrid` component of the list returned by `tune_grid()` contains these estimation results. Here, the $k$-nearest neighbors model with $k = 1$ achieves the best performance across all forecast accuracy metrics, which is driven by the increasing variance of the series. The optimal parameter combination according to RMSE is then used to generate the final time series forecast: ```{r} r$best r$forecast ``` Below, we plot the $k$ values against their estimated forecast accuracy in terms of RMSE: ```{r} plot(r$tuneGrid$k, r$tuneGrid$RMSE, type = "o", pch = 19, xlab = "k (number of nearest neighbors)", ylab = "RMSE", main = "Estimated accuracy") ``` The following example illustrates tuning a model with multiple hyperparameters. Here, we evaluate a random forest model using three tuning parameters: ```{r} m <- create_model(UKDriverDeaths, lags = 1:12, method = "rf") r <- tune_grid(m, h = 12, tuneGrid = expand.grid(num.trees = c(200, 500), replace = c(TRUE, FALSE), mtry = c(4, 8)), type = "normal", size = 12 ) r$tuneGrid ``` # Preprocessings and transformations Model forecast accuracy can often be improved by applying transformations or preprocessing techniques to the target time series. Currently, **utsf** focuses on preprocessing methods designed for trended time series. This is because most models currently integrated into the package (such as regression trees and $k$-nearest neighbors) predict averages of the target values in the training set (i.e., historical observations). When a series exhibits a trend, these historical averages typically fall outside the range of future values. The following example demonstrates forecasting a trended series (the annual revenue passenger miles flown by commercial airlines in the United States) using a random forest model: ```{r} m <- create_model(airmiles, lags = 1:4, method = "rf", trend = "none") f <- forecast(m, h = 4) autoplot(f) ``` In this case, no preprocessing is applied to handle the trend (`trend = "none"`). As shown below, the forecast fails to capture the upward trajectory of the series because regression tree models can only predict average values of the training targets. The following subsections explain how to handle trended series using three different transformations: ## Differencing A common approach for handling trended series is first-differencing (computing the differences between consecutive observations). The model is then trained on the differenced series, and the resulting forecasts are back-transformed to the original scale. Below is an example: ```{r} m <- create_model(airmiles, lags = 1:4, method = "rf", trend = "differences", nfd = 1) f <- forecast(m, h = 4) autoplot(f) ``` In this example, the `trend` parameter is set to `"differences"` to apply first-differencing. The order of differencing is specified via `nfd` (which typically defaults to 1). Setting `nfd = -1` automatically estimates the required order of differencing using the `ndiffs()` function from the **forecast** package. In this case, if the estimated order is 0, no differencing is applied. ```{r} # The order of first differences is estimated using the ndiffs function m <- create_model(airmiles, lags = 1:4, method = "rf", trend = "differences", nfd = -1) m$differences ``` ## The additive transformation This transformation is also used for dealing with trended series in other R packages, such as **tsfknn** or **tsfgrnn**. The additive transformation modifies the targets of the training set as follows: * The target associated with a feature vector is transformed by subtracting the mean of that feature vector. * To back-transform a prediction, the mean of the corresponding input vector is added to it. The effect of the additive transformation on the training examples can be easily inspected using the package API. For instance, consider the training set of a model fit without preprocessing: ```{r} timeS <- ts(c(1, 3, 7, 9, 10, 12)) m <- create_model(timeS, lags = 1:2, trend = "none") cbind(m$features, Targets = m$targets) ``` Now, consider the effect of the additive transformation on the training examples: ```{r} timeS <- ts(c(1, 3, 7, 9, 10, 12)) m <- create_model(timeS, lags = 1:2, trend = "additive", transform_features = FALSE) cbind(m$features, Targets = m$targets) ``` In addition to transforming the targets, the features can also be transformed. A feature vector is transformed by subtracting its mean from each element. The effect of this transformation is shown below: ```{r} timeS <- ts(c(1, 3, 7, 9, 10, 12)) m <- create_model(timeS, lags = 1:2, trend = "additive", transform_features = TRUE) cbind(m$features, Targets = m$targets) ``` Transforming the features removes level effects from the feature space, which empirically improves forecast accuracy. Below, we forecast the `airmiles` time series using the additive transformation with a random forest model: ```{r} m <- create_model(airmiles, lags = 1:4, method = "rf") f <- forecast(m, h = 4) autoplot(f) ``` By default, the additive transformation is applied to both targets and features. ## The multiplicative transformation The multiplicative transformation is similar to the additive transformation, but it is designed for series with exponential trends, whereas the additive transformation is better suited for linear trends. Under the multiplicative transformation, training examples are processed as follows: * The target associated with a feature vector is transformed by dividing it by the mean of that feature vector. * To back-transform a prediction, it is multiplied by the mean of the corresponding input vector. * Optionally, each element of a feature vector is transformed by dividing it by the vector's mean. Below is an example using a synthetic time series with an exponential trend, comparing forecasts generated with the additive and multiplicative transformations: ```{r} t <- ts(10 * 1.05^(1:20)) m_m <- create_model(t, lags = 1:3, method = "rf", trend = "multiplicative") f_m <- forecast(m_m, h = 4) m_a <- create_model(t, lags = 1:3, method = "rf", trend = "additive") f_a <- forecast(m_a, h = 4) library(vctsfr) plot_predictions(t, predictions = list(Multiplicative = f_m$pred, Additive = f_a$pred)) ``` In this case, the forecast generated with the multiplicative transformation captures the exponential trend perfectly. As an additional example, we compare the additive and multiplicative transformations applied to the `AirPassengers` series (a monthly series exhibiting an exponential trend): ```{r, fig.width= 6} m_m <- create_model(AirPassengers, lags = 1:12, method = "knn", trend = "multiplicative") f_m <- forecast(m_m, h = 36) m_a <- create_model(AirPassengers, lags = 1:12, method = "knn", trend = "additive") f_a <- forecast(m_a, h = 36) library(vctsfr) plot_predictions(AirPassengers, predictions = list(Multiplicative = f_m$pred, Additive = f_a$pred)) ``` # Default parameters The `create_model()` function requires only one mandatory parameter: the time series used to train the model. Below, we describe the default settings for the remaining parameters: * `lags`: An integer vector defining the autoregressive lags. If `frequency(ts) == f` (where `ts` is the input time series and $f > 1$), the default lags are `1:f` (e.g., `1:4` for quarterly data and `1:12` for monthly data), capturing potential seasonal patterns. If `frequency(ts) == 1`: * Significant lags are automatically selected based on the partial autocorrelation function (`stats::pacf()`). * If no lag is statistically significant, lags `1:5` are selected. * If only a single lag is significant and the additive or multiplicative transformation is applied to both targets and features, lags `1:5` are used instead. This prevents applying feature-level transformations when only one lag is present, as scaling a single feature by its own mean is uninformative. * `method`: Defaults to the $k$-nearest neighbors algorithm. * `trend`: Defaults to the additive transformation (`"additive"`). This method has proven effective for trended time series and generally improves performance across most series types. * `transform_features`: Defaults to `TRUE`, ensuring that both features and targets are transformed when using the additive or multiplicative transformation. # More on using your custom models In previous sections we have explained how to leverage the package's features to apply your own models for time series forecasting. In this section we describe how to take advantage of other package's features. We will employ the random forest model implemented in the **randomForest** package used previously. First, let us examine the computation of prediction intervals: ```{r} library(randomForest) myModel <- function(X, y, ...) { param <- list(...) args <- list(x = X, y = y, ntree = 200) args <- args[!(names(args) %in% names(param))] args <- c(args, param) do.call(randomForest::randomForest, args = args) } m <- create_model(USAccDeaths, method = myModel) f <- forecast(m, h = 12, level = 90, PI = TRUE) f ``` As can be seen, all these features are used in the same way as with the built-in models. For example, let us estimate the forecast accuracy of the previous model for a forecasting horizon of 5: ```{r} efa(m, h = 5, type = "normal", size = 8)$per_horizon ``` By default, the additive transformation is applied. Let us select the multiplicative transformation to estimate the forecast accuracy under this setting: ```{r} m <- create_model(USAccDeaths, lags = 1:12, method = myModel, trend = "multiplicative") efa(m, h = 5, type = "normal", size = 8)$per_horizon ``` Finally, we perform parameter tuning over different values for the `mtry` and `ntree` arguments: ```{r} m <- create_model(USAccDeaths, method = myModel) r <- tune_grid(m, h = 6, tuneGrid = expand.grid(mtry = 3:6, ntree = c(200, 500)), type = "normal" ) r$tuneGrid ```