Time series forecasting using amplitude-frequency analysis of STL components

The problem of economic time series analysis and forecasting using amplitude-frequency analysis of STL decomposition is considered. An amplitude-phase operator was chosen as an apparatus for extraction the series harmonic components, the advantages of which (compared to the Fourier transform) are: calculations speed, result accuracy, simplicity and interpretability of software implementation. The forecast quality was carried out using the MAPE metric. Significantly higher prediction quality was achieved compared to Facebook Prophet forecasting package.


Introduction
In this work the subject of analysis is economic time series, such as: exchange rates, the cost of precious metals and oil. The series are difficult to predict by classical methods, since they have a pronounced residual component. For example, the ARIMA and ARMA models [1,2], artificial neural network [3][4][5][6] do not take these differences of time series into account and, as a rule, they are able to predict only the trend and the seasonal components, which leads to a low forecast quality.
It is proposed to investigate the features of the STL decomposition [7,8] of the considered time series and use the amplitude-phase operator as an apparatus for amplitude-frequency analysis.

Amplitude and phase operators
A significant feature of the proposed approach is amplitude-phase operator (AFO) for the amplitudefrequency analysis of time series STL components [9][10][11]. AFO has several advantages over classical spectral analysis tools such as Fast Fourier Transform (FFT). Computational speed is especially important for data processing tasks. Compared to FFT, the number of real-valued operations in processing the signal spectrum is tens of times less for one-dimensional data (e.g., a time series) and hundreds of times less for processing data defined on a two-dimensional grid of nodes (e.g., digital images) [12].
AFO in the classical form processes a continuous signal, but there is also a discrete interpretation of the method [13]. If the time series consists of N measurements, then, according to the Kotelnikov theorem [14], harmonics τ μ = a μ cos µt + b μ sin µt with numbers µ = 1,…,[N/2] can be extracted. Further, for the AFO method, it is necessary to determine the order n of the trigonometric polynomial T n (t), which approximates the time series. On the one hand, µ ≤ n ≤ [N/2], on the other hand, as n grows, the So, let s = min{r: rµ − 1 ≥ n}, m = (s + 1)µ, signal defined on uniform grid of m nodes t k . Then values of harmonic τ μ in points t k are calculated by the formula When implementing AFO, the following features were taken into account. First, since the grid of nodes in the data is fixed, and when implementing AFO, it can change (insignificantly) depending on extracted harmonic number, sometimes it was necessary to choose unevenly located points as nodes, but this fact does not lead to a large calculation error if the series has a high sampling frequency ( is about several tens or more). Secondly, in order to save RAM, it was decided to store not a whole harmonic, but its Fourier coefficients, which are calculated using AFO by the formulas: Thirdly, the discrepancy of the nodes grids was taken into account: AFO parameters are calculated for the grid of nodes located at [0,2π), and the time series measurements are either in natural points or in datetime format points (Python).

Features of STL decomposition for analyzed time series
The time series, forecast of which is the purpose of this work, have typical features of STL decomposition (figure 1). Namely, the trend has a structure close to piecewise linear. The amplitude of a seasonal component is insignificant (about 1-2% of time series average value). The residual component, on the contrary, makes a significant contribution to the STL decomposition and has a specific structure (for more details in Section 3.3). Residual component amplitude is about 5-10% of time series average value.
Thus, the analysis of STL time series decomposition leads to the following conclusion. Subject to the availability of expert data on the behavior of the time series trend at the forecast period, it is possible to achieve high forecast accuracy if the residual component will be qualitatively investigated. If such data are not available, then classical methods of trend forecasting should be considered. Next, models for trend and residual component predicting will be proposed. Seasonality in examples will be continued in periodicity for prediction period.

Trend forecasting models
So, in considered time series, trend has a structure close to piecewise linear, and linearity often persists over several months ( figure 1 (b)). This allows in some cases to make a qualitative forecast by linear extrapolation of a period (month) preceding to forecast. Namely, the trend curve is interpolated by a straight line at two boundary points of the previous period and continues for one period into the future.
An alternative method of trend forecasting is extrapolation from known values using an Exponentially Weighted Moving Average (EWMA) [15,16]. Based on the known W values of a trend, using an EWMA with the α parameter, a forecast is built for one point. Then the window W is shifted 3 one unit to the right and a following prediction is obtained in the same way. The optimal EWMA parameters were selected: W = 2 periods = 60 days, α = 0.001.

Amplitude-frequency analysis of residual component
During the analysis of residual component, it was noted that there is a distinct structure: some harmonics have significant amplitude and qualitatively approximate the series. If they dominate in training sample, they are more likely to dominate in the future. Figure 2 (a) shows random component of time series (EUR/RUB rate, January -November 2020) and its approximation by harmonics with numbers [6,4,11,17,5], which have amplitudes [0.49, 0.39, 0.35, 0.35, 0.31] respectively.
The question arises: how to find the required dominant harmonics. The following decision was taken. Let find the percentile Q(γ) of the γ level for sample obtained from residual absolute values. In other words, we find the "typical amplitude" for the series. For example, in the given example, this value should be of the order of 1.9 ( figure 2 (b)). For γ = 0.95, we obtain Q(γ) = 1.86. Then find the numbers of series harmonics that have the largest amplitude and their sum is close enough to Q(γ).

Examples and comparative analysis of the method
During numerical experiments for forecasting time series, the following optimal algorithm macroparameters were selected. Data sampling was one day: this is a maximum sampling frequency of considered time series. One month was chosen as time series period. For the training sample, we take 5-11 periods preceding one predicted month. Since computational experiments have not shown an improvement in prediction quality when selecting macroparameters on a validation set, its use is impractical.
MAPE (Mean Absolute Percentage Error) metric was chosen to assess a forecast quality, because, firstly, it is absolute and allows you to evaluate a result independently data scale; secondly, the predicted series obviously do not contain zero or close to zero values (in the denominator), therefore, it is possible not to modify the metric [17,18].
The prediction of all STL decomposition components is compared with the prediction obtained using the Facebook Prophet library [19,20]. Figure 3 shows examples of forecasting time series: (a) USD/RUB exchange rateforecast for August 2020 based on data for January -July 2020; (b) GBR/RUB exchange rateforecast for September 2020 based on data for January -August 2020.

Conclusion
The time series considered in this work have a specific structure and cannot be predicted by classical methods. This is due to the fact that usually the residual component is excluded from a forecast, and this should not be done in the analyzed time series, because this component has a significant amplitude and makes a considerable contribution to STL decomposition.
The proposed method of amplitude-frequency analysis of a residual component made it possible to achieve a high forecast quality. In addition, the use of AFO method instead of FFT significantly increases the speed of calculations, which is especially important when analyzing a big data.