Discussed in #80

To simplify and reuse reconciler's code, it is convenient to declare a reconciler class with common methods like:

• fit
• predict
• fit_predict
• sample

In its current version some methods like HierarchicalForecast.reconcile need the Y_hat_df pd.DataFrame input to be indexed by 'unique_id' to operate correctly.

This makes running the methods challenging, in particular since StatsForecast.predict output is not indexed by 'unique_id'.
I suggest to make the needed set_index('unique_id') operation implicit in the reconcile method.

The example (How to use) in README.md loads data with

Y_df, S, tags = HierarchicalData.load('./data', 'TourismLarge')

However the HierarchicalData.load method returns a tuple of only 2, causing not enough values to unpack. (source here)

Python error message:

Traceback (most recent call last):
  File "hierarchicalforecast/example.py", line 11, in <module>
    Y_df, S, tags = HierarchicalData.load('./data', 'TourismLarge')
ValueError: not enough values to unpack (expected 3, got 2)

It would also be great if more examples can be provided.

Thanks!

A user reported the following error,

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
/var/folders/01/4mysj8cx1bjg_w8rbw097jwr0000gp/T/ipykernel_22347/4254128446.py in <module>
      7 hrec = HierarchicalReconciliation(reconcilers=reconcilers)
      8 #Y_rec_df = hrec.reconcile(Y_hat_df, Y_df_train, S, tags)
----> 9 Y_rec_df = hrec.reconcile(Y_hat_df=Y_hat_df, Y_df=Y_df_train, S=S, tags=tags)

~/opt/anaconda3/envs/ts_recon/lib/python3.9/site-packages/hierarchicalforecast/core.py in reconcile(self, Y_hat_df, S, tags, Y_df, level, bootstrap)
    148                 kwargs = {key: common_vals[key] for key in kwargs}
    149                 fcsts_model = reconcile_fn(y_hat=y_hat_model, **kwargs)
--> 150                 fcsts[f'{model_name}/{reconcile_fn_name}'] = fcsts_model['mean'].flatten()
    151                 if (pi and has_level and level is not None) or (bootstrap and level is not None):
    152                     for lv in level:

~/opt/anaconda3/envs/ts_recon/lib/python3.9/site-packages/pandas/core/frame.py in __setitem__(self, key, value)
   3653         else:
   3654             # set column
-> 3655             self._set_item(key, value)
   3656 
   3657     def _setitem_slice(self, key: slice, value):

~/opt/anaconda3/envs/ts_recon/lib/python3.9/site-packages/pandas/core/frame.py in _set_item(self, key, value)
   3830         ensure homogeneity.
   3831         """
-> 3832         value = self._sanitize_column(value)
   3833 
   3834         if (

~/opt/anaconda3/envs/ts_recon/lib/python3.9/site-packages/pandas/core/frame.py in _sanitize_column(self, value)
   4533 
   4534         if is_list_like(value):
-> 4535             com.require_length_match(value, self.index)
   4536         return sanitize_array(value, self.index, copy=True, allow_2d=True)
   4537 

~/opt/anaconda3/envs/ts_recon/lib/python3.9/site-packages/pandas/core/common.py in require_length_match(data, index)
    555     """
    556     if len(data) != len(index):
--> 557         raise ValueError(
    558             "Length of values "
    559             f"({len(data)}) "

ValueError: Length of values (4704) does not match length of index (2688)

The error is related to a difference between the number of base forecasts expected by S.

• Change back the crown icon to "Hierarchical Crown Forecast", or better use favicon to change the world image from the explorer tab towards the library icon.
• Change PERMBU description in the nbs/index.ipynb and add links to probabilistic methods in the description.

Forecasts are calculated correctly, but if the forecasts dataframe Y_h is not sorted according to S, forecasts are assigned to the wrong place:

If a repository contributor adds a comment #self-assign to any issue, that will assign the issue to the user. This will unburden repo-maintainers from having to assign issues to people.

Trying to implement on forecasts that i already have (both aggregations and low level)

I made my own summing matrix but when running the tutorial to see what the y_hat_df was formatted like, ran into this error:

"""
TypeError Traceback (most recent call last)
in ()
2 models=[(ets, 4, 'ZZA')],
3 freq='QS', n_jobs=-1)
----> 4 Y_hat_df = fcst.forecast(h=8, fitted=True)
5 Y_fitted_df = fcst.forecast_fitted_values()

3 frames
/usr/lib/python3.7/multiprocessing/pool.py in get(self, timeout)
655 return self._value
656 else:
--> 657 raise self._value
658
659 def _set(self, i, obj):

TypeError: seasonal_decompose() got an unexpected keyword argument 'period'

"""

Which is odd bc I use seasonal decompose and know it has a period argument..

This happens on all tutorials or examples when opening in collab..

Do you have an example of implementation with forecasts/tags you already have or could you please fix the examples?

For now, the evaluation is performed for each hierarchy considering all observations regardless of the time series. In metrics that are not "distributive" for several series, such as rmse, the performance can be confusing.

The following line is not general. Not all datasets contain the unique_id "['Total']".

The total unique_id could be recovered using tags.

Hi the google colab notebook is not available.. Kindly update it

The 'Nixtla/hierarchicalforecast' repository doesn't contain the 'examples/TourismSmall.ipynb' path in 'main'.

PERMBU throws an error when n_samples is different than y_insample.shape[1].

In the presence of non stationarity ERM might degrade in performance.
It would be a good idea to add a test, and a warning for the user if the series are non-stationary.

From the search of other available reconciliation methods. It seems that we are still lacking the ones from
"George Athanasopoulos, Rob J. Hyndmanm Roman Ahmed, and Lin Shang. Optimal combination forecasts for hierarchical".

I would like to confirm.

The section Non-negative MinTrace is not as informative as Non-Negative Hierarchical Reconciliation, especially for potential users unfamiliar with the MinTrace method.

Kindly update the example in the "How to use" section at the PyPI website, as well as the link to the Colab Notebook.

I've been stuck in that non-working example, but here in the repo, I found the correct one. You have developed a fantastic package, but folks like me who first face the PyPI documentation might find it not working while it is (with updated code).

The main issue I found at PyPI is using auto_arima and naive instead of AutoARIMA and Naive.

Thanks for your contribution =D

See conda-forge/hierarchicalforecast-feedstock#4.

A clever way to simplify lower level reconciliation functions is to assume 'hierarchically ordered series'.
It can be achieved with the following lines at the beginning of the reconcile method on the core.HierarchicalForecast with these lines:

Y_df.unique_id = Y_df.unique_id.astype('category')
Y_df.unique_id = Y_df.unique_id.cat.set_categories(S_df.index)

If we decide to operate like this, we should add in the docstrings of the methods that are order dependent the assumption.
Currently the PERMBU method operates like this.

These two issues are closely related: Challenging HierarchicalForecast method input's indexing.

Users have reported issues trying to run MinTrace(method=='mint_shrink').

I suspect that it is related to datasets with time series of distinct lengths (in our examples all time series have the same length and the method works fine). In particular, I think that line 441 below is not working correctly. It should ignore the NAS in the residuals matrix generated by the difference in lengths, but perhaps it is not ignoring them.

From Linear-GitHub Sync

Both Bootstrap and Normality probabilistic reconcilers require inputs that depend on a previous run of mean hierarchical reconciliation methods.

Currently the instantiation of probabilistic reconcilers in HierarchicalForecast.reconcile is complicated due to the effort of sending samplers as part of the reconcile_fn arguments.

I believe that separating the responsabilities between reconcile_fn and prob_reconcile_fn, used after the use of the reconcile_fn can improve code's readability and will help to simplify the instantiation of the samplers.

##########################Getting this error#####################

RemoteTraceback Traceback (most recent call last)
RemoteTraceback:
"""
Traceback (most recent call last):
File "/usr/lib/python3.7/multiprocessing/pool.py", line 121, in worker
result = (True, func(*args, **kwds))
File "/usr/local/lib/python3.7/dist-packages/statsforecast/core.py", line 154, in forecast
raise error
File "/usr/local/lib/python3.7/dist-packages/statsforecast/core.py", line 149, in forecast
res_i = model.forecast(h=h, y=y_train, X=X_train, X_future=X_f, fitted=fitted, **kwargs)
File "/usr/local/lib/python3.7/dist-packages/statsforecast/models.py", line 481, in forecast
mod = ets_f(y, m=self.season_length, model=self.model)
File "/usr/local/lib/python3.7/dist-packages/statsforecast/ets.py", line 968, in ets_f
raise NotImplementedError('tiny datasets')
NotImplementedError: tiny datasets
"""

The above exception was the direct cause of the following exception:

NotImplementedError Traceback (most recent call last)
in
2 models=[ETS(season_length=4, model='ZAA')],
3 freq='M', n_jobs=-1)
----> 4 Y_hat_df = fcst.forecast(h=4, fitted=True)
5 Y_fitted_df = fcst.forecast_fitted_values()

Some methods like TopDown, BottomUp, ERM benefit a lot from sparse matrix operations.

Luckily the numpy's and scipy.sparse's matrix multiplications can both be expressed with .dot method.

It would be convenient to make the HierarchicalForecast library agnostic to the method's sparsity by changing all the corresponding matrix multiplications.

It would be convenient to have prefitted ARIMA base forecasts for medium and large datasets on S3.

• We guarantee replicability.
• We can save a lot of time for potential users and other researchers.
• Circle CI tests that check time efficiency can greatly benefit from this dataset.

After going through the package and paper found that supports only point forecasts reconciliation. I would like to know if it supports probabilistic forecast reconciliation. I appreciate your quick response to my inquiry for information.

The Base predictions and Reconciliation sections don't have a complete explanation of the hierarchical reconciliation methodology. Ideas to improve those sections:

Base predictions:

• Define what we mean by "base predictions".
• Use a different model. By construction, the Naive model is coherent, so this model does not well motivate the use of hierarchical reconciliation. In particular, RandomWalkWithDrift could work.

Reconciliation

• Show that the base forecasts are not coherent, using plot_hierarchical_predictions_gap (fix #85 first).
• Explain the methods used for reconciliation. In this case, the explanation for BottomUp is missing.

Let's reevaluate the bootstrap boolean parameter from the core.HierarchicalForecast main class to prob_method that accepts strings to be able to use PERMBU, Bootstrap or MinTrace-normality.

We need to think of a way to make the methods compatible between one-another or to declare a new ProbabilisticReconciler class.

In the .evaluate of the HierarchicalEvaluation class, it sets h = len(Y_h.loc[Y_h.index[0]]) but when the look ahead period is only 1 period ahead, Y_h.loc[Y_h.index[0]] returns values rather than a dataframe.

The length of these values is then 2 because of the ds column and the forecast column. This then returns an error during the reshape of Y_test here: y_test_cats = Y_test.loc[cats, 'y'].values.reshape(-1, h) because it tries to reshape to 2 even though there is only 1 record.

Would changing h = len(Y_h.loc[Y_h.index[0]]) to h = len(Y_h.loc[[Y_h.index[0]]]) with the double brackets to ensure a dataframe is returned work?

Let me know if you need any additional detail or a reproducible example as this is my first issue statement

In some cases, utils.aggregate returns a different order of unique_id in Y_df and S. Even though core.HierarchicalReconciliation ensures the same order, it could reduce a good user experience.

Although probabilistic hierarchical forecasting is still a niche topic, there is a clear need to extend hierarchical forecasting available methods beyond point/mean forecasting.

This paper has an R implementation that can be improved, along with pointers to other methods we can extend towards their predictive distributions.

By recommendation of Shanika, we can change the current matrix:
$P_{MinT} = (SW^{−1} S)^{-1}S^{t} W^{−1}$

To avoid the double computation of the inverses that scale cubically.
Her paper contains the alternative that uses a single inverse with lower dimensionality.

Hello, when I call pip install hierarchicalforecast and then try to import the OptimalCombination class I get an import error.
I was able to fix this issue but calling pip install git+https://github.com/Nixtla/hierarchicalforecast to import the newest version of the repo, so I think the pip is not up to date with the latest version of the repo.

Perhaps the readme can be updated to include:

https://arxiv.org/abs/2207.03517

In the method core.HierarchicalForecast.reconcile
I suggest we change the name of the S dataframe in to S_df

• It homogeneizes inputs across libraries.
• It makes a distinction between method's internal matrix S.
• Suffix _df helps our user's intuition on the type of input expected.

It would be convenient to add continuous integration tests to monitor performance and fitting times of the reconciliation methods.
Automated tests will hasten the development and improvements of the code.

Some methods don't use the in-sample values of the time series (for example BottomUp, MinTrace(method='ols'), TopDown(method='forecast_proportions')). Therefore, the Y_df argument should be optional to avoid overcomplexity using the reconcile method and should be requested only when necessary.

Dear colleagues,

I am working with a problem where I do need to predict the performance of a financial portfolio formed by different assets. As of today I am approaching this as a combination of a classification and regression problem but now I have been task to do this by all of the portfolios that our firm has.

I can replicate the same method as before as I am predicting at asset level but I have been thinking to frame this as a hierarchical forecasting problem.

The only thing is that the number of portfolios grow over time and the number of assests is different from portfolio to portfolio and I am bit unsure if there is a way to reconcialiate forecasts if the hierachy size is changing over time without retraining the models.

Looking to the notation that R.Hyndman uses:

It seems that the matrixes they need to be mxn compatible.

Is this the right interpreation or is there any way to reforecast at hierarchical level without retraining again the models.

The fable/hts replication 2 notebook, has different base predictions compared to the
Australian Prison Population example from this hts example.

Two ideas come to my mind:

• Different data sources
• Different ETS model execution.

I followed the tutorial, and everything worked until I reached the following part:

reconcilers = [ BottomUp(), TopDown(method='forecast_proportions'), MiddleOut(middle_level='filial/categoria', top_down_method='forecast_proportions')] hrec = HierarchicalReconciliation(reconcilers=reconcilers) Y_rec_df = hrec.reconcile(Y_hat_df=Y_hat_df, Y_df=Y_train_df, S=S, tags=tags)

It works fine in small datasets, but I'm trying to use this in a dataset with the following structure:

Filial (36 unique)
----Categoria (50 unique)
--------Master(~30k unique)

Each series have 36 time steps.

The code is currently running the reconcile() function, it is stuck in this part for the last 2 hours. It would be nice if there was an indication of how long it will probably take or at least some sort of verbose showing of what is going on.

My computer has an i5-12500H, 16Gb of RAM, and an SSD for storage.

Discussed in #47

Kindly ensure you follow best academic citation practice in your academic papers, e.g., the accompanying paper here:
https://arxiv.org/abs/2207.03517

I'm aware that in industry it might not be common or expected to comment on prior art (or even considered counterproductive for marketing reasons), and that Nixtla is primarily a commercial venture. But, since you are publishing in the arXiv you are at least making the claim of adhering to basic scientific standards, so I think you should at least give it a try.

A proper "literature review" or "prior art" section is the minimum in any scientific paper, and that means commenting on the specific context of your contribution, not just providing a generic list of packages in an appendix.

As past contributors to sktime, you are probably aware of its hierarchical framework functionality, and the fact that the algorithms you present are also implemented there. See, e.g., this conference presentation from April:
https://github.com/sktime/sktime-tutorial-pydata-berlin-2022

Perhaps more important even are the earlier python packages when it comes to hierarchical forecasting, which have developed pertinent designs, and which you don't even cite, e.g.,:
https://github.com/carlomazzaferro/scikit-hts, FYI @carlomazzaferro
https://github.com/AngelPone/pyhts, FYI @AngelPone

And there is, of course, even more in R.

You can't claim to be a "reference" (in the title of your paper!) without making reference to what has come before. Pun intended.

From your mission statement: "We intend to continue maintaining and increasing the repository, promoting collaboration across the forecasting community."

I do hope to see that, giving due credit in your scientific papers would be a good start.

Let me know if you have any questions on best scientific practice.

It would be good to specify in the code the type of covariance estimator used by the Normality class.
Further in the future it would be good to add the capability to switch types of covariance estimator to the Normality reconciliator with a covariance_type input.

That assumption is too restrictive for the reconcile method.

It would be beneficial for complex models like MinT-shr and ERM-lasso to inherit the BaseEstimator sklearn class to help with the automated search of hyperparameters.

When I open the documentation and click on the source code links, it leads to a 404 webpage error on github.

Top-down methods don't work with non-strictly hierarchies, (for example groups).

When a time series only has one observation, the HierarchicalPlot.plot_hierarchically_linked_series fails. The error arises because of the following line,

In that case, Y_df.loc[series] returns a pd.Series object instead of a pd.DataFrame object and therefore the index cannot be set.

Some methods and functionality (bootstrapped prediction intervals) need the fitted values of the models. But the following iteration omits to request them.

The above leads to unintelligible bugs.