- Notifications
You must be signed in to change notification settings - Fork0
{ffc} R 📦 to fit dynamic functional time series models and produce functional forecasts
License
Unknown, MIT licenses found
Licenses found
nicholasjclark/ffc
Folders and files
| Name | Name | Last commit message | Last commit date | |
|---|---|---|---|---|
Repository files navigation
FunctionalForeCasting
The goal of theffc 📦 is to forecast complex, time-changingfunctional relationships by integrating Generalized Additive Models withdynamic factor functional basis expansions.
Key benefits:
- Model functional responses that change shape over time (not justmagnitude)
- Forecast entire curves into the future, not just single values
- Handle complex multivariate time series with functional structure
- Seamless integration with the powerful
mgcvandfableecosystems
The package introduces dynamic functional predictors using the newfts() term, which decomposes functional time series into time-varyingbasis coefficients that can be forecasted using either independent timeseries models from thefable package or with efficient dynamic factormodels using precompiled Stan models.
You can install the development version offfc fromGitHub with:
# install.packages("pak")pak::pak("nicholasjclark/ffc")library(ffc)# Fit a model with time-varying coefficientsmod <- ffc_gam( response ~ fts(predictor, time_k = 10), data = your_data, time = "time_column", family = gaussian())# Forecast the functional coefficientsfc <- forecast(mod, newdata = future_data, model = "ARDF")library(fable)library(tsibble)library(tidyverse)library(ggplot2); theme_set(theme_bw(base_size = 12))Our aim here is to forecast the number of domestic visitors toMelbourne, Australia. The data can be found in thetsibble::tourismdata set. For now we need to explicitly add thequarter andtimevariables to the data, but in future this will be done automatically forseamless integration with thetsibbleverse
tourism_melb <- tourism |> filter( Region == "Melbourne", Purpose == "Visiting" ) |> mutate( quarter = lubridate::quarter(Quarter), time = dplyr::row_number() )tourism_melb#> # A tsibble: 80 x 7 [1Q]#> # Key: Region, State, Purpose [1]#> Quarter Region State Purpose Trips quarter time#> <qtr> <chr> <chr> <chr> <dbl> <int> <int>#> 1 1998 Q1 Melbourne Victoria Visiting 666. 1 1#> 2 1998 Q2 Melbourne Victoria Visiting 601. 2 2#> 3 1998 Q3 Melbourne Victoria Visiting 529. 3 3#> 4 1998 Q4 Melbourne Victoria Visiting 575. 4 4#> 5 1999 Q1 Melbourne Victoria Visiting 623. 1 5#> 6 1999 Q2 Melbourne Victoria Visiting 530. 2 6#> 7 1999 Q3 Melbourne Victoria Visiting 479. 3 7#> 8 1999 Q4 Melbourne Victoria Visiting 538. 4 8#> 9 2000 Q1 Melbourne Victoria Visiting 618. 1 9#> 10 2000 Q2 Melbourne Victoria Visiting 549. 2 10#> # ℹ 70 more rowsSplit into training and testing folds. We wil aim to forecast the last 5quarters of the data
train <- tourism_melb |> dplyr::slice_head(n = 75)test <- tourism_melb |> dplyr::slice_tail(n = 5)Now fit anffc_gam. We use time-varying level and time-varyingseasonality components, together with a Tweedie observation model(because our outcome,Trips, consists of non-negative real values).This model is simpler so we use the'gam' engine for fitting:
mod <- ffc_gam( Trips ~ # Use mean_only = TRUE to model a time-varying mean fts( time, mean_only = TRUE, time_k = 50, time_m = 1 ) + # Time-varying seasonality fts( quarter, k = 4, time_k = 15, time_m = 1 ), time = "time", data = train, family = tw(), engine = "gam")Theautoplot() method is handy for viewing the time-varying basiscoefficients fromffc_gam() models. Here we draw 10 realisations fromthe estimated coefficient distributions and plot them.
fts_coefs(mod, times = 10, summary = FALSE) |> autoplot()
We can also draw the time-varying basis coefficients using support fromthegratia package, whichhas helpful functions for plotting smooth effects:
gratia::draw(mod)
These plots show the time-varyingcoefficients for a set of basisfunctions. We have one function representing the mean of the series(essentially a constant) as well as three basis functions representingthe quarterly seasonality. Each of these has a coefficient that canchange through time, allowing the entire functional series to changeshape over time. We can compute forecast distribution by fitting thebasis coefficient forecast models in parallel (which is automaticallysupported within thefable package). Here we fit independentexponential smoothing models to each coefficient time series
fc <- forecast( object = mod, newdata = test, model = "ETS", summary = FALSE)We can also convert resulting forecasts to afable object forautomatic plotting and/or scoring of forecasts
# Using the new as_fable method for seamless conversionfc_ffc <- as_fable(mod, newdata = test, forecasts = fc)fc_ffc#> # A fable: 5 x 10 [1Q]#> # Key: Region, State, Purpose [1]#> Quarter Region State Purpose Trips quarter time .dist .mean .model#> <qtr> <chr> <chr> <chr> <dbl> <int> <int> <dist> <dbl> <chr> #> 1 2016 Q4 Melbourne Victor… Visiti… 804. 4 76 sample[200] 839. FFC_E…#> 2 2017 Q1 Melbourne Victor… Visiti… 734. 1 77 sample[200] 760. FFC_E…#> 3 2017 Q2 Melbourne Victor… Visiti… 670. 2 78 sample[200] 774. FFC_E…#> 4 2017 Q3 Melbourne Victor… Visiti… 824. 3 79 sample[200] 747. FFC_E…#> 5 2017 Q4 Melbourne Victor… Visiti… 985. 4 80 sample[200] 836. FFC_E…Leverage the fabletools ecosystem for forecast analysis
# Calculate accuracy metricsaccuracy(fc_ffc, test)#> # A tibble: 1 × 13#> .model Region State Purpose .type ME RMSE MAE MPE MAPE MASE RMSSE#> <chr> <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>#> 1 FFC_ETS Melbour… Vict… Visiti… Test 11.9 90.6 78.4 0.165 9.60 NaN NaN#> # ℹ 1 more variable: ACF1 <dbl># Generate prediction intervals fc_intervals <- hilo(fc_ffc, level = c(80, 95))fc_intervals#> # A tsibble: 5 x 12 [1Q]#> # Key: Region, State, Purpose [1]#> Quarter Region State Purpose Trips quarter time .dist .mean .model#> <qtr> <chr> <chr> <chr> <dbl> <int> <int> <dist> <dbl> <chr> #> 1 2016 Q4 Melbourne Victor… Visiti… 804. 4 76 sample[200] 839. FFC_E…#> 2 2017 Q1 Melbourne Victor… Visiti… 734. 1 77 sample[200] 760. FFC_E…#> 3 2017 Q2 Melbourne Victor… Visiti… 670. 2 78 sample[200] 774. FFC_E…#> 4 2017 Q3 Melbourne Victor… Visiti… 824. 3 79 sample[200] 747. FFC_E…#> 5 2017 Q4 Melbourne Victor… Visiti… 985. 4 80 sample[200] 836. FFC_E…#> # ℹ 2 more variables: `80%` <hilo>, `95%` <hilo># Distribution summariesfc_summary <- fc_ffc |> summarise( mean_forecast = mean(.dist), median_forecast = median(.dist), q25 = quantile(.dist, 0.25), q75 = quantile(.dist, 0.75) )fc_summary#> # A tsibble: 5 x 5 [1Q]#> Quarter mean_forecast median_forecast q25 q75#> <qtr> <dbl> <dbl> <dbl> <dbl>#> 1 2016 Q4 839. 843. 778. 896.#> 2 2017 Q1 760. 758. 711. 815.#> 3 2017 Q2 774. 770. 730. 822.#> 4 2017 Q3 747. 743. 695. 802.#> 5 2017 Q4 836. 831. 785. 898.Next we can explore how to compare forecasts fromffc models totraditional time series models by again leveraging the simplicity andpower of thefable ecosystem
# Generate FFC forecasts with different modelsfc_ffc_arima <- as_fable(mod, newdata = test, model = "ARIMA")fc_ffc_ets <- as_fable(mod, newdata = test, model = "ETS")# Generate traditional model forecastsfc_traditional <- train |> model( ARIMA = ARIMA(Trips), ETS = ETS(Trips) ) |> forecast(h = 5)# Calculate accuracy for all modelsacc_ffc_arima <- accuracy(fc_ffc_arima, test)acc_ffc_ets <- accuracy(fc_ffc_ets, test)acc_traditional <- accuracy(fc_traditional, test)# Extract MAPE values for titlesmape_ffc_arima <- round(acc_ffc_arima$MAPE, 1)mape_ffc_ets <- round(acc_ffc_ets$MAPE, 1)mape_arima <- round(acc_traditional$MAPE[acc_traditional$.model == "ARIMA"], 1)mape_ets <- round(acc_traditional$MAPE[acc_traditional$.model == "ETS"], 1)# Create comparison plotslibrary(patchwork)p1 <- autoplot(fc_ffc_arima, train) + geom_line(data = test, aes(y = Trips), color = "black") + ggtitle(paste0("FFC with ARIMA (MAPE: ", mape_ffc_arima, "%)"))p2 <- autoplot(fc_ffc_ets, train) + geom_line(data = test, aes(y = Trips), color = "black") + ggtitle(paste0("FFC with ETS (MAPE: ", mape_ffc_ets, "%)"))p3 <- autoplot(filter(fc_traditional, .model == "ARIMA"), train) + geom_line(data = test, aes(y = Trips), color = "black") + ggtitle(paste0("Traditional ARIMA (MAPE: ", mape_arima, "%)"))p4 <- autoplot(filter(fc_traditional, .model == "ETS"), train) + geom_line(data = test, aes(y = Trips), color = "black") + ggtitle(paste0("Traditional ETS (MAPE: ", mape_ets, "%)"))(p1 | p2) / (p3 | p4)
These plots illustrate how theffc models outperform traditionalforecasting models for this forecasting experiment by thinking aboutthis as a functional time series problem.
If you encounter a clear bug, please file an issue with a minimalreproducible example onGitHub
Contributions are very welcome, but please see ourCode ofConductwhen you are considering changes that you would like to make.
Theffc project is licensed under anMIT open source license
About
{ffc} R 📦 to fit dynamic functional time series models and produce functional forecasts
Topics
Resources
License
Unknown, MIT licenses found
Licenses found
Code of conduct
Uh oh!
There was an error while loading.Please reload this page.
Stars
Watchers
Forks
Releases
Packages0
Contributors2
Uh oh!
There was an error while loading.Please reload this page.