FPGA-Based Implementation of a Multilayer Perceptron Suitable for Chaotic Time Series Prediction

Abstract

Many biological systems and natural phenomena exhibit chaotic behaviors that are saved in time series data. This article uses time series that are generated by chaotic oscillators with different values of the maximum Lyapunov exponent (MLE) to predict their future behavior. Three prediction techniques are compared, namely: artificial neural networks (ANNs), the adaptive neuro-fuzzy inference system (ANFIS) and least-squares support vector machines (SVM). The experimental results show that ANNs provide the lowest root mean squared error. That way, we introduce a multilayer perceptron that is implemented using a field-programmable gate array (FPGA) to predict experimental chaotic time series.

1. Introduction

Chaotic time series prediction techniques have the challenge of forecasting long horizons of observed data, which must guarantee low error that can be measured as root mean squared error (RMSE). Among the kinds of chaotic signals that require attention to predict their future behavior are for example those associated to human diseases and pathologies, namely: epilepsy, pancreatic beta cell, electroencephalogram, electrocardiogram, and so on. Besides, predicting chaotic time series is a challenge because they are unpredictable data. Nowadays, many researchers are investigating the appropriate prediction technique for this kind of time series. In the recent work [1], the authors showed the prediction capabilities of three techniques, namely: artificial neural networks (ANNs), the adaptive neuro-fuzzy inference system (ANFIS) and least-squares support vector machines (LS-SVM). From the experiments, SVM and ANFIS were overcome by ANNs. Those experimental results are used herein as the base to highlight that ANNs are a good option for chaotic time series prediction. There may exist other time series prediction techniques, all of them focused on predicting future events using past information. However, only a few techniques can be successfully applied in the long-term prediction of chaotic time series, due to problems such as: accumulation of errors, lack of information and sensitivity to initial conditions.

In the state of the art, one can find several techniques in order to solve the problem on chaotic time series prediction considering short-term and long-term prediction horizons, e.g., ANNs [2,3,4], fuzzy systems [5,6,7], SVM [8,9], hybrid systems [10,11,12], etc. ANNs have been employed in [2] to predict the spontaneous fluctuation in vascular caliber; in that work, the authors basically used a multilayer perceptron (MLP). In [13], the authors showed that the dynamics of a chaotic system are present in the time series analysis of the river flow, and by combining phase-space reconstruction and ANN techniques, they predicted the daily river flow. In [14], it was shown that financial time series have chaotic behavior, and they can be predicted using an MLP. On the other hand, SVM has been widely used in regression and classification problems, and its principles can easily be extended to time series prediction, as shown in [8]. The authors in [15] predicted the global incident solar radiation applying SVM, and also in [16], it was applied to model river water pollution. The other prediction technique applied herein is ANFIS, which has emerged as a synergic hybrid intelligent system because it combines the human-like reasoning style of a fuzzy logic system with the learning and computational capabilities of ANNs [17]. In [18] the authors reported the prediction of chaotic time series of Mackey and Glass, the Dow Jones Company and the Mexican exchange stock.

This article considers as cases of study different chaotic time series generated by a chaotic oscillator and with different values of their maximum Lyapunov exponent (MLE). In chaos theory, it is well known that the higher the MLE value, the higher the unpredictability of a chaotic time series. We show preliminary results already provided in [1] by applying three techniques: ANNs, LS-SVM and ANFIS. As ANNs showed the lowest RMSE for predicting the same chaotic time series, then we introduce an MLP implemented in a field-programmable gate array (FPGA). This electronic system can find application in healthcare informatics [19].

Section 2 describes the generation of chaotic time series from a chaotic oscillator implemented in an FPGA, and their analysis to evaluate the MLE is performed by using the free and online available time series analysis (TISEAN) tool. Section 3 describes three time series prediction techniques: ANNs, SVM and ANFIS. Their comparison is given in Section 3.4. Section 4 introduces an MLP implemented in an FPGA and shows experimental prediction results. Finally, the conclusions are listed in Section 5.

2. Generating Chaotic Time Series

The time series data that are used herein to predict their future behavior are generated by the chaotic oscillator described by:

$$\begin{array}{c}\dot{x}=y\\ \dot{y}=z\\ \dot{z}=-ax-by-cz+d_{1}f\left(x\right)\end{array}$$

(1)

where $x,y,z$ are the state variables; $a,b,c,d_1$ are real and positive coefficients; and $f\left(x\right)$ is a nonlinear function described by (2) with $s\left(x\right)=60.606$ as the slope, $\alpha =0.0165$ the break points and $k=1$ the saturation level.

$$\begin{array}{c}\hfill f\left(x\right)=\left\{\begin{array}{cc}\hfill k& x<\alpha \hfill \\ \hfill s\left(x\right)& -\alpha \le x\le \alpha \hfill \\ \hfill -k& x<-\alpha \hfill \end{array}\right.\end{array}$$

(2)

This chaotic oscillator was implemented using an FPGA, as already shown in [20]; their coefficients $a,b,c,d_1$ were varied as shown in [21,22,23] in order to maximize MLE, which is related to the unpredictability of the chaotic time series data. For instance,

Table 1. Optimized maximum Lyapunov exponent (MLE) for a two-scroll attractor and their associated ($a,b,c,d_1$) coefficient values.

Case	a	b	c	${\mathit{d}}_1$	MLE
1	1.0000	1.0000	0.4997	1.0000	0.3761
2	1.0000	0.7884	0.6435	0.6665	0.3713
3	0.8661	1.0000	0.3934	0.9903	0.3607
4	0.7746	0.6588	0.5846	0.4931	0.3460
5	1.0000	0.7000	0.6780	0.1069	0.3437
6	1.0000	0.7000	0.7000	0.2542	0.3425
7	0.5610	0.9470	0.3460	0.6810	0.2225
8	0.7000	0.7000	0.7000	0.7000	0.1117

lists the coefficient values associated with eight different MLE values. Figure 1 shows examples of chaotic time series corresponding to the state variable x and with different values of MLE. In this work, the time series data of this state variable x are used to predict future behavior of six and 12 steps ahead. The amplitudes of all the chaotic time series data are normalized in the range $[\pm 1]$.

3. Chaotic Time Series Prediction Techniques

Time series prediction techniques can be as simple as predicting one future value, known as single step ahead prediction, or complex so as to predict more than one future value, which can be described by $\widehat{y}(t+k)$ to save a specific input vector at any time t. That way, we are interested in predicting horizons of six $k=6$ and twelve $k=12$ steps ahead. In this manner, the format of training data is described by (3) and (4), respectively. In both cases, the prediction of six and 12 steps ahead is verified by evaluating the root mean squared error (RMSE) described by (5), which has been used to measure the prediction performance of ANNs, ANFIS and SVM techniques in [1].

The chaotic time series data are divided into two sets: training and validation sets. In those cases, we used 2000 samples for training the prediction technique (ANN, ANFIS and SVM) and 2000 samples for validation. In the training phase, the parameters of the prediction techniques are adjusted. In the validation phase, the prediction errors are analyzed when new data are used.

$$\begin{array}{ccc}\hfill \widehat{y}(t+6)& =& x_1\left(t\right),x(t-6),x(t-12),x(t-18)\hfill \end{array}$$

(3)

$$\begin{array}{ccc}\hfill \widehat{y}(t+12)& =& x_1\left(t\right),x(t-12),x(t-24),x(t-36)\hfill \end{array}$$

(4)

$$RMSE=\sqrt{\frac{1}{n}\sum _{i=1}^n{[y_i-{\widehat{y}}_i]}^2}$$

(5)

3.1. Artificial Neural Networks

An ANN is a set of elementary processing units called neurons whose processing capability is stored in the connections by synaptic weights and whose adaptation depends on learning. Figure 2 shows the structure of an artificial neuron where $x_{1,2,\dots ,j}$ represents the input, which has an associated weight w. The artificial neuron has a bias ($\gamma$), and f denotes the activation function. That way, the state of an artificial neuron is evaluated by (6).

$$f\left(u\right)=\sum _{i=1}^n{x}_j{w}_j+b$$

(6)

A particular case of an ANN is the multilayer perceptron (MLP), which consists of an input and an output layer, and one or more hidden layers. In an MLP, the input signals propagate through the network layer-by-layer [24]. In an MLP, the numbers of neurons and layers can be the difference between success and failure [25]. For instance, one guideline to choose the number of hidden neurons is the application of the geometric pyramid rule given in [25], where the number of neurons follows a pyramid shape. First, one must set the number of neurons in the input layer n and the number of neurons in the output layer m.

Table 2. Computation of the number of hidden layer neurons using the geometric pyramid rule.

Hidden Layers	Geometric Pyramid Rule
1	$r={(n·m)}^{1/2}$	$h_1=r$
2	$r={(n/m)}^{1/3}$	$h_1=m·r^2,h_2=m·r$
3	$r={(n/m)}^{1/4}$	$h_1=m·r^3,h_2=m·r^2,h_3=m·r$
⋮	⋮	⋮
i	$r={(n/m)}^{1/i+1}$	$h_1=m·r^i,h_2=m·r^{i-1},...,h_i=m·r$

lists several equations associated with different numbers of hidden layers ($h_i$). As one can see, for a three-layer network with n input neurons and m output neurons, the hidden layer has $\sqrt{n·m}$ neurons. A similar rule applies to networks with i hidden layers; however, the evaluation of the number of hidden layer neurons is slightly more complex.

In this work, we designed two MLPs from the application of the geometric pyramid rule. The first one consists of five layers, with sixteen neurons in the input layer ($n=16$), three hidden layers ($h_1,h_2,h_3$) and one neuron in the output layer $m=1$. Therefore, from Table 2: $r=2$ and $h_1=8,h_2=4,h_3=2$. The second MLP consists of six layers with 32 neurons in the input layer ($n=32$), four hidden layers and one neuron in the output layer $m=1$, so that $r=2$ and $h_1=16,h_2=8,h_3=4,h_4=2$.

A third MLP introduced in [26] is also used herein. It consists of six layers, with different numbers of neurons.

Table 3. Characteristics of three MLPs: the one from [26] and the others designed by applying the geometric pyramid rule with five and six layers.

Characteristics	MLP from [26]	MLP Applying the Geometric Pyramid Rule
Number of layers	6	5	6
Number of neurons	$4\times 7\times 4\times 8\times 5\times 1$	$16\times 8\times 4\times 2\times 1$	$32\times 16\times 8\times 4\times 2\times 1$
Activation function	$h\times h\times h\times h\times h\times h\times l$	$h\times h\times h\times h\times l$	$h\times h\times h\times h\times h\times l$
Learning algorithm	Levenberg–Marquardt algorithm

lists the characteristics of the three MLPs that are used herein to predict chaotic time series, where the neurons of the input and hidden layers have a hyperbolic tangent function (h) and the neuron of the output layer has a linear function (l).

3.2. Adaptive Neuro-Fuzzy Inference System

ANFIS is a type of system that incorporates ANNs and fuzzy logic. In this case, the neurons of the ANN incorporate TSK (Takagi–Sugeno–Kang) so that the system is adaptive, fault tolerant, distribute, learns from experience and presents generalizability. Assuming that the fuzzy control system has an output f and two inputs (x; y), thus the two fuzzy rules are represented by (7).

$$\begin{array}{c}\hfill R^1:\mathrm{IF}x\mathrm{is}{A}_1\mathrm{AND}y\mathrm{is}{B}_1,\mathrm{THEN}{f}_1=p_{1}x+q_{1}y+r_1\\ \hfill R^2:\mathrm{IF}x\mathrm{is}{A}_2\mathrm{AND}y\mathrm{is}{B}_2,\mathrm{THEN}{f}_2=p_{2}x+q_{2}y+r_2\end{array}$$

(7)

Figure 3 shows the general structure of ANFIS, where the function of nodes in each layer is:

• Layer 1 contains input nodes passing external signals to the next layer.
• Layer 2 performs the parameter adjustment of the input membership function.
• Layers 3 and 4 perform fuzzy operations.
• Layers 5 and 6 weight and provide the output of the system.

In ANFIS, four inputs and 2${}^4$ = 16 rules are used. Furthermore, genfis1is used to generate an initial membership function matrix (generalized bell function) from the training data.

3.3. Support Vector Machine

The prediction of time series using least-squares support vector machines (LS-SVM) is done using a set of historical values of time series as inputs and a single output as the target value. When the prediction of the time series is developed by SVMs, given time series data $y_t$ composed by a set of consecutive real values $y_t=\{y_t\in \mathbb{R}\mid y_1,y_2,\dots ,y_n\}$, the prediction is represented by dividing $y_t=\{y_1,y_2,\dots ,y_n\}$ in windows $w=(y_t,\dots ,y_{t+p-1})$ of size p, in order to find a function $f(y_t,\dots ,y_{t+p-1})=y_{t+p}$. The radial basis function is used as a kernel function because it provides excellent performance under assumptions of general softness [27]. In the space of characteristics, LS-SVM has the form $f\left(\mathit{x}\right)={\mathit{w}}^T\phi \left(\mathit{x}\right)+b$, where $\mathit{w}$ is a vector of weights. Figure 4 shows a functional structure detailed in [28].

3.4. Chaotic Time Series Prediction by ANN, ANFIS and SVM

The three MLPs from Table 3 are compared to ANFIS and SVM for the prediction of the eight chaotic time series data listed in Table 1.

Table 4. RMSE obtained during the validation stage with a six steps ahead prediction horizon.

MLE	ANN			ANFIS	LS-SVM
	5 Layers	6 Layers	[26]
0.1117	7.31 × 10${}^{-4}$	3.91 × 10${}^{-4}$	1.70 × 10${}^{-3}$	0.0153	0.0120
0.2225	3.00 × 10${}^{-3}$	1.80 × 10${}^{-3}$	5.49 × 10${}^{-3}$	0.0460	0.0180
0.3425	1.52 × 10${}^{-3}$	2.98 × 10${}^{-4}$	1.83 × 10${}^{-3}$	0.0276	0.0012
0.3437	8.78 × 10${}^{-4}$	3.59 × 10${}^{-4}$	1.80 × 10${}^{-3}$	0.0244	0.0041
0.3460	1.09 × 10${}^{-3}$	5.14 × 10${}^{-4}$	2.09 × 10${}^{-3}$	0.0180	0.0103
0.3607	3.22 × 10${}^{-3}$	2.29 × 10${}^{-3}$	5.03 × 10${}^{-3}$	0.0446	0.0816
0.3713	1.03 × 10${}^{-3}$	9.63 × 10${}^{-4}$	2.17 × 10${}^{-3}$	0.0315	0.0037
0.3761	1.93 × 10${}^{-3}$	8.14 × 10${}^{-4}$	2.87 × 10${}^{-3}$	0.0373	0.0062

lists the prediction error of the eight chaotic time series associated with different MLE values. The RMSE was obtained by using a prediction horizon of six steps ahead. Figure 5 shows the performance comparison of the prediction techniques with respect to their RMSE value obtained during the validation stage. As one sees, the three MLPs produce lower error. Comparing the MLPs listed in Table 3, it is observed that the six-layer and five-layer topologies, both designed from the geometric pyramid rule, provide the least prediction error.

Considering the previous results, the MLPs with five layers and six layers are selected to predict the eight chaotic time series, but now to achieve a twelve steps ahead horizon. Figure 6 shows the RMSE obtained during the validation stage, again for the eight MLE values associated with the chaotic time series from Table 1.

4. FPGA Implementation of an MLP

This section shows the FPGA implementation of the MLP with five layers shown in Figure 7, which follows the methodology already given in [20]. In this MLP, the information flows in one direction, so that the prediction depends largely on the input signals. The MLP is interfaced with a personal computer using the serial communication protocol to feed the input data, as described in [20]. The MLP is implemented using the FPGA Ciclone IV GX FPGA DE2i-150 from Altera. The descriptions of the digital blocks are performed using fixed-point notation, which provides a high speed and a low cost of hardware resources. In this manner, $N=32$ bits are used to represent a fixed point format, where the most significant bit represents the sign, three bits are used to represent the integer part and 28 bits are used to represent the fractional part. This allows representing numbers over the range $\{-2^{N-1},2^{N-1}-1\}$, with a numerical resolution equal to 3.7253 × 10${}^{-9}$, enough to represent the weight and bias of the neurons in the MLP.

From Figure 2, one can identify the basic building blocks of a neuron, namely: multiplication, sum and a block dedicated to implementing the hyperbolic tangent function. The VHDL block descriptions of the adder and multiplier are shown in Figure 8; they perform operations with 2’scomplement, and they have a clock (CLK) and a reset (RST) pin to make the system synchronous.

The hyperbolic tangent activation function can be approached by using the piecewise-linear (PWL) approximation technique developed by [29] by using segments of order two, as shown in (8), where $\beta$ and $\theta$ determine the slope and gain $H_{s1}\left(z\right)$ between $-L\le z\le L$. In this manner, the hyperbolic tangent is approached herein by (9), where $\theta =0.25$, $L=2$ and $\beta =1$.

$$\begin{array}{c}\hfill H_{s1}\left(w\right)=\left\{\begin{array}{cc}\hfill w(\beta -\theta w)& \mathrm{for}0\le w\le L\hfill \\ \hfill w(\beta +\theta w)& -L\le w\le 0\hfill \end{array}\right.\end{array}$$

(8)

$$\begin{array}{c}\hfill G_{s1}\left(w\right)=\left\{\begin{array}{cc}\hfill 1& \mathrm{for}L\le w\hfill \\ \hfill H_{s1}& \mathrm{for}-L<w<L\hfill \\ \hfill -1& w\le -L\hfill \end{array}\right.\end{array}$$

(9)

Figure 9 shows the block description to implement the hyperbolic tangent activation function using PWL functions described by (9). The value of $H_{s1}\left(w\right)$ is obtained by simplifying operations; therefore, it can be expressed as $w(\beta \pm \theta w)$, where an extra multiplier allows one to perform the conditions described by (8). That is, when the product $\theta w$ is positive, then $H_{s1}\left(w\right)=w(\beta -\theta w)$, and when the product $\theta w$ is negative, then $H_{s1}\left(w\right)=w(\beta +\theta w)$. On the other hand, the value of $G_{s1}\left(w\right)$ is obtained using a comparator block.

Figure 10 shows the synthesis of an artificial neuron. Recall that the design of each neuron depends on the number of inputs, weights, biases and activation function. That way, the MLP shown in Figure 7 has 30 neurons between the first and fourth layers, which are synthesized by using the neuron from Figure 10, because all 30 neurons have an activation function of a hyperbolic type, as listed in Table 3. The last layer has an activation function of a linear type, so that the neuron is modified at $f\left(u\right)$ from Figure 10.

Figure 11 shows the experimental results of the five-layer MLP, using the validation dataset and six steps ahead prediction with the chaotic time series with MLE = 0.3761. The MLP achieves an RMSE = 0.019.

Figure 12 shows the experimental results of the five-layer MLP, using the validation dataset and twelve steps ahead prediction with the chaotic time series with MLE = 0.3761. The MLP achieves an RMSE = 0.0266.

5. Conclusions

This article showed the performance comparison among three MLPs, LS-SVM and ANFIS to predict chaotic time series that were generated from a chaotic oscillator implemented in an FPGA. Eight chaotic time series were generated with different MLE values, and the prediction techniques were compared with respect to the RMSE. As MLPs showed the least RMSE, we introduced the FPGA implementation of an MLP with five layers designed by applying the geometric pyramid rule. The experimental results were performed for a prediction horizon of six and twelve steps ahead. The MLP predicted the chaotic time series with a very low RMSE. This confirmed the usefulness of the FPGA to implement an MLP to predict chaotic time series data.