remove(list=ls())
library(fpp3)
library(fredr)
library(tidyverse)
library(patchwork)
library(knitr)
library(kableExtra)
library(gridExtra)
library(writexl)
library(tseries)
library(tsibbledata)
library(forecast)
library(quantmod)
library(feasts)
library(broom)
Sys.setlocale("LC_TIME", "en_US.UTF-8")Module 10 Discussion_HTS>S
Part I: Theoretical Explanation of HTS and GTS
1. Differences between HTS and GTS
The primary distinction between Hierarchical Time Series (HTS) and Grouped Time Series (GTS) lies in the structural relationship between the levels.
- HTS (Hierarchical): Levels are organized in a strict nested (parent-child) structure. A higher level is the direct aggregation of its unique sub-categories.
- Example: Geography (Country > State > Region) or Pharmaceutical classification (ATC1 > ATC2).
- GTS (Grouped): Levels are defined by crossing factors where categories do not have a unique hierarchical nesting. Different grouping attributes can overlap.
- Example: Tourism data grouped by “purpose of travel” (business, holiday) and “region”. A “Business” trip can occur in any state, meaning the categories are independent dimensions.
2. Forecast Method vs. Reconciliation Method
It is crucial to distinguish between the act of forecasting and the act of ensuring consistency across the system:
- Forecast Method: This refers to generating base forecasts (\(\hat{y}_h\)) for each individual series independently. These forecasts are usually “incoherent,” meaning the sum of the bottom-level forecasts does not match the top-level forecast.
- Reconciliation Method: This is the process of adjusting the base forecasts to ensure aggregation consistency. It uses a mapping matrix to resolve discrepancies so that the reconciled forecasts (\(\tilde{y}_h\)) are “coherent” (the hierarchy sums up correctly).
3. Comparison of the Four Reconciliation Methods
| Method | Definition | Strengths | Weaknesses |
|---|---|---|---|
| Top-down | Forecasts the total (top level) and disaggregates it to the bottom levels based on historical proportions. | Simple to implement; works well when top-level data is stable with less noise. | Ignores specific trends or seasonality present only in the sub-series. |
| Bottom-up | Forecasts all bottom-level series and aggregates them upward to get the totals. | Captures individual dynamics of each sub-series. | Bottom-level data is often noisy; ignores correlations between levels. |
| Middle-out | Starts at an intermediate level. Uses bottom-up for levels above it and top-down for levels below it. | A compromise that balances the stability of higher levels with the detail of lower levels. | Requires subjective judgment to select the “middle” level. |
| MinT (Minimum Trace) | An optimal reconciliation approach that uses a weighted least squares framework based on the forecast error covariance matrix. | Most accurate & robust; minimizes the total forecast variance and utilizes information from all levels. | Computationally more intensive as it requires estimating the covariance matrix. |
Part II:
Step 1: Data Exploration
For this analysis, I have selected the aus_livestock dataset from the tsibbledata package, which records meat production in Australia for human consumption. The dataset is structured by 8 states and 7 livestock types.
By performing a diagnostic check using a cross-tabulatio, I discovered that some combinations contain structural zeros. From a statistical perspective, if an entire time series contains no data or zero variance, it will result in a singular covariance matrix. This leads to a numerical error during the MinT (Minimum Trace) reconciliation process. Therefore, identifying and handling these zeros is a critical first step.
# Check original distribution
table(tsibbledata::aus_livestock$State, tsibbledata::aus_livestock$Animal)
Bulls, bullocks and steers Calves
Australian Capital Territory 510 558
New South Wales 510 558
Northern Territory 486 558
Queensland 510 558
South Australia 510 558
Tasmania 510 558
Victoria 510 558
Western Australia 510 558
Cattle (excl. calves) Cows and heifers Lambs
Australian Capital Territory 558 0 558
New South Wales 558 510 558
Northern Territory 558 486 0
Queensland 558 510 558
South Australia 558 510 558
Tasmania 558 510 558
Victoria 558 510 558
Western Australia 558 510 558
Pigs Sheep
Australian Capital Territory 558 558
New South Wales 558 558
Northern Territory 558 558
Queensland 558 558
South Australia 558 558
Tasmania 558 558
Victoria 558 558
Western Australia 558 558Step 2: Data Cleaning & GTS Setup
In this stage, I limited the study period to start from January 2000 to ensure data consistency. The last two years, 2017 and 2018, are reserved as the test set to evaluate the model’s forecasting performance.
During the data cleaning process, I identified that 9 out of the original 56 combinations (8 states × 7 livestock types) contained no historical data at all within our chosen timeframe. Specifically: Australian Capital Territory (ACT): All 7 livestock categories are zero (as there are no active abattoirs in this territory). Northern Territory (NT): No records exist for Sheep and Lambs.
By excluding these structural zeros, the total number of series was reduced from 56 to 47. Filtering out these inactive series is essential to maintain numerical stability during the reconciliation process and to ensure the models focus on active production trends.
# Fill gaps and align start time
livestock_all <- tsibbledata::aus_livestock %>%
fill_gaps(Count = 0) %>%
filter_index("2000 Jan" ~ .)
# Calculate variance and flag status for all combinations in one go
analysis <- livestock_all %>%
filter_index(. ~ "2016 Dec") %>%
as_tibble() %>%
group_by(State, Animal) %>%
summarise(v = var(Count), total_count = sum(Count), .groups = "drop") %>%
mutate(Status = if_else(!is.na(v) & v > 0, "Valid", "Excluded"))
# Separate valid keys and excluded keys
valid_keys <- analysis %>% filter(Status == "Valid")
excluded_list <- analysis %>% filter(Status == "Excluded")
# Show the excluded list for your report
print(excluded_list)# A tibble: 7 × 5
State Animal v total_count Status
<fct> <fct> <dbl> <dbl> <chr>
1 Australian Capital Territory Bulls, bullocks and ste… 0 0 Exclu…
2 Australian Capital Territory Calves 0 0 Exclu…
3 Australian Capital Territory Cattle (excl. calves) 0 0 Exclu…
4 Australian Capital Territory Lambs 0 0 Exclu…
5 Australian Capital Territory Pigs 0 0 Exclu…
6 Australian Capital Territory Sheep 0 0 Exclu…
7 Northern Territory Sheep 0 0 Exclu…# Build the final GTS using only valid combinations
livestock_gts <- livestock_all %>%
semi_join(valid_keys, by = c("State", "Animal")) %>%
aggregate_key(Animal * State, Count = sum(Count))Step 3: Modeling & Forecasting
Since this study utilizes a Grouped Time Series (GTS) framework, the fitted models primarily employ the Bottom-up and MinT (Minimum Trace) reconciliation methods.
Methods such as Top-down and Middle-out are specifically designed for Hierarchical Time Series (HTS), where levels follow a strict, unique nested structure (parent-child relationship). In a GTS context, since factors like state and livestock type are crossed rather than nested, applying Top-down or Middle-out would require an arbitrary decision on which factor to prioritize. To maintain statistical integrity and account for the multi-dimensional correlations without forcing a hierarchy, I have focused on Bottom-up (aggregating base signals) and the optimal MinT approach (utilizing shrinkage estimation for error covariance).
# Split into Train (2000-2016) and Test (2017-2018)
train <- livestock_gts %>% filter_index(. ~ "2016 Dec")
test <- livestock_gts %>% filter_index("2017 Jan" ~ "2018 Dec")
# Build ETS models and reconcile with BU and MinT (Shrinkage)
fit <- train %>%
model(ETS = ETS(Count)) %>%
reconcile(
BU = bottom_up(ETS),
MINT = min_trace(ETS, method = "mint_shrink")
)
# Forecast for the year 2018
fc <- fit %>% forecast(h = "2 year")Step 4: Visualization
It is difficult to visually distinguish the fitting differences among the three models as they all track the historical data closely. While the models successfully capture the overall trend and strong seasonality of the training period, they tend to struggle with accurately reflecting large, sudden fluctuations in the actual production counts.
augment(fit) %>%
filter(is_aggregated(Animal), is_aggregated(State)) %>%
as_tibble() %>%
ggplot(aes(x = Month, y = .fitted, color = .model)) +
geom_line(aes(y = Count), color = "gray80") +
geom_line() +
facet_wrap(~ .model, ncol = 1) +
theme_minimal() + theme(legend.position = "none") +
labs(title = "Total Level: Model Fit (2000-2016)", y = "Count"
)
Regarding the forecast intervals, the Bottom-Up (BU) model appears to be the most “stable” or precise, as it produces the narrowest confidence intervals. However, this precision comes with a trade-off: because the BU intervals are so tight, the actual values frequently fall outside the 95% range during volatile periods. In contrast, the ETS (Base) model provides wider confidence intervals that more successfully encompass the actual data points when fluctuations occur.
# Forecast with Confidence Interval
fc %>%
filter(is_aggregated(Animal), is_aggregated(State)) %>%
autoplot(livestock_gts %>% filter_index("2015 Jan" ~ .), level = 95) +
facet_wrap(~ .model, ncol = 1) +
theme_minimal() + theme(legend.position = "bottom") +
labs(title = "Total Level: 2-Year Forecast (2017-2018)", y = "Count")
The accuracy table confirms that BU is the best-performing model for this dataset, followed by ETS and MinT.
accuracy_total <- fc %>%
filter(is_aggregated(Animal), is_aggregated(State)) %>%
accuracy(livestock_gts) %>%
select(.model, RMSE, MAE, MASE) %>%
arrange(RMSE)
accuracy_total