Forecasting Bitcoin closing price series using linear regression and neural networks models

Time series analysis

Collected data

We tested our algorithms on six daily price series. Three of them are stock market series, all the data were extracted from the ’Historical Data’ available on (http://www.finance.yahoo.com) website; the other ones are cryptocurrencies, namely Bitcoin, Ethereum and Litecoin price daily series, all the data were extracted from (http://www.coinmarketcap.com) website.

• Daily stock market prices for Microsoft Corporation (MSFT), from 9/12/2007 to 11/11/2011.

• Daily stock market prices for Intel Corporation (INTC), from 9/12/2007 to 11/11/2010.

• Daily stock market prices for National Bankshares Inc. (NKSH), from 6/27/2008 to 8/29/2011.

• Daily Bitcoin, Ethereum and Litecoin price series, from 15/11/2015 to 12/03/2020.

We state once more that we choose these price series and the related time intervals as benchmark to compare our results with well known literature results obtained by using other methods.

Specifically, we have chosen for the stock market series the same time intervals chosen in Kazem et al. (2013). The choice of Bitcoin as criptocurrency is quite natural since it represents about 60 % of the Total Market Capitalization. We chose Ethereum and Litecoin since they are among the most important and well-known cryptocurrencies. It is worth noting that, for the stock market series we used the same data of the work we compare to, whereas for the cryptocurrencies we used all the available data to have more significant results.

The dataset was divided into two sets, a training part and a testing part. After some empirical test the partition of the data which lead us to optimal solutions was 80 % of the daily data for the training dataset and the remaining for the testing dataset.

Data pre-processing

For both models we prepared our dataset in order to have a set of inputs (X) and outputs (Y) with temporal dependence. We performed a one-step ahead forecast: our output Y is the value from the next (future) point of time while the inputs X are one or several values from the past, i.e., the so called lagged values. From now on we identify the number of used lagged values with the lag parameter. In the Linear Regression and Univariate LSTM models the dataset includes only the daily closing price series, hence there is only one single lag parameter for the close feature. On the contrary, in the Multiple Linear Regression and Multivariate LSTM models the dataset includes both close and volume (USD) series, hence we use two different lag parameters, one for the close and one for the volume feature. In both cases, we attempted to optimize the predictive performance of the models by varying the lag from 1 to 10.

Univariate versus multivariate forecasting

A univariate forecast consists of predicting time series made by observations belonging to a single feature recorded over time, in our case the closing price of the series considered. A multivariate forecast is a forecast in which the dataset consists of the observations of several features. In our case we used:

• for BTC, ETH and LTC series all the features provided by Coinmarketcap website: Open, High, Low, Close, Volume.

• for MSFT, INTC, NKSH series all the features provided by Yahoofinance website: Date, Open, High, Low, Close, Volume.

We observed that adding features to the dataset did not lead to better predictions, but performance and sometimes also results worsened. For this reason, we decided to use in the multivariate analysis only the close and volume features, that provided the best results.

Statistical analysis

As a first step we carried out a statistical analysis in order to check for non-stationarity in the time series. We used the augmented Dickey-Fuller test and autocorrelation plots (Banerjee et al., 1993; Box & Jenkins, 1976). A stochastic process with a unit root is non-stationary, namely shows statistical properties that change over time, including mean, variance and covariance, and can cause problems in predictability of time series models. A common process with unit root is the random walk. Often price time series show some characteristics which makes them indistinguishable from a random walk. The presence of such a process can be tested using a unit root test.

The ADF test is a statistical test that can be used to test for a unit root in a univariate process, such as time series samples. The null hypothesis H₀ of the ADF test is that there is a unit root, with the alternative H_a that there is no unit root. The most significant results provided by this test are the observed test statistic, the Mackinnon’s approximate p-value and the critical values at the 1%, 5% and 10% levels.

The test statistic is simply the value provided by the ADF test for a given time series. Once this value is computed it can be compared to the relevant critical value for the Dickey-Fuller Test.

Critical values, usually referred to as α levels, are an error rate defined in the hypothesis test. They give the probability to reject the null hypothesis H₀. So if the observed test statistic is less than the critical value (keep in mind that ADF statistic values are always negative (Banerjee et al., 1993), then the null hypothesis H₀ is rejected and no unit root is present.

The p-value is instead the probability to get a ”more extreme” test statistic than the one observed, based on the assumed statistical hypothesis H₀, and its mathematical definition is shown in Eq. (3). (3)${\displaystyle p_{value}=P\left(\right.t\ge t_{observed}\left|\right.H_0\left)\right.}$

The p-value is sometimes called significance, actually meaning the closeness of the p-value to zero: the lower the p-value, the higher the significance.

In our analysis we performed this test using the adfuller() function provided by the statsmodels Python library, and we chose a significance level of 5%.

Furthermore, the autocorrelation plot, also known as correlogram, allowed us to calculate the correlation between each observation and the observations at previous time steps, called lag values. In our case we employed the autocorrelation_plot() function provided by the python Pandas library (Mckinney, 2011).

Forecasting

We decided to follow two different approaches: the first uses two well-known statistical methods: Linear Regression (LR) and Multiple Linear Regression (MLR). The second uses two very common neural networks (NN): Multilayer Perceptron (MLP) NN and Long Short-Term Memory (LSTM) NN. The reasons of this choices are explained below.

LSTM networks

Long Short-Term Memory networks are nothing more than a prominent variations of Recurrent Neural Network (RNN). RNN’s are a class of artificial neural network with a specific architecture oriented at recognizing patterns in sequences of data of various kinds: texts, genomes, handwriting, the spoken word, or numerical time series data emanating from sensors, markets or other sources (Hochreiter & Schmidhuber, 1997). Simple recurrent neural networks are proven to perform well only for short-term memory and are unable to capture long-term dependencies in a sequence. On the contrary, LSTM networks are a special kind of RNN, able at learning long-term dependencies. The model is organized in cells which include several operations. LSTM hold an internal state variable, which is passed from one cell to another and modified by Operation Gates (forget gate, input gate, output gate). These gates control how much of the internal state is passed to the output and work in a similar way to other gates. These three gates have independent weights and biases, hence the network will learn how much of the past output and of the current input to retain and how much of the internal state to send out to the output.

In our case the inputs are the closing prices of the previous days and the number of values considered depends on the lag parameter. The output is the forecast price. We used the Keras framework for deep learning. Our model consists of one stacked LSTM layer with 64 units each and the densely connected output layer with one neuron. We used Adam optimizer and MSE (mean squared error) as a loss.

We optimized our LSTM model searching for the best set of epochs and batch size ”hyperparameters” values. These hyperparameters strongly depend on the number of observations available for the experiment. Due to the recently birth of the cryptocurrency markets, the dimensions of our datasets are quite limited (around 1000 observations), therefore we decided to vary the epochs hyperparameter from 300 to 800 with a step of 100. The Keras LSTM algorithm we used sets as default value for batch size 32. So, for each fixed epoch, we trained the model varying the batch size within the interval $\left[\right.22,82\left]\right.$ with a step of 10. We did not take into account values less than 300 epochs, nor greater than 800 in order to avoid underfitting and overfitting problems. Furthermore, we did not consider batch size values less than 22, since they would lead to extremely long training times. Similarly, batch size values greater than 82 would not allow to find a good local minimum point of the chosen loss function during the learning procedure. The results obtained during the hyperparameters tuning are shown in Fig. 2.

This figure shows the MAPE error as a function of the batch size hyperparameter, for each fixed epoch. As can be seen from the figure, we considered the batch size equal to 72 to be the optimal value. In fact, it is an excellent compromise, having a low MAPE value, which is also practically the same for all tested epochs. The optimal choice for the epochs hyperparameter is 600, which is the one that minimizes the MAPE error for batch size equal to 72, and is consistently among the best choices for almost all batch sizes considered. Therefore, the best set of epochs and batch size ”hyperparameters” values we chose is 600 and 72, respectively.

Time regimes

The time series considered are found to be indistinguishable from a random walk. This peculiarity is common for time series of financial markets, and in our case is confirmed by the predictions of the models, in which the best result is obtained considering only the price of the previous day.

The purpose is to find an approach that allow us to avoid time series differencing technique, in view of the fact that we are interested in prices and not in price variations represented by integrated series of d-order. For this reason, each time series was segmented into short partially overlapping sequences, in order to find if shorter time regimes are present, where the series do not resemble a random walk. Finally, to continue with the forecasting procedure, a train and a test set were identified within each time regime.

For each regime we always sampled 200 observations - namely 200 daily prices. The beginning of the next regime is obtained with a shift of 120 points from the previous one. Thus, every regime is 200 points wide and has 80 points in common with the following one.

We chose a regime length of 200 days because, in this way, we obtain at least 5 regimes (from 5 to 12) for each time series to test the effectiveness of our algorithms, without excessively reducing the number of samples needed for training and testing. The choice was determined also according to the following: we performed the augmented Dickey-Fuller test on subsets of the data, starting from the whole set and progressively reducing the data window and sliding it through the data. The first subset of data that does not behave as random walks appears at time interval of 230 days, which we rounded to 200.

Since the time series considered have different lengths, the partition in regimes has generated:

• Bitcoin, Ethereum and Litecoin: 12 regimes

• Microsoft: 8 regimes

• Intel and National Bankshares: 5 regimes

From a mathematical point of view, the used approach can be described as follows.

Let us target a vector $\overrightarrow{OA}$ along the t axis, with length 200. This vector is identified by the points O(1, 0), A(a, 0) ≡ (200, 0). The length of this vector represents the width of each time regime.

Let $\overrightarrow{OH}$ be a fixed translation vector along the t axis, identified by the points O(1, 0) and H(h, 0) ≡ (120, 0). The length of $\overrightarrow{OH}$ represents the translation size.

For the sake of simplicity, let us label the $\overrightarrow{OA}$ and $\overrightarrow{OH}$ vectors with $\overrightarrow{A}$ and $\overrightarrow{H}$.

Let $\overrightarrow{A^{\prime }}$ be the vector $\overrightarrow{A}$ shifted by $\overrightarrow{H}$ and ${\overrightarrow{A}}^n$ the vector $\overrightarrow{A}$ shifted by n times $\overrightarrow{H}$.

Therefore, the vector that identifies the nth sequence to be sampled along the series is given by: (6)${\displaystyle {\overrightarrow{A}}^n=\overrightarrow{A}+n\overrightarrow{H}}$ where $n\in \left[\right.0,\frac{D-A}{h}\left]\right.$, being D the dimension of the sampling space, A the time regimes width and h the translation size.

So the nth time regime is given by: (7)${\displaystyle R^n=f\left(\right.{\overrightarrow{A}}^n\left)\right.=f\left(\right.\overrightarrow{A}+n\overrightarrow{H}\left)\right.}$ where f is the function that maps the values along the t axis (dates) to the respective regimes y values (actual prices).

Performance Measures

To evaluate the effectiveness of different approaches, we used the relative Root Mean Square Error (rRMSE) and the Mean Absolute Percentage Error (MAPE), defined respectively as: (8)${\displaystyle relativeRMSE=\sqrt{\frac{1}{N}\sum _{i=1}^N\left(\right.\frac{y_i-f_i}{y_i}{\left)\right.}^2}}$ (9)${\displaystyle MAPE=\frac{1}{N}\sum _{i=1}^N\left|\right.\frac{y_i-f_i}{y_i}\left|\right.}$

In both formulas y_i and f_i represent the actual and forecast values, and N is the number of forecasting periods. These are scale free performance measures, so that they are well appropriate to compare model performance results across series with different orders of magnitude, as in our study.

Time series analysis

In Fig. 3 we report the decomposition of Bitcoin (Fig. 3A-3D) and Microsoft (Fig 3E-3H) time series, for comparison purposes, as obtained using the seasonal_decompose() method, provided by the Python statsmodels library (Skipper & Perktold., 2010).

The seasonal_decompose() method requires to specify whether the model is additive or multiplicative. In the Bitcoin time series, the trend of increase at the beginning is almost absent (from around 2016-04 to 2017-02); in later years, the frequency and the amplitude of the cycle appears to change over time. The Microsoft time series shows a non-linear seasonality over the whole period, with frequency and amplitude of the cycles changing over time. These considerations suggest that the model is multiplicative. Furthermore, if we look at the residuals, they look quite random, in agreement with their definitions. The Bitcoin residuals are likewise meaningful, showing periods of high variability in the later years of the series.

It is also possible to group the data at seasonal intervals, observing how the values are distributed and how they evolve over time. In our work we grouped the data of the same month over the years we considered. This is achieved with the ’Box plot’ of month-wide distribution, shown in Fig. 4 (Fig. 4A: Bitcoin; Fig. 4B: Microsoft). The Box plot is a standardized way of displaying the distribution of data based on five numbers summary: minimum, first quartile, median, third quartile and maximum. The box of the plot is a rectangle which encloses the middle half of the sample, with an end at each quartile. The length of the box is thus the inter-quartile range of the sample. The other dimension of the box has no meaning. A line is drawn across the box at the sample median. Whiskers sprout from the two ends of the box defining the outliers range. The box length gives an indication of the sample variability, and for the Bitcoin samples shows a large variance, in almost all months, except for April, September and October. Not surprisingly, bitcoin volatility is much higher than Microsoft one. The line crossing the box shows where the sample is centred, i.e., the median. The position of the box in its whiskers and the position of the line in the box also tell us whether the sample is symmetric or skewed, either to the right or to the left. The plot shows that the Bitcoin monthly samples are therefore skewed to the right. The top whiskers is much longer than the bottom whiskers and the median is gravitating towards the bottom of the box. This is due to the very high prices that Bitcoin reached throughout the period between 2017 and 2018. These large values tend to skew the sample statistics. In Microsoft, an alternation between samples skewed to the left and samples skewed to the right occurs, except for the sample of October that shows a symmetric distribution. Lack of symmetry entails one tail being longer than the other, distinguishing between heavy-tailed or light-tailed populations. In the Bitcoin case we can state that the majority of the samples are left skewed populations with short tails. Microsoft shows an alternation between heavy-tailed and light-tailed distributions. We can see that some Microsoft samples, particularly those with long tails, present outliers, representing anomalous values. This is due to the fact that heavy tailed distributions tend to have many outliers with very high values. The heavier the tail, the larger the probability that you will get one or more disproportionate values in a sample.

Tables 1 and 2 show the statistics calculated for each time series and for each short time regime. The unit of measurement of the values in the tables is the US dollar ($). In Table 1 we can observe that the only series for which the trimmed mean, obtained with trim_mean() method provided by the Python scipy library (Jones, Oliphant & Peterson, 2001), with a cut-off percentage of 10%, is significantly different from the mean are Bitcoin, Ethereum and Litecoin. In particular the trimmed mean decreased. This is due to the fact that these cryptocurrencies, for a long period of time, registered a large price increment and this implies a shift of the mean to the right (i.e., to highest prices). This confirms that cryptocurrencies distribution is right-skewed. Table 2 shows that stock market series time regimes present a lower σ than BTC, ETH and LTC ones, namely that cryptocurrencies distribution has higher variance.

Figures 5 and 6 show the autocorrelation plots of BTC and MSFT series. The others stock market series are not presented because they show the same features of the MSFT series. Both autocorrelation plots (Figs. 5C and 6C) show a strong autocorrelation between the current price and the closest previous observations and a linear fall-off from there to the first few hundred lag values. We then tried to make the series stationary by taking the first difference. The autocorrelation plots of the ’differences series’ (Figs. 5C and 5D) show no significant relationship between the lagged observations. All correlations are small, close to zero and below the 95% and 99% confidence levels.

As regards the augmented Dickey-Fuller results, shown in Table 3, looking at the observed test statistics, we can state that all the series follows a unit root process. We remind that the null hypothesis H₀ of the ADF test is that there is a unit root. In particular, all the observed test statistics are greater than those associated to all significance levels. This implies that we can not reject the null hypothesis H₀, but does not imply that the null hypothesis is true.

Observing the p-values, we notice that for the stock market series we have a low probability to get a “more extreme” test statistic than the one observed under the null hypothesis H₀. Precisely, for both MSFT and INTC we got a probability of 29%, for NKSH a probability of 25%. The same considerations also apply to the Bitcoin, Ethereum and Litecoin cryptocurrency time series. We conclude that H₀ can not be rejected and so each time series present a unit root process.

We conclude that all the considered series show the statistical characteristics typical of a random walk.

Time series forecasting

Tables 4 and 5 show the best results, in terms of MAPE and rRMSE, obtained with the different algorithms applied to the entire series. From now on, let us label the closing and the volume features lag parameters with k_p and k_v respectively. In particular, Table 4 reports the results obtained using the Linear Regression algorithm for univariate series forecast, using only closing prices, and the Multiple Linear Regression model for multivariate series, using both price and volume data.

Table 5 shows the results obtained with the LSTM neural network, distinguishing between univariate LSTM, using only closing prices, and multivariate LSTM, using both price and volume data.

Small values of the MAPE and rRMSE evaluation metrics suggest accurate predictions and good performance of the considered model.

From the analysis of the series in their totality, it appears that linear models outperforms neural networks. However, for both models, the majority of best results are obtained for a lag of 1,thus confirming our hypothesis that the series are indistinguishable from a random walk.

In order to perform the time series forecasting, we also implemented a Multi-Layer Perceptron model. Since the LSTM network outperforms the MLP one, we decided to show only the LSTM results. This is probably due to the particular architecture of the LSTM network, that is able to capture long-term dependencies in a sequence.

It should be noted that better predictions are obtained for stock market series rather than for the cryptocurrencies one. In particular, the best result is obtained for Microsoft series, with a MAPE of 0, 011 and k_p equal to 1. This is probably due to the high price fluctuations that Bitcoin and the other cryptocurrencies have suffered during the investigated time interval. This is confirmed by the statistics shown in Table 1. It must be noted that the addition of the volume feature to the dataset does not improve the predictions.

In order to perform prices forecast we changed the approach and decided to split the time series analysis using shorter time windows of 200 points, shifting the windows by 120 points, with the aim of finding local time regimes where the series do not follow the global random walk pattern.

Tables 6 and 7 show the results obtained with our approach of partitioning the series into shorter sequences. Let us label the moving step forward with h. Particularly, in Table  6 are presented the results obtained using the Linear Regression algorithm for univariate series forecast, using only closing prices, and the Multiple Linear Regression model for multivariate series, using both price and volume data. This approach, has the advantage of being simple to implement and requires low computational complexity. Nevertheless, has led to good results, similar to those present in the literature, if not better as in the Microsoft, Bitcoin and National Bankshares cases, where the MAPE error is lower that 1%.

Table 7 shows the results obtained with the LSTM neural network, distinguishing between univariate LSTM, using only closing prices, and multivariate LSTM, using both price and volume data. For each time regimes we show the best results obtained on a specific time window defined by the k_p and k_v values reported in Tables 6 and 7. Note that we highlighted the best results in bold. In particular, it is worth noting that introducing the time regimes, the best result is obtained for the Bitcoin time series, outperforming also the financial ones.

These results show how such innovative partitioning approach allowed us to avoid the ”random walk problem”, finding that best results are obtained using more than one previous price. Furthermore, this method leads to a significant improvement in predictions. It is worth noting that, from this analysis the best result arise from the Bitcoin series, with a MAPE error of 0.007, a temporal window k_p of 7 and a translation step h of 120, obtained using both regression models and LSTM network.

Another interesting consideration that arises from the results is that, as stated previously in the analysis of the series in their entirety, the linear regression models generally outperform the neural networks ones, while in the short-time regimes approach the different models yielded to similar results.

For a direct feedback we report in Table 8 the best results obtained in the papers we compared to and our best ones. In the event that the best MAPE error results from different models, we consider the model whose computational complexity is the least as best. It is noticeable that our results outperform those obtained in the benchmark papers, providing notable contribution to the literature.