Two deep learning approaches to forecasting disaggregated freight flows: convolutional and encoder–decoder recurrent

Two deep learning approaches to forecasting disaggregated freight flows: convolutional and encoder–decoder recurrent Isidro Lloret1 · José A. Troyano2 · Fernando Enríquez2 · Juan-José González-de-la-Rosa3 Abstract Time series forecasting of disaggregated freight flow is a key issue in decision-making by port authorities. For this purpose and to test new deep learning techniques we have selected seven time series of imported goods from Morocco to Spain through the port of Algeciras, and we have tested two forecasting deep neural networks models: dilated causal convolutional and encoder–decoder recurrent. We have experimented with four different granularities for each series: quarterly, monthly, weekly and daily. The results show that our neural network models can manage these raw series without first removing seasonality or trend. We also highlight the ability of neural models to work with a fixed input size of one year, being able to make good predictions using the same input size for all granularities. The two deep learning models have globally improved the benchmarks of the M4 Competition of forecasting. Each neural network model obtains its best results under different circumstances: the recurrent one with daily granularity and intermittent series, and the convolutional one with weekly and monthly granularities. Keywords Disaggregated freight transport · Time-series forecasting · Machine learning · Deep learning · Dilated convolutional neural network · Encoder–decoder recurrent neural network 1 Introduction Seaports are major infrastructure elements of transportation networks in supply chains. Port authorities are responsible BIsidro Lloret [email protected] José A. Troyano [email protected] Fernando Enríquez [email protected] Juan-José González-de-la-Rosa [email protected] 1Department of Computer Engineering, Engineering School of Algeciras, University of Cádiz, S/N Ramón Puyol Av., 11202 Algeciras, Spain 2Department of Languages and Computer Systems, Computer Engineering School, University of Seville, S/N Reina Mercedes Av., 41012 Sevilla, Spain 3Research Unit in Computational Instrumentation and Industrial Electronics PAIDI-TIC-168, Area of Electronics, Engineering School of Algeciras, University of Cádiz, S/N Ramón Puyol Av., 11202 Algeciras, Spain for making decisions regarding these infrastructures and their facilities in light of freight transport demand forecasts (Systematics 1997; Chopra 2019). Disaggregated freight flow allows for a segmented analysis of future freight transport demand in order to help prioritise transport investments, transport policy development and the growth of the industry of logistics (Havenga 2013). The aim of this work is to evaluate the usefulness of machine learning techniques for the prediction of the flow of goods through seaports. Machine learning techniques base their predictions on inferred models from data. There are multiple machine learning applications in many domains, as healthcare, education and sports, but has only recently been applied to commercial transactions, so its application to the study of international trade trends is limited (Batarseh et al. 2019). We are especially interested in deep learning models. Deep learning is a specific kind of machine learning that allows a machine to learn complicated concepts by building them from simpler concepts. The different components of this type of model are organized into a deep and layered structure (Goodfellow et al. 2016). Our experiments seek to give insights into the applicability of the new deep learning techniques to the forecasting of the flow of goods. Specifically, we have used two deep learning models based on the following techniques: –Convolutional neural networks A type of artificial neural network in which neurons are organized into hierarchical blocks, each of which specializes in enhancing certain features of the input data. –Recurrent neural networks Another type of neural network composed of nodes that, in addition to receiving information from the previous layer, can also obtain it from themselves, thus allowing the sequential processing of data. 2.1 Data description and forecast horizons Our dataset consists of seven disaggregated time series of goods importations from Morocco to Spain through the port of Algeciras in the period 2000–2016, corresponding to the European tariff chapters shown in Table 1. The seven series represent a great variety of seasonality and trend, which will allow us to evaluate the behavior of our forecasting techniques under very different conditions. Maritime traffic time series are updated in the databases with the current timestamp for every port operation (loading or unloading). From the original series, which have a oneminute granularity, we have resampled them into four levels of granularity (quarterly, monthly, weekly and daily), and for each we have considered a multi-step forecast horizon as defined in the M4 Competition (Makridakis et al. 2020). This way we obtain a set of 28 series (available at, Xxxxx 2020)for short, medium and long term forecasting. And therefore we can support very different types of decisions, like strategic, planning and operational decisions. Table 2shows for each granularity level, the forecast horizon, the seasonality period and the time interval of the series. Daily series are the ones presenting two types of seasonality (annual and weekly). For these type of series we have shortened the time interval to the last two and a half years, which is long enough for our prediction purposes. Figure 1shows the seven series with monthly granularity. Validation and test periods are plotted in green and orange, respectively. We describe how we have done the validation and test splitting in Sect. 2.2. As we can see, there are highly seasonal series (like chapter 7), intermittent series (like chapter 10) and series with a mix of trend, seasonality and noise (like chapter 3). Series with small or zero values can cause problems in the calculation of error measurements (Sect. 4.2). To avoid this, we have treated all the series in the same way, adding a constant value of 1 to them. Figure 2shows the autocorrelation plots of the daily series in period 2000–2016. These plots shed light on which series have weekly, annual or both seasonalities. As we can observe, the seven series show different behavioural patterns. The weekly seasonality seems to be present in all the series, even in chapters 7 and 10 where it appears in an overlapping way to the annual seasonality. Chapter 7, shows a well defined annual seasonality. And chapter 10 seems to integrate an annual seasonality with a semiannual one. 2.2 Machine learning approach to forecasting A time series is a sequence of observations sttaken from a variable at successive equally spaced points tin time. For the purpose of multi-step forecasting, we represent a series by The data we have selected for the experimentation consists of seven disaggregated time series of goods importations from Morocco to Spain through the port of Algeciras in the period 2000–2016. The series present a great variety of seasonality and trend, which will allow us to evaluate the behavior of our forecasting techniques under very different conditions. We have also resampled the original series into four granularities (quarterly, monthly, weekly and daily), each associated with a specific forecast horizon (8 quarters, 18 months, 13 weeks, and 14 days), thus covering short, medium and long term forecasting. With this variety we wanted to avoid one of the most common shortcomings of many time series forecasting works, in which conclusions are based on a few or a single time series, as well as on onestep or short-term predictions (Makridakis et al. 2018). We consider this wide range of experiments to be an additional guarantee for the conclusions we draw from our results. Our experiments show that both deep learning models globally beat the overall weighted average (OWA) accuracy measure of all baselines. We use as baselines all the methods included in the M4 Competition benchmarks. By levels of granularity, the recurrent model obtain the best results for daily data and the convolutional model for weekly data. This paper is organized in six sections. The second section deals with defining the specific details of our forecasting task. Section 3 describes the two deep learning models we have used in our experiments. In Sect. 4 we present our experimental methodology. We discuss the results in Sect. 5 and extract conclusions in Sect. 6. 2 The task In this section we describe the time series forecasting task. First, we present the dataset used in our experiments, and then we explain the machine learning approach to forecasting. Table 1European tariff chapter of selected series Chapter Description 3 Fish and crustaceans, molluscs and other aquatic invertebrates 7 Edible vegetables and certain roots and tubers 10 Cereals 16 Preparations of meat, of fish or of crustaceans, molluscs or other aquatic invertebrates 62 Articles of apparel and clothing accessories, not knitted or crocheted 69 Ceramic products 85 Electrical machinery and equipment and parts thereof; sound recorders and reproducers, television image and sound recorders and reproducers, and parts and accessories of such articles Table 2 Time series characteristics by granularity Granularity level Seasonality period Forecast horizon Time interval Quarterly Year 8 2000–2016 Monthly Year 18 2000–2016 Weekly Year 13 2000–2016 Daily Year/Week 14 2014–07 to 2016–12 Fig. 1 Monthly series. Green represents the validation data and orange the test data Fig. 2 Autocorrelation plot of daily series in period 2000–2016 S={s1,s2,...,sT,sT+1,...,sT+h,...},(1) x s T- n+1 s T- n+2 …sT 12…n y sT+1 sT+2 …sT+h 12…h Multi output forecasting model x s T- n+1 s T- n+2 …sT 12…n y sT+1 1 Recurrent forecasting model (a) (b) Fig. 3 Machine learning vector-to-vector models for time series forecasting. The input vector xcontains the nvalues previous to the prediction. The output vector ycontains the predicted values. aThe multiple output model obtains the complete h(forecast horizon) values of the forecast. bThe recurrent model obtains a one-step prediction value that has to be used iteratively to obtain the hvalues so that it is temporarily ordered. We iterate in the same way until we obtain the hone-step predictions for the vector y. Once the mappings of our model are defined, we need a dataset of examples so that the machine learning algorithm can experience through a supervised learning scheme. The where sT is the last observation previous to the prediction of the observed values {sT +1,...,sT +h}, being h the forecast horizon. A machine learning algorithm is able to learn from experience to perform a task by combining an optimization procedure, a cost function, a model and a dataset. Time series forecasting can be performed as a regression task, where a machine learning algorithm is asked to make a prediction y for the next h time steps, conditioned by the n preceding observations x. To address this, we have considered two forecasting strategies (Fig. 3), namely the multiple output and the recursive strategies (Bontempi et al. 2012). In the multiple output strategy, our model defines a mapping from Rn → Rh to get y directly from x, however, in the recursive strategy our model defines a mapping from Rn → R to obtain a one-step prediction y from x. For the next prediction, we modify the vector x by removing the oldest observation and adding y Fig. 4Training sets for multiple output and recursive forecasting strategies X= ⎡ ⎢ ⎢ ⎢ ⎣ s1s2... s n s2s3... s n+1 . . . . . ..... . . sT−h−n+1 sT−h−n+2 ... s T−h ⎤ ⎥ ⎥ ⎥ ⎦ Y= ⎡ ⎢ ⎢ ⎢ ⎣ sn+1 sn+2 ... s n+h sn+2 sn+3 ... s n+1+h . . . . . ..... . . sT−h+1 sT−h+2 ... s T ⎤ ⎥ ⎥ ⎥ ⎦ X= ⎡ ⎢ ⎢ ⎢ ⎣ s1s2... s n s2s3... s n+1 . . . . . ..... . . sT−nsT−n+1 ... s T−1 ⎤ ⎥ ⎥ ⎥ ⎦ Y= ⎡ ⎢ ⎢ ⎢ ⎣ sn+1 sn+2 . . . sT ⎤ ⎥ ⎥ ⎥ ⎦ (a) Multiple output (b) Recursive training sets for our two forecasting strategies are created from the observed values {s1,...,sT}, and consists of matrices Xand Y(Fig. 4), where each training pair is provided by the ith row of Xand Yrespectively. The goal for our model is to learn from the example pairs how to predict y from xwhen using the multiple output strategy, and how to predict the single output yfrom xin the recursive case. The strength of the multiple output strategy is its ability to capture the dependence between the predicted values. In the case of the recursive strategy, every predicted value depends on the accuracy of the previous ones, accumulating errors, so it may deteriorate. The advantage of the recursive strategy is that it results in a training set larger that the direct one, increased by h−1 pairs. In addition, it uses data closer to the prediction as input for the training. In our experiments we have reached better results using the direct strategy for daily and weekly granularities, and using the recursive one for monthly and quarterly granularities. This makes sense because, for monthly and quarterly granularities we have less data points for training. During the learning process, the optimization procedure aims to minimize a cost function, that provides a measure of the error in the model forecasts. This is achieved by adjusting certain internal parameters of the model, represented by the vector θ, that define the mapping y=f(x,θ)of the model; from now, when the forecasting strategy is recursive, we will consider ya one-dimensional vector. Among the possible cost functions, we have chosen the mean absolute error (MAE), defined over the complete dataset of mtraining examples (x,y)by: J(θ)=1 m m  i=1 y−f(x,θ)1,(2) where .1denotes L1norm. Minimizing the MEA cost function leads to a function f(x,θ)that estimates the median value of yfor each x; assuming that the training pairs are a subset of the infinitely many samples from the true data generating distribution (Goodfellow et al. 2016). We have used other regression-specific cost functions and obtained worse results. This is consistent with the recent work of Chai et al. (2019), where they concluded that the MAE cost function can achieve a better robustness and generalization capability than the conventional mean square error (MSE); and the investigation of Qi et al. (2020) about the advantages of the MAE in our context of deep neural network-based vector-to-vector regression. Due to the complexity of neural network models, they are usually optimized using an iterative method, such as gradient descent, that allows controlling the overfitting during optimization. To do this, a percentage of training pairs is separated, so that they do not influence the optimization, although the value of the cost function J(θ)is calculated for these pairs, saving the parameters θeach time a new minimum is reached for them. Once the optimization is complete, the saved values are used to define the model. It is a common practice to split the time series into training and test portions (Hyndman and Athanasopoulos 2018). The test data are used to evaluate the accuracy of a model trained with previous data (Fig. 5a). It is necessary to adjust the model’s hyperparameters before the test forecasting are made. To achieve this, the training data are split into a validation and a new training part (Fig. 5b). Then, the model is trained repeatedly until the hyperparameters are tuned to get the best prediction regarding the validation data. A length for the test and validation data equal to the forecast horizon is just the right size for our task. 3 Neural network architectures The capability of machine learning models in the effective determination of seasonal effect and trend in time series is a controversial issue (Makridakis et al. 2018). In our deep learning models, we have found that it is possible, achieving the best results with an input data size of one year, agreeing with the conclusions of Hamzaçebi (2008), who proposes Fig. 5 Model training, validation and test datasets during: athe test process, bthe validation process dimensional input vector xand a learnable vector k, called kernel or filter. A one-unit convolutional layer outputs a vector ygiven by yt=(xk)t= i xt+i−1ki,(3) where denotes a cross-correlation operator. Each ytcontains a similarity measure, obtained through the dot product of a region of the input x(starting at xt) and the filter k. The output yis sometimes referred to as the feature map. When forecasting time series, the value of ytcannot depend on future values of x, so the causal convolution is used to address it: yt= i xt−i+1ki,(4) where an implicit zero padding of xis assumed in the implementation, so that both xand yare n-dimensional. Figure 6 depics this operation. A causal convolutional layer with m units would provide a n×mfeature map. It is very common to work with a stack of convolutional layers, so they must be able to process a n×mmatrix Xto get a n×mmatrix Y, in which each feature Yt,fis calculated using the 2D-slice Kf,:,:of the filter tensor and a region of the matrix Xending at row t: Yt,f= i,j Xt−i+1,jKf,i,j(5) To deal with long-range time dependencies in the time series, one solution is to increase the input region affecting the features (receptive field). The dilated causal convolution applies the filter over a region larger than its size by skipping input values with a certain rate r: an input size equal to the length of the seasonal period for improving forecasts with seasonal series. 3.1 Convolutional networks A convolutional neural network (Aloysius and Geetha 2017) (CNN) is a type of artificial neural network in which neurons are organized in blocks, resulting in a neuron specialization similar to that present in the neurons of the visual cortex of a biological brain. Although the original concepts of the CNNs are inspired by image recognition, the model is applied very successfully to many other tasks. The hierarchical structure of CNNs makes them less dense, and less prone to overfitting, than the general (fully connected) neural network model. CNNs usually include two types of layers: convolutional and pooling layers. Convolutional layers apply filters by means of operations that involve a region of the input data. In image processing, for example, these filters allow certain aspects of the images to be highlighted, such as borders, vertical lines, etc. Pooling layers are usually applied after a convolution and they summarize, by reducing the dimensionality, the features extracted by the convolutional filters. It is usual to complete the architecture of a CNN with a final fully connected layer that performs the high-level reasoning of the network, using as input the features extracted by the previous convolutional and pooling layers. CNNs designed for image processing work with 2D data, but they can be defined for any other dimension. 1D CNNs, for example, can process sequential data such as natural language texts or time series (Lim et al. 2017). 3.1.1 The dilated causal convolutional layer A 1D-convolutional layer is composed of a number of units that perform the convolution operation. Usually in practice, each unit computes the cross-correlation operation of an n- Input vector 12 345678 12 345678 0 -5 Output vector (features) 13231122 Kernel / Filter vector 5 510 -5 5-10 050 -5 5 From timestep 1 to 8 Dot product Fig. 6 Causal convolutional unit example. The dot product of the kernel vector and a slice of the input is calculated for each time step (the stride value is 1), padding with zeros at the beginning; and the features are stored in the output vector. Causal convolution means that features are not based on future time steps h x1 nx1nxmnxm nxm x (input) 1 2 3 4 5 6 7 n=8 Flatten to (nxm) x1 Dense (h units) 1 2 … h y (output) Feature maps resulting from a stack of 3 convolutional layers with mfilters Fig. 7 DCCNN forecasting model. Without loss of generality, we have assumed an input vector xof length n=8. Each convolutional layer (omitted for clarity) contains mfilters of 2 rows that get a new feature map. The yellow lines show how an increasing dilation rate in powers of 2 enlarge the receptive field affecting the features. Finally, the last feature map is resized so that the dense layer can predict the output y Yt,f= i,j Xt−ir+r,jKf,i,j(6) Figure 7shows the feature maps obtained by a stack of dilated causal convolutional layers, and how a dilation rate that increases exponentially in each layer, results in an exponential growth of the receptive field. To complete the convolutional layer, a learnable bias bfis added per filter, and also an activation function: Yt,f=φ⎛ ⎝ i,j Xt−ir+r,jKf,i,j+bf⎞ ⎠(7) where φis usually the rectified linear activation function (ReLU). 3.1.2 Dilated causal convolutional neural network (DCCNN) In this deep learning approach, we have leveraged new architectures from the work of Oord et al. (2016), based on the use of a stack of dilated causal convolutional layers to deal with long-term dependencies. Our DCCNN model defines a vector-to-vector mapping y=f(x),fromthenprevious values of a time series to the next hforecasted values, and is trained to estimate the median value of yconditioned to x, as explained in Sect. 2.2. The model is implemented using a stack of dilated causal convolutional layers as defined in Eq. (7). The stack has l layers, each consisting of munits (filters), and a dilation rate rstarting at 1 and doubling with each successive layer. The first layer maps the input vector xto a n×mfeature map. Successive layers take as input the feature map of the previous layer to produce a new one, in which each feature is influenced by a receptive field that is consequently duplicated each time, starting from a length of 2 (Fig. 7). We have determined the number of layers lso that the size of the receptive field does not exceed the dimensionality of the input vector x: l=log2n,(8) where .denotes the integer part. Once the input vector xis processed by the learned stack of convolutional layers, the resulting meaningful feature map needs to be reshaped to a vector. Then, a learned dense layer with hunits maps the features vector to the forecasting output vector y. Dilation significantly reduce the number of convolutional layers required, not being necessary the addition of pooling layers. The number of filters mof the convolutional layers is a hyperparameter tuned during model validation. 3.2 Recurrent networks Recurrent neural networks (Mandic and Chambers 2001)is another major architecture in the current catalogue of deep neural network models. They are basically characterized by having layers that, in addition to receiving information from the previous layer of the network, also receive it from the state of the layer itself in the previous iteration, i.e. they have cycles and therefore state memory. This makes them dependent on time and causes the order in which the data are provided to directly influence the results obtained, something appropriate when working with time series. 3.2.1 From the dense layer to the GRU layer A typical dense (fully connected) layer contains mneural units that map from a d-dimensional input vector xto an m-dimensional output vector y, and is defined by y=φ(Wx +b), (9) where φis the activation function that makes nonlinearity possible, Wis a learnable m×dweight matrix, and ba learnable m-dimensional bias vector. If we feedback the output to the input, we get a recurrent neural network (RNN) layer, that can operate on a sequence of inputs x1,x2,...,xτ into a sequence of “hidden” outputs h1,h2,...,hτ, with each output htdepending on the input xt(at the same time step) and the output ht−1(at the previous time step): ht=φ(Wxt+Uht−1+b), (10) where φis usually the tanh activation function and Uis a learnable m×mweight matrix. Considering the RNN layer as a dynamical system, htwill represent the state at time step t. The RNN layer can use the state to store information about activations triggered by previously processed data. Learning about long sequences, using a gradient descent method, presents the main problem that the flow of the back-propagation error tends to disappear. Hochreiter and Schmidhuber (1997) provided a remedy to this problem with the Long Short-Term Memory (LSTM) architecture, in which they achieve a constant error flow through the use of gating mechanisms to control the flow of information. Later, Cho et al. (2014) proposed the Gated Recurrent Unit (GRU) architecture, based on this same concept and structurally similar, but simpler and usually faster than an LSTM, yielding similar results (Chung et al. 2014). In a GRU layer, the output state htis a mixture of the previous state ht−1and the proposed estate ˜ ht, governed in each dimension by an update gate ut, ht=ut˜ ht+(1−ut)ht−1,(11) where denotes the Hadamard product, allowing that a value close to 1 in a dimension of the gate utcauses the update of this dimension in the output state htto the value for this dimension in the proposed state ˜ ht. Similarly, a value close to 0 causes an update in the output state from the previous state, allowing to remember useful information for the future. The update gate utis a mapping from the current input xtand the previous output estate ht−1to the range 0 to 1, ut=σ(Wuxt+Uuht−1+bu), (12) where σis the logistic sigmoid function, Wuand Uuare learnable weight matrices, and buis a learnable bias vector. The output candidate ˜ htis calculated similarly to the output state in an RNN layer, but using a reset gate rtto selectively erase part of the information contained in the previous state ht−1: ˜ ht=φ(Wxt+U(rtht−1)+b). (13) The reset gate rtis computed in a similar way as the update gate: rt=σ(Wrxt+Urht−1+br). (14) As a whole, the GRU layer is capable of preserving relevant information in future, forgetting what is no longer of interest. 3.2.2 Encoder–decoder recurrent neural network (EDRNN) This deep learning approach is based on the encoder–decoder architecture proposed by Cho et al. (2014) and Sutskever et al. (2014), to map sequences to sequences. Our EDRNN model defines a vector-to-vector mapping y=f(x),fromthenprevious values of a time series to the next hforecasted values, which is trained to estimate the median value of yconditioned to x, as explained in Sect. 2.2. The model is implemented using GRU layers as defined in Eqs. (11–14) for our particular case of scalar sequence of inputs xt, and consists of two separate neural networks (encoder and decoder) jointly trained. The encoder–decoder architecture (Fig. 8) allows processing input and output sequences of different lengths. Its design and operation are described below: (a) Encoder A learned GRU layer with munits processes the input vector x=(x1,...,xn)in temporal order, updating the m-dimensional output state htfor every input value so the final state ccodifies the whole history of the input. The number of units mis a hyperparameter tuned during model validation. 1 2 … m GRU (m units) n…21 GRU (m units) c(final state) Dense (1 unit) Decoder output 1 1 State Decoder input (a) Encoder (b) Decoder 1 2 … m 12…h x(input vector) y(output vector) Fig. 8 EDRNN forecasting model. aEncoder. The input vector xis summarized in a state vector c.bDecoder. The state of the GRU layer is initialized to the final state cof the encoder, and the decoder input to xn. In orange we highlight the prediction loop where, for each iteration, the new predicted value is stored in the next component of the output vector y, and the GRU and dense layer outputs are fed back Table 3Benchmarks of the M4 Competition used as baseline Statistical benchmarks without removing the seasonality or trend 1 Naïve 1 All forecasts are equal to the last observation prior to the prediction 2 Seasonal Naïve (Naïve S) All forecasts for the same time step within the seasonal period are equal to the last known observation of the same time step. This means that the last observed seasonal period is repeated consecutively in the prediction Statistical benchmarks removing the seasonality 3 Naïve 2 Like Naïve 1 4 SES Simple exponential smoothing 5 Holt Holt’s exponential smoothing 6 Damped Damped exponential smoothing 7 Combination (Comb) The arithmetic average of methods SES, Holt and Damped 8 Theta It was the winner method of the M3 Competition. It decomposes a series into two or more Theta lines that are extrapolated separately and combined to form a forecast (Assimakopoulos and Nikolopoulos 2000) Machine learning benchmarks removing the seasonality and linear trend 9 MLP A multilayer perceptron 10 RNN A basic recurrent neural network (b) Decoder It consists of a GRU layer with munits and a dense layer with 1 unit. At each time step tthe dense layer maps the state htof the GRU layer to the decoder output yt, and the GRU layer maps the previous state ht−1and the previous decoder output yt−1to the current state ht. Thus, the decoder model can be defined by yt=g(ht−1,yt−1). (15) The decoder is trained according to this mapping, with the initial ht−1set to the final state cof the encoder, and the initial yt−1set to the last encoder input xn.To generate the first prediction y1of the decoder, we start by initializing it in the same way as before. Then, we iterate by feeding back each decoder output as a new decoder input, and each output state of the GRU layer as the new input state, until we complete hconsecutive predictions. 4 Methodology In this section, we present our experimental methodology. First, we briefly describe the baseline methods. In the second subsection, we present the different evaluation metrics we have used. Finally, in the third subsection we specify the main aspects of our experimental scenarios. 4.1 Baseline We have taken as baseline the ten benchmarks of the M4 Competition (Makridakis et al. 2020) listed in Table 3and whose code is available online.1 1https://github.com/M4Competition/M4-methods. In some of these methods, the seasonality of the series is removed when passing a 90% autocorrelation test, and in others the linear trend may also be removed. In these cases, once the forecasts are made and before the results are evaluated, both are restored. 4.2 Evaluation metrics Among the forecast accuracy metrics typically used, we have adopted those of the M4 Competition. These metrics are expressed in terms of: the training data {s1,s2,...,sT}; the test data {sT+1,sT+2,...,sT+h}, where his the forecast horizon; and the forecast data {sT+1,sT+2,...,sT+h}. Detailed information on these metrics is given below: – The symmetric mean absolute percentage error (sMAPE) metric was originally defined by Armstrong (1985) and later modified: sMAPE =2 h T+h  t=T+1 |st−st| |st|+|st|·100(%)(16) This percentage metric provides a result between 0 and 200%. It has the disadvantage of being undefined for intermittent series when any observed value stand its prediction stare equal to zero. Also problematic are small values of st(Hyndman and Koehler 2006). Moreira- Matias et al. (2013) encountered these problems and added a constant to the denominator to produce more accurate statistics, based on the work of Jaynes (2003). We have adopted a similar solution, adding 1 to all the series before making the forecasts, and subtracting it after Hochreiter S, Schmidhuber J (1997) Long short-term memory. Neural Comput 9(8):1735–1780 Hyndman RJ, Athanasopoulos G (2018) Forecasting: principles and practice. OTexts Hyndman RJ, Koehler AB (2006) Another look at measures of forecast accuracy. Int J Forecast 22(4):679–688 Jaynes, ET (2003) In: Bretthorst GL (ed) Probability theory: the logic of science. Cambridge University Press, Cambridge. https://doi. org/10.1017/CBO9780511790423 Lim H, Park J, Han Y (2017) Rare sound event detection using 1d convolutional recurrent neural networks. In: Proceedings of the detection and classification of acoustic scenes and events 2017 workshop, pp 80–84 Lloret I (2020) Imported goods from morocco to Spain through the port of algeciras (2000–2016). Mendeley Data, V1. https://doi.org/10. 17632/gnh24w3mrg.1 Makridakis S, Spiliotis E, Assimakopoulos V (2018) Statistical and Machine Learning forecasting methods: concerns and ways forward. PLoS ONE 13(3):e0194889. https://doi.org/10.1371/ journal.pone.0194889 Makridakis S, Spiliotis E, Assimakopoulos V (2020) The m4 competition: 100,000 time series and 61 forecasting methods. Int J Forecast 36(1):54–74 Mandic DP, Chambers J (2001) Recurrent neural networks for prediction: learning algorithms, architectures and stability. Wiley, Hoboken Moreira-Matias L, Gama J, Ferreira M, Mendes-Moreira J, Damas L (2013) Predicting taxi-passenger demand using streaming data. IEEE Trans Intell Transp Syst 14(3):1393–1402 Oord Avd, Dieleman S, Zen H, Simonyan K, Vinyals O, Graves A, Kalchbrenner N, Senior A, Kavukcuoglu K (2016) Wavenet: a generative model for raw audio. arXiv preprint arXiv:1609.03499 Qi J, Du J, Siniscalchi SM, Ma X, Lee CH (2020) On mean absolute error for deep neural network based vector-to-vector regression. IEEE Signal Process Lett 27:1485–1489 Sutskever I, Vinyals O, Le QV (2014) Sequence to sequence learning with neural networks. In: Advances in neural information processing systems, pp 3104–3112 Systematics C (1997) Inc. NCHRP report 388: a guidebook for forecasting freight transportation demand. Transportation Research Board, Washington DC