thfrec {FoReco} | R Documentation |
Forecast reconciliation through temporal hierarchies (temporal reconciliation)
Description
Forecast reconciliation of one time series through temporal hierarchies (Athanasopoulos et al., 2017). The reconciled forecasts are calculated either through a projection approach (Byron, 1978), or the equivalent structural approach by Hyndman et al. (2011). Moreover, the classic bottom-up approach is available.
Usage
thfrec(basef, m, comb, res, mse = TRUE, corpcor = FALSE,
type = "M", sol = "direct", keep = "list", v = NULL, nn = FALSE,
nn_type = "osqp", settings = list(), bounds = NULL, Omega = NULL)
Arguments
basef |
(\(h(k^\ast + m) \times 1\)) vector of base forecasts to be reconciled, containing the forecasts at all the needed temporal frequencies ordered as [lowest_freq' ... highest_freq']'. |
m |
Highest available sampling frequency per seasonal cycle (max. order of temporal aggregation, \(m\)), or a subset of \(p\) factors of \(m\). |
comb |
Type of the reconciliation. Except for bottom up, all other options correspond to a different (\((k^\ast + m) \times (k^\ast + m)\)) covariance matrix, \(k^\ast\) is the sum of (\(p-1\)) factors of \(m\) (excluding \(m\)):
|
res |
vector containing the in-sample residuals at all the temporal
frequencies ordered as |
mse |
Logical value: |
corpcor |
Logical value: |
type |
Approach used to compute the reconciled forecasts: |
sol |
Solution technique for the reconciliation problem: either
|
keep |
Return a list object of the reconciled forecasts at all levels (if keep = "list") or only the reconciled forecasts matrix (if keep = "recf"). |
v |
vector index of the fixed base forecast (\(\mbox{min}(v) > 0\) and \(\mbox{max}(v) < (k^\ast + m)\)). |
nn |
Logical value: |
nn_type |
"osqp" (default), "KAnn" (only |
settings |
Settings for osqp (object |
bounds |
(\((k^\ast + m) \times 2\)) matrix with temporal bounds: the first column is the lower bound, and the second column is the upper bound. |
Omega |
This option permits to directly enter the covariance matrix:
|
Details
Let \(m\) be the highest available sampling frequency per seasonal cycle, and denote \({\cal K} = \left\lbrace k_m, k_{p-1}, \ldots, k_{2}, k_1\right\rbrace\) the \(p\) factors of \(m\), in descending order, where \(k_p=m\), and \(k_1=1\). Define \(\mathbf{K}\) the \(\left( k^\ast \times m\right)\) temporal aggregation matrix converting the high-frequency observations into lower-frequency (temporally aggregated) ones: \[ \mathbf{K} = \left[\begin{array}{c} \mathbf{1}_m' \cr \mathbf{I}_{\frac{m}{k_{p-1}}} \otimes \mathbf{1}_{k_{p-1}}' \cr \vdots \cr \mathbf{I}_{\frac{m}{k_{2}}} \otimes \mathbf{1}_{k_{2}}' \cr \end{array}\right].\] Denote \(\mathbf{R} = \left[\begin{array}{c} \mathbf{K} \cr \mathbf{I}_m \end{array}\right]\) the \(\left[(k^\ast+m) \times m \right]\) temporal summing matrix, and \(\mathbf{Z}' = \left[ \mathbf{I}_{k^\ast} \; -\mathbf{K} \right]\) the zero constraints kernel matrix.
Suppose we have the \(\left[(k^\ast+m) \times 1\right]\) vector
\(\widehat{\mathbf{y}}\) of unbiased base forecasts for the
\(p\) temporal aggregates of a single time series \(Y\)
within a complete time cycle, i.e. at the forecast horizon \(h=1\)
for the lowest (most aggregated) time frequency. If the base forecasts
have been independently obtained, generally they do not fulfill the
temporal aggregation constraints, i.e. \(\mathbf{Z}'
\widehat{\mathbf{y}} \ne \mathbf{0}_{(k^\ast \times 1)}\).
By adapting the general point forecast reconciliation according to
the projection approach (type = "M"
),
the vector of temporally reconciled forecasts
is given by:
\[\widetilde{\mathbf{y}} = \widehat{\mathbf{y}} -
\mathbf{\Omega}\mathbf{Z}\left(\mathbf{Z}'\mathbf{\Omega}
\mathbf{Z}\right)^{-1}\mathbf{Z}'\widehat{\mathbf{y}},\]
where \(\mathbf{\Omega}\) is a \(\left[(k^\ast+m)
\times (k^\ast+m)\right]\) p.d. matrix, assumed known. The alternative
equivalent solution (type = "S"
) following the
structural reconciliation approach by Athanasopoulos et al. (2017) is given by:
\[\widetilde{\mathbf{y}} = \mathbf{R}\left(\mathbf{R}'
\mathbf{\Omega}^{-1}\mathbf{R}\right)^{-1}\mathbf{R}'
\mathbf{\Omega}^{-1}\widehat{\mathbf{y}}.\]
Bounds on the reconciled forecasts
When the reconciliation makes use of the optimization package osqp,
the user may impose bounds on the reconciled forecasts.
The parameter bounds
permits to consider lower (\(\mathbf{a}\)) and
upper (\(\mathbf{b}\)) bounds like \(\mathbf{a} \leq
\widetilde{\mathbf{y}} \leq \mathbf{b}\) such that:
\[ \begin{array}{c}
a_1 \leq \widetilde{y}_1 \leq b_1 \cr
\dots \cr
a_{(k^\ast + m)} \leq \widetilde{y}_{(k^\ast + m)} \leq b_{(k^\ast + m)} \cr
\end{array} \Rightarrow
\mbox{bounds} = [\mathbf{a} \; \mathbf{b}] =
\left[\begin{array}{cc}
a_1 & b_1 \cr
\vdots & \vdots \cr
a_{(k^\ast + m)} & b_{(k^\ast + m)} \cr
\end{array}\right],\]
where \(a_i \in [- \infty, + \infty]\) and \(b_i \in [- \infty, + \infty]\).
If \(y_i\) is unbounded, the \(i\)-th row of bounds
would be equal
to c(-Inf, +Inf)
.
Notice that if the bounds
parameter is used, sol = "osqp"
must be used.
This is not true in the case of non-negativity constraints:
-
sol = "direct"
: first the base forecasts are reconciled without non-negativity constraints, then, if negative reconciled values are present, the"osqp"
solver is used; -
sol = "osqp"
: the base forecasts are reconciled using the"osqp"
solver.
In this case it is not necessary to build a matrix containing
the bounds, and it is sufficient to set nn = "TRUE"
.
Non-negative reconciled forecasts may be obtained by setting nn_type
alternatively as:
-
nn_type = "sntz"
("set-negative-to-zero") -
nn_type = "osqp"
(Stellato et al., 2020)
Value
If the parameter keep
is equal to "recf"
, then the function
returns only the (\(h(k^\ast + m) \times 1\)) reconciled forecasts vector, otherwise (keep="all"
)
it returns a list that mainly depends on what type of representation (type
)
and solution technique (sol
) have been used:
recf |
(\(h(k^\ast + m) \times 1\)) reconciled forecasts vector, \(\widetilde{\mathbf{y}}\). |
Omega |
Covariance matrix used for forecast reconciliation, \(\mathbf{\Omega}\). |
nn_check |
Number of negative values (if zero, there are no values below zero). |
rec_check |
Logical value: |
varf (type="direct") |
(\((k^\ast + m) \times 1\)) reconciled forecasts variance vector for \(h=1\), \(\mbox{diag}(\mathbf{MW}\)). |
M (type="direct") |
Projection matrix, \(\mathbf{M}\) (projection approach). |
G (type="S" and type="direct") |
Projection matrix, \(\mathbf{G}\) (structural approach, \(\mathbf{M}=\mathbf{R}\mathbf{G}\)). |
S (type="S" and type="direct") |
Temporal summing matrix, \(\mathbf{R}\). |
info (type="osqp") |
matrix with information in columns
for each forecast horizon \(h\) (rows): run time ( |
Only if comb = "bu"
, the function returns recf
, R
and M
.
References
Athanasopoulos, G., Hyndman, R.J., Kourentzes, N., Petropoulos, F. (2017), Forecasting with Temporal Hierarchies, European Journal of Operational Research, 262, 1, 60-74.
Byron, R.P. (1978), The estimation of large social accounts matrices, Journal of the Royal Statistical Society A, 141, 3, 359-367.
Di Fonzo, T., and Girolimetto, D. (2023), Cross-temporal forecast reconciliation: Optimal combination method and heuristic alternatives, International Journal of Forecasting, 39(1), 39-57.
Hyndman, R.J., Ahmed, R.A., Athanasopoulos, G., Shang, H.L. (2011), Optimal combination forecasts for hierarchical time series, Computational Statistics & Data Analysis, 55, 9, 2579-2589.
Nystrup, P., Lindström, E., Pinson, P., Madsen, H. (2020), Temporal hierarchies with autocorrelation for load forecasting, European Journal of Operational Research, 280, 1, 876-888.
Schäfer, J.L., Opgen-Rhein, R., Zuber, V., Ahdesmaki, M., Duarte Silva, A.P., Strimmer, K. (2017), Package ‘corpcor’, R package version 1.6.9 (April 1, 2017), https://CRAN.R-project.org/package= corpcor.
Schäfer, J.L., Strimmer, K. (2005), A Shrinkage Approach to Large-Scale Covariance Matrix Estimation and Implications for Functional Genomics, Statistical Applications in Genetics and Molecular Biology, 4, 1.
Stellato, B., Banjac, G., Goulart, P., Bemporad, A., Boyd, S. (2020). OSQP: An Operator Splitting Solver for Quadratic Programs, Mathematical Programming Computation, 12, 4, 637-672.
Stellato, B., Banjac, G., Goulart, P., Boyd, S., Anderson, E. (2019), OSQP: Quadratic Programming Solver using the 'OSQP' Library, R package version 0.6.0.3 (October 10, 2019), https://CRAN.R-project.org/package=osqp.
See Also
Other reconciliation procedures:
cstrec()
,
ctbu()
,
htsrec()
,
iterec()
,
lccrec()
,
octrec()
,
tcsrec()
,
tdrec()
Examples
data(FoReco_data)
# top ts base forecasts ([lowest_freq' ... highest_freq']')
topbase <- FoReco_data$base[1, ]
# top ts residuals ([lowest_freq' ... highest_freq']')
topres <- FoReco_data$res[1, ]
obj <- thfrec(topbase, m = 12, comb = "acov", res = topres)