<figure><img src="image_1.png" /></figure>

HAL
open science

# Time series forecasting using Machine Learning: some results of two research projects

Gerardo Rubino

► To cite this version:

Gerardo Rubino. Time series forecasting using Machine Learning: some results of two research projects. CORE 2024 - 23th International Conference on Computer Science Congress, IPN, Mexico, Sep 2024, Mexico DF, Mexico. hal-04756220

HAL Id: hal-04756220

https://inria.hal.science/hal-04756220v1

Submitted on 28 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

Public Domain
# Time series forecasting using Machine Learning:
some results of two research projects

G. Rubino, INRIA, France

CORE 2024, Mexico DF, sept. 2024

<figure><img src="image_2.png" /></figure>
# The general problem

* The problem is: we have all the values (for instance, real numbers) of a time series x, in a time interval $[T_0, T]$, and we want to forecast some information related to the future of this signal (that is, related to the numbers x(t), for $t > T$).

* This a classic **statistical** problem.

* The emergence of **Machine Learning** (and its current evolution) made it possible to attack it using these new tools. This is the topic of this talk.

* We will be concerned here, inside Machine Learning, with Supervised Learning.

* We just finished a project with people from Uruguay, Argentina, Chile and France, *ClimateDL*, where we looked at this type of problem for a specific signal, *temperature*. I will provide some descriptions about the project and about some of the obtained results.

* Inside *ClimateDL*, we developed a subproject about a specific methodology for these problems, which I will present first.

<figure><img src="image_3.png" /></figure>
Outline

I/ — The EvoESN tool (with S. Basterrech, Copenaghue)

II/ — ClimateDL: some (other) results (mainly with P. Rodríguez, Montevideo)

III/ — Conclusions

<figure><img src="image_4.png" /></figure>
# Outline

I/ - EvoESN

II/ — ClimateDL: some (other) results

III/ – Some conclusions and perspectives

<figure><img src="image_5.png" /></figure>
An independent subproject

* EvoESN = Evolutionary Echo State Network

* This is the main output of a subproject of ClimateDL that emerged as an independent methodological contribution of ClimateDL.

* The problem was the exploration of a family of tools designed for the general problem described at the beginning and then, to try to improve some known weaknesses of the approach.

<figure><img src="image_6.png" /></figure>
# Models with temporal information

* The important feature when using Machine Learning to forecast some aspect of the future of a time series is the way we will handle the time dimension.

* Some general remarks follow, concerning the behavior of fundamental neural architectures face to this problem.

* Basically, a layered Neural Network (NN) is an easy to train tool but it performs poorly when used for time series forecasting.

* At the opposite side of the spectrum, a recurrent NN captures more naturally and more efficiently the time dimension information. However, it is harder to train, leading to pretty long training times.

* The idea of the approach I will discuss here (called Reservoir Computing, RC) is the design of a model trying to combine the positive sides of both architectures and to avoid the negative ones.

<figure><img src="image_7.png" /></figure>
# Layered Neural Networks

Directed graph without circuits:

* Three types of neuron layers: *input*, *hidden* and *output* layer.

* The input layer contains neurons whose inputs are the input pattern data.

* The hidden layers receive inputs from each unit in the preceding layers.

* The neurons in the output layer are the only neurons that produce the output values of the whole network.

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>
Recurrent connectivities

Directed graph with circuits:

* Class of network where there is **at least one circuit** in the graph.

* Connections among hidden units are allowed with a time delay.

* Fundamental difference with feedforward networks, a network with recurrences defines a dynamical system.

* This creates an **internal state** which makes the method a discrete time state-space model.

<figure><img src="image_10.png" /></figure>
# Reservoir Computing (RC)

* Recently (the 2000s), the Reservoir Computing (RC) approach was proposed to overcome the important limitation of training procedures of RNNs while introducing no significant disadvantages for time series processing applications.

* RC idea:

* (i) a NN where the input neurons are connected to a much larger set of neurons having circuits (the “reservoir”), which in turn connects to a set (usually small) of output neurons (the “readout”), together with

* (ii) the only part of the weights that learn is the set of weights of the connections between the reservoir and the output neurons.

* In some cases, RC outperforms classic method for forecasting the behavior of time series.

* Nowadays other possibilities have appeared to better train standard recurrent structures (not covered here).

<figure><img src="image_11.png" /></figure>
ESN: Echo State Network

* Consider a neural network with *d* input neurons, a reservoir with *N* ≫ *d* units, and a set of *o* output neurons. The input at time *t* is vector **u**(*t*) ∈ ℝ<sup>*d*</sup>, and the output at time *t* is vector **y**(*t*) ∈ ℝ<sup>*o*</sup>. At time *t*, the state of the network is a vector **x**(*t*) ∈ ℝ<sup>*N*</sup>. The time variable *t* is an integer.

* The weights from input to reservoir are grouped into the *d* × *N* matrix **W**<sup>in</sup>; The weights inside the reservoir are grouped into the *N* × *N* square matrix **W**<sup>res</sup>; the weights from the reservoir to the outputs are grouped into the *N* × *o* matrix **W**<sup>out</sup>.

* The weights in matrices **W**<sup>in</sup> and **W**<sup>res</sup> are initialized and then, **frozen**. The weights in **W**<sup>out</sup> evolve during the learning phase.

<figure><img src="image_12.png" /></figure>
Basically, a typical implementation of an ESN means, formally, choosing

• $x(t) = f(\mathbf{W}^{\text{in}}\mathbf{u}(t) + \mathbf{W}^{\text{res}}x(t-1))$

• and $y(t) = \mathbf{W}^{\text{out}}x(t)$,

• where $f(x) = \tanh(x)$ coordinate-wise.

• For simplicity, bias terms are omitted (as usual, we replace them with dummy terms).

Learning means then to find a good matrix $\mathbf{W}^{\text{out}}$ using any possible approach (such as a gradient descent). Matrices $\mathbf{W}^{\text{in}}$ and $\mathbf{W}^{\text{res}}$ are randomly initialized and then, frozen.

<figure><img src="image_13.png" /></figure>
ESP

* Echo State Property (ESP): asymptotically, the network state $x(t)$ converges as $t \to \infty$, whatever the initialization of the frozen weights.

* The ESP holds under reasonable technical conditions, the most important one being the fact that the spectral radius of $W^{\text{res}}$ is < 1.

* Main question: even if theoretical results say that the behavior of the tool is **asymptotically** independent of the initial weights in the reservoir, what happens in practice?

* We explored this issue, and attacked the problem posed by the fact that, in practice, the initial values of the reservoir's weights matter.

<figure><img src="image_14.png" /></figure>
# Our proposal

* For technical reasons, supported by experimental data, we can take the $W^{\text{res}}$ matrix (by far, the largest of the 3) sparse. Asume that we will have $M \ll N^2$ entries a priori non-zero, that is, $N^2 - M$ frozen zeros in $W^{\text{res}}$.

* We store the $M$ non-zero values of the matrix in a vector $v$ and the position in the matrix of value $v[k]$ using another vector $p$, in value $p[k]:$ if $p[k] = (i, j)$, then $W^{\text{res}}[i, j] = v[k]$.

* In our method, we modify the main paradigm of RC by touching the coefficients of $W^{\text{res}}$, that is, the values of vector $v$. The idea is to improve the learning process by modifying those weights, in spite of the theoretical results saying that *in the limit*, the initial values don't matter.

* We first compress vector $v$ using a Fourier-like transformation. For specific technical reasons, we use the Discrete Cosine Transform (DCT), recalled in next slide.

<figure><img src="image_15.png" /></figure>
Discrete Cosine Transform (DCT)

* The DCT has properties that make it suitable for encoding large matrices.

* We use one of the most common versions, the normalized DCT of Type II: given an arbitrary sequence of H reals s = (s₀, s₁, ..., sH-1), its DCT is DCT(s) = S = (S₀, S₁, ..., S_H-1), where

* S₀ = $\sqrt{\frac{1}{H}} \sum_{h=0}^{H-1} s_h$,

* and for $1 \le l \le H-1$, $S_l = \sqrt{\frac{2}{H}} \sum_{h=0}^{H-1} s_h \cos(\frac{\pi l(2h+1)}{2H})$.

* This transformation is bijective; its inverse is denoted IDCT here, that is, with previous notation, IDCT(S) = s. We then have

$$s_l = \frac{S_0}{\sqrt{H}} + \sqrt{\frac{2}{H}} \sum_{h=1}^{H-1} S_h \cos\left(\frac{\pi h (2l + 1)}{2H}\right).$$

<figure><img src="image_16.png" /></figure>
Our method

* For some comfort in discussions, we usually call *chromosome* a vector such as DCT(v). We call *phenotype* a vector of weights such as v. We also call *frequency space* the universe of chromosomes.

* Assume that we store the first *C* elements of DCT(v) (the low frequency values) in a vector *α* having size *M*, completing the rest of *α* with zeros (*M* − *C* values); this last operation is called *padding*.

* Our method (called EvoESN for Evolutionary ESN) is described here:

* We start by randomly initializing a vector (a chromosome) with *C* components.

* Denote by $γ_0$ the result of padding this vector.

* Then, we build (the phenotype) $v_0 = \text{IDCT}(\gamma_0)$.

* We train our ESN with the **weights** in $v_0$.

* If the training doesn't pass the testing phase, we perform a "genetic modification" (context: optimization using genetic algorithms) of $γ_0$ producing a new chromosome $γ_1$ and we repeat the process.

* When the training iteration passes the test, the procedure is stopped.

<figure><img src="image_17.png" /></figure>
<figure><img src="image_18.png" /><figcaption>Visualization of the genotype-phenotype mapping of EvoESN.</figcaption></figure>

<figure><img src="image_19.png" /></figure>
Evaluation

* This technique was implemented as described here and compared favorably with standard implementations of ESNs on usual benchmarks (for instance, on chaotic systems of differential equations).

* We studied its behavior regarding the different parameters previously introduced (such as C/M, M/N²), of the number of iterations, etc. Part of them are mentioned in the paper we published on the tool (see next slide).

* The tool is now being implemented inside a package distributed by INRIA.

<figure><img src="image_20.png" /></figure>
* A preliminary version of our technique appears in

**Evolutionary Echo State Network: evolving reservoirs in the Fourier space**, Sebastián Basterrech and Gerardo Rubino, in IEEE International Joint Conference on Neural Networks (IJCNN'22), Padua, Italy, July 18-23, 2022.

* An improved method, with a more complete description and more extensive results illustrating our approach, was published in

**Evolutionary Echo State Network: A neuroevolutionary framework for time series prediction**, Sebastián Basterrech and Gerardo Rubino, in ASOC (Applied Soft Computing), No. 144, p. 110463, 2023.

<figure><img src="image_21.png" /></figure>
Illustrations

* See paper.

<figure><img src="image_22.png" /></figure>
# Outline

I/ — EvoESN

II/ — ClimateDL: some (other) results

III/ — Some conclusions and perspectives

<figure><img src="image_23.png" /></figure>
Some extensions to the original problem

First of all, a couple of generalizations of our basic problem which are relevant for
our project.

* 1st: instead of having a real-valued signal x, we consider the case where $x \in \mathbb{R}^d$ for some $d \ge 1$, a vector time series.

* In the temperature example, assume time is an integer, the time unit is a day, and $x_1(t)$ is the maximal temperature at day $t$.

* Then, $x_2(t)$ could be the cumulated precipitations at day $t$ (starting at some reference initial day), $x_3(t)$ a measure of atmospheric pressure, etc.

* 2nd: instead of consider a single signal x, sometimes we have

* a set of signals $x^{(1)}, \dots, x^{(K)}$, where $x^{(k)}$ corresponds, for instance, to geographical location $k$,

* possibly some information about some kind of inter-dependencies between the set of signals;

* in this type on context, we typically want to forecast some future information concerning all/some of them.

In our project, we consider both extensions.

<figure><img src="image_24.png" /></figure>
The ClimateDL problem

* We are interested in forecasting events related to extreme weather in a future period, and in a given large area. More specifically, we want to explore the problem of forecasting how hot will be next summer, for example. Our area is a large part of South America. (According to climatologists, from a physical viewpoint, forecasting how cold will be next winter is harder).

* Our data is (supposed to be) rich: in many geographical points in the area, we have more than 30 years of histories of dozens of parameters, plus many other a priori relevant data for the problem at hand.

* In the project, we consider the two previously described extensions of the general problem.

* ClimateDL is a joint research project between scientists from Uruguay, Argentina, Chile and France, composed of data scientists, climatologists, meteorologists, and mathematicians. The project ended in the first quarter of 2024. The talk concerns basically the work of part of the data scientists.

<figure><img src="image_25.png" /></figure>
In a nutshell

* Today, the operational sources of the type of information we consider (some characteristics of a future large period, say next summer) are a few (less than 20 centers in the world). Around high temperatures, they basically provide a kind of probability that next summer will be “warmer than usual”, “as hot as usual” or “cooler than usual”.

* Since we were going to rely on Machine Learning (ML) methods, at the beginning we decided to try to predict properties of the random variable “number of days of next summer where the average temperature is very high” (appropriately defined), following similar paths than in other (very different) application areas.

* All the methods we tried failed (and then, we understood why the available forecasts target that simple 3-valued qualitative output).

<figure><img src="image_26.png" /></figure>
* We moved to target that 3-valued variable, and took as a rule to use ML techniques that could run on a "standard machine" (while having a good level of RAM and a few GPU), hoping that the rich data we had could lead to something interesting.

* Keeping the same constraint of using only "modest" resources, we explored also some related issues (see below).

* Basically, we got state-of-the-art accuracy levels when using ML, in spite of using standard tools and "normal" hardware.

* The performance was, in several cases, better than that of statistical tools.

* We are happy with these results, but they confirm that these problems are really hard.

<figure><img src="image_27.png" /></figure>
About today's production systems

* There are less than 20 centers in the world publishing a forecasting report about next summer with the 3-states qualitative output (example: IRI, at Columbia).

* They all provide similar reports. Importantly, **they don't agree**. Then, some of them build then a new predictor by combining the outputs of all the centers.

* Some examples of key components in these forecasting centers:

* NMME (North American Multi-Model Ensemble) – several NPW (Numerical Weather Prediction) models running on large computing platforms, from several institutions

* CPT (Climat Predictability Tool) – statistical techniques, from Columbia

* SMN-CPT.NMME: a variant of CPT combining CPT with NWP models used inside NMME

* ...

* Using, say, 15 predictors in the 3-classes outputs, they can count a "score" for each class and publish a kind of "probability" for each one; the accuracy is obviously very low.

<figure><img src="image_28.png" /></figure>
Our results

Basically, we explored 4 problems in this project. One is the RC area (first part of the talk).

In the remaining 3, we attacked similar or related problems as the original one, in order to explore the general area, and we also looked at the mathematical foundations associated with using ML in exploring the future of time series.

Results:

1/ We studied a "meteorology" version where we looked at the target maximal *temperature of a day* and tried to forecast its value up to 10 days in the future. The result was published as follows:

**Maximal Temperature forecasting under spatio-temporal interrelations using Machine Learning**, Gonzalo Marco, María Eugenia Miranda, Pablo Rodríguez-Bocca and Gerardo Rubino, in ECML'23 (European Conference on Machine Learning and Principles and Practice of Knowledge Discovery in Databases), Track: 2nd Workshop on Machine Learning for Irregular Time Series: Advances in Generative Models, Global Models and Self-Supervised Learning (ML4ITS'23), Torino, September 22, 2023.

<figure><img src="image_29.png" /></figure>
2/ We also explored the target *mean temperature* over a future segment of, for example, 90 days. Result:

Seasonal average temperature forecast with the AutoGluonTS modern autoML tool, Diego Kiedanski, Pablo Rodríguez-Bocca and Gerardo Rubino, in ECML'23 (European Conference on Machine Learning and Principles and Practice of Knowledge Discovery in Databases), Track: MACLEAN: MACHine Learning for EArth ObservatioN (MACLEAN'23), Torino, September 18, 2023.

3/ The results concerning the original problem of ClimateDL (class of next summer) isn't submitted yet, but the research work is done. We are now writing. The authors are Soledad Collazo, Guillermo Pereira, Pablo Rodríguez-Bocca and Gerardo Rubino.

<figure><img src="image_30.png" /></figure>
Some details about our data

* An important part of the data concerns the Southern part of South America (Argentina, Chile, Uruguay, plus small regions in Paraguay and Brazil), and was produced by climatologists in Argentina, obtained from more than 200 weather measuring stations in the area (see next slide).

* To this set we added data coming from other public sources, including data concerning the situation at three oceans (Pacific, Atlantic and also Indian). In particular, we used ERA5 data here<fnref n="1" />.

* An important engineering work was done to clean the data, to compress other external data (specially data concerning the state of the oceans, including atmospheric associated parameters), etc.

* The result was a huge dataframe (Pandas) with more than 30 years of daily information, plus some preliminary pre-computations.

<figure><img src="image_31.png" /></figure>

<footnote n="1">¹ERA5 is the fifth ECMWF atmospheric project that provides hourly estimates of many atmospheric, land and oceanic climate variables (ECMWF is the European Centre for Medium-Range Weather Forecasts)</footnote>
The locations of the stations

<figure><img src="image_32.png" /></figure>

<figure><img src="image_33.png" /></figure>
About result 1/

We studied a “meteorology” version where we looked at the target *maximal temperature of a day* and tried to forecast its value up to 10 days in the future.

* Problem: forecast the max daily temperatures at a given geographical place, for the next *j* days, where $1 \le j \le 10$ (a meteorological question).

* Tools used: **XGBoost** and **Graph Wavenet**.

* XGBoost uses *decision trees*, it is very popular, with a nice records of good results. It needed to use an important engineering effort (*feature engineering*) to use the available spatial data.

* Graph Wavenet belongs to the Graph Neural Networks family. It handles the spatial interdependencies in the data in the usual ML black box way, and it needed almost no engineering effort to be used.

* Results: very similar with both tools. At the state-of-the-art level, but running on modest hardware. When *j* is close to 1, Graph Wavenet wins; when it gets close to 10, the contrary happens.

<figure><img src="image_34.png" /></figure>
A comment/hope

Consider the problem of predicting the close future (meteorology case).

* Here, extremely heavy platforms are the only operational sources of forecasts.

* In November 2022 appeared "Pangu-Weather: A 3D High-Resolution System for Fast and Accurate Global Weather Forecast", by 6 researchers from Huawei Cloud Computing (arXiv:2211.02556v).

* Using ML and a powerful computing platform, they reached (at least) the level of current production systems, using a very rich data set (source: ERA5), several Deep Learning procedures capable of predicting many basic parameters (temperature, wind, humidity, etc.), on a fine spatial grid. Some contributions on the ML side too.

* And in December 2022 appeared "GraphCast: Learning skillful medium-range global weather forecasting", by a large team (composed of 18 authors) from Google/Deepmind (arXiv:2212.12794v1).

* Similar characteristics (rich data, powerful hardware facilities), techniques based on the Graph Neural Network paradigm, and even better accuracy levels than previous reference. With also contributions in ML.

* The problem is different and, somehow, simpler, and still very hard (chaotic aspects are critical); these two breakthroughs are then promising.

<figure><img src="image_35.png" /></figure>
About result 2/

* Here, the metric studied is the *mean temperature* over a future period of any number of consecutive days (90 days for next summer). Easier than our main problem. We finally relied entirely on AutoGluon, an AutoML tool.

* **AutoML** frameworks primarily focus on automatic model/hyperparameter selection. A recent member of this family is **AutoGluon-Tabular**. It also offers the automatic combination of many models into a single output better than any of them, by means of a novel way of combining models.

* It also adds strategies for choosing the appropriate models in order to offer a much better use of the available computing time, feature engineering and data processing services, etc. It has been reported to be faster and much more accurate than its competitors, and it is very easy to use.

* We used an extension, **AutoGluonTS** (TS = Time Series), actually a generalization, that can process data coming from time series (plus external variables, etc.). It efficiently transforms the problem into a tabular one. Then, it passes the job to AutoGluon-Tabular. All this is done automatically.

<figure><img src="image_36.png" /></figure>
Some details about the results

* The basic available data is the min and max temperatures per day, at any of the 137 selected weather stations; as everyone else in the same situation, we “define” the average temperature at each day as their arithmetic mean at that day. We know this value for, say, 40 years, from 1981 to 2020. We used the 30 years of data 1981..2010 for training our tools.

* Our time variable *t* (in days) starts at 1, which corresponds to Jan 1, 1981. Denote by *year(t)* the year of *t* (*year(1)* = 1981 = *year(200)*, *year(400)* = 1982, etc.).

* Consider any period of 90 days ending at time *t* with *year(t) ≥ 2011*. To find its class,

* we look at all times $t' = t - 365k$, $k \in \mathbb{N}$, for which we have data.

* We also add $t' = t \pm j - 365k$ with $j = 1$ or 2.

* At all these $t'$ dates, we considered the mean temperature on the 90-days period ending at $t'$. This set of values defines an empirical distribution.

* We then compute the terciles of that empirical distribution built from all those values, which we denote by $\tau_{1,t}$ and $\tau_{2,t}$. We then classify the period ending at *t* using the two terciles.

<figure><img src="image_37.png" /></figure>
An example

Station: 83964, Period: 2017-07-01 - 2017-10-01

<figure><img src="image_38.png" /></figure>

<figure><img src="image_39.png" /></figure>
Another example

<figure><img src="image_40.png" /></figure>

<figure><img src="image_41.png" /></figure>
## About the results

* Basically, we reached the state-of-the-art with standard ML techniques, using existing tools and without any need of new methods and powerful computation platforms.

* There are a few other ML-based works addressing the same problem, using more powerful hardware and developing specific new techniques, but obtaining in general similar performances (in some cases, our tools performed better, never worse).

* For the evaluations themselves, we needed some *feature engineering* effort, that we qualify as moderate.

* However, an important engineering effort was needed on the preparation of the data. We had to study and then to transform, in some cases, several of the parameters that we used in the experiments.

* See our papers for the precise definition of the different relevant quantities and how they are measured.

<figure><img src="image_42.png" /></figure>
Outline

I/ – EvoESN

II/ — ClimateDL: some (other) results

III/ – Some conclusions and perspectives

<figure><img src="image_43.png" /></figure>
* Our main conclusion is that, given that modern ML tools can compete, in this specific type of problem, with today's super-heavy production models, a more extensive and intensive experimental analysis must be performed, with the goal of looking for better accuracy and more insight.

* Having said this, the absolute accuracy is low. Something intuitively clear, now, is the fact that the data we had isn't enough to make good predictions, because of the complexity (including the chaotic side) of the targeted physical system.

* Another idea to continue our project is to look at the tools developed in works such as the two recent papers published on the short-term forecasting (≤ 10 days), at the ML techniques used there, and at the methodology followed to prepare the data. All this must be analyzed to see if something similar can be added to our initial methods.

* A last perspective to this work is to look at more detailed possible outputs than the 3-valued one output.

<figure><img src="image_44.png" /></figure>
Thank you.

gerardo.rubino@inria.fr

<figure><img src="image_45.png" /></figure>