Forecasting Uranium Resource Price Prediction by Extreme Learning Machine with Empirical Mode Decomposition and Phase Space Reconstruction

. A hybrid forecasting approach combining empirical mode decomposition (EMD), phase space reconstruction (PSR), and extreme learning machine (ELM) for international uranium resource prices is proposed. In the first stage, the original uranium resource price series are first decomposed into a finite number of independent intrinsic mode functions (IMFs), with different frequencies. In the second stage, the IMFs are composed into three subseries based on the fine-to-coarse reconstruction rule. In the third stage, based on phase space reconstruction, different ELM models are used to model and forecast the three subseries, respectively, according to the intrinsic characteristic time scales. Finally, in the foruth stage, these forecasting results are combined to output the ultimate forecasting result. Experimental results from real uranium resource price data demonstrate that the proposed hybrid forecasting method outperforms RBF neural network (RBFNN) and single ELM in terms of RMSE, MAE, and DS.


Introduction
Uranium resource products have been widely used in economics, military, social lives, and so on and have a revolutionary effect on the real world in many areas.Uranium resource is both the material basis for the development of nuclear energy and a kind of strategic resource.Therefore, more accurate forecasts for international uranium resource prices play an increasingly important role in the development planning and utilization of nuclear energy.
Forecasting the uranium resource prices is one of the most important and challenging tasks due to its inherent nonlinearity and nonstationary characteristics.In the past decades, price prediction has attracted increasing attention by lots of academic researchers.The forecasting approaches used in the literature can be classified into two categories: statistical models and artificial intelligence models [1,2].However, the statistical models cannot effectively capture nonlinear patterns hidden in price time series owing to the fact that these models are developed based on the underlying assumption that the time series being forecasted are linear and stationary [3].In order to overcome this limitation of statistical models, a good deal of nonlinear models have been proposed, among which the artificial neural network (ANN) has attracted a growing interest by researchers due to their excellent nonlinear modeling capability [4][5][6][7][8].Many studies conclude that the ANN model outperforms various conventional statistical models.However, ANN suffers from local minimum traps and the difficulty of determining the hidden layer size and learning rate [9].A new learning algorithm for the single-hidden-layer feedforward neural network (SLFN) called the extreme learning machine (ELM) has been proposed recently and overcomes the aforementioned disadvantages [10,11].In the learning process of ELM, the input weights and hidden biases are randomly chosen, and the output weights are analytically determined by using the Moore-Penrose generalized inverse.ELM can learn much faster with a higher generalization performance than the traditional gradient-based learning algorithms and solves the problem of stopping criteria, learning rate, learning epochs, and local minima [10][11][12][13].In recent years, ELM has attracted a lot of attention and has become an important method in nonlinear modeling [11][12][13].
When building intelligent prediction models directly using original values, it is difficult to obtain satisfactory forecast results due to the high-frequency, nonstationary, and chaotic properties of uranium resource price data.Hence, in order to further improve the prediction performance, before constructing a forecasting model, recent research efforts on modeling time series with complex nonlinearity, dynamic variation, and high irregularity are to initially utilize information extraction techniques to extract features hidden in the data, then use these extracted characteristics to construct the forecasting model [14][15][16][17][18].That is to say, by some means of suitable feature extractions or signal processing methods, the useful or interesting information which may not be observed directly from the original data can be revealed in the extracted features.Thereby, an effective forecasting model possessing better prediction precision will be developed.
Empirical mode decomposition (EMD), based on Hilbert-Huang Transform (HHT), is very suitable for decomposing nonlinear and nonstationary time series, which adaptively represents the local characteristic of the given signal [19,20].By using EMD, any complicated signal can be decomposed into a finite and often small number of Intrinsic Mode Functions (IMFs), which have simpler frequency components and stronger correlations, thus are easier and more accurate to forecast [8].Recently, the EMD has been widely used in many fields, such as in the analysis of the atmosphere time series [21], river water turbidity forecasting [22], crude oil price prediction [23], short-term wind power prediction, and so forth, [8,13,14,16,24].
In this study, we propose a hybrid uranium resource prices forecasting model by integrating EMD, phase space reconstruction (PSR), and extreme learning machine (ELM).Firstly, the original uranium resource price series are first decomposed into a finite number of independent intrinsic mode functions (IMFs), with different frequencies.Secondly, the IMFs are composed into three subseries based on the fine-to-coarse reconstruction rule.And then, based on phase space reconstruction, different ELM models are used to model and forecast the three subseries, respectively, according to the intrinsic characteristic time scales.Finally, these forecasting results are combined with the ultimate forecasting result output.Moreover, experimental results from real uranium resource price data demonstrate that the proposed hybrid forecasting method outperforms RBF neural network (RBFNN) and single ELM in terms of RMSE, MAE, and DS.
The rest of this paper is organized as follows.Section 2 gives brief overviews of EMD, PSR, and ELM.The proposed model is described in Section 3. Section 4 compares the experimental results obtained by the proposed hybrid approach and those by RBF neural network and single ELM, and this paper is concluded in Section 5.

Empirical Mode Decomposition (EMD).
EMD is a new signal processing technique.Unlike wavelet decomposition, EMD is not required to determine a filter base function before decomposition.The main idea of EMD is to decompose original time series data into a sum of oscillatory functions, namely, intrinsic mode functions (IMFs).In the EMD, the IMFs must satisfy two conditions: (1) the number of extrema (sum of maxima and minima) and the number of zero crossings must either equal or differ at most by one, and (2) at any point, the mean value of the envelope defined by the local maxima and the envelope defined by the local minima is zero.
The essence of EMD is the sifting process which extracts IMFs from the original data.The algorithm of EMD is described as follows [19,20,25].
Step 1. Identify all the local extrema including the minimum values and maximum values in time series data ().
Step 2. Generate the upper and lower envelopes  max () and  min () by a cubic spline line.
Step 3. Calculate the mean value  1 () from the upper and lower envelopes and then generate the mean envelope as Step 4. Calculate the difference between the time series data () and the mean value  1 ().The first difference ℎ 1 () is designed as a protointrinsic mode function Step 5. Check whether the protointrinsic mode function ℎ 1 () satisfies the properties of IMF or not.Ideally, ℎ 1 () should be an IMF.However, it may generate a new extremum and shift or exaggerate the existing extrema in the sifting process.
If properties of ℎ 1 () satisfy all the requirements of an IMF, ℎ 1 () is denoted as the th IMF   () and substitutes the residue  1 () for the original time series data (); that is, Otherwise, ℎ 1 () is not an IMF.Then, it substitutes ℎ 1 () for the original time series ().
Step 6. Repeat from Step 1 to Step 5.The sifting process stops when the residue satisfies one of the termination criteria.First, the residue or the th component is smaller than the predetermined threshold or becomes a monotonic function such that no more IMF can be extracted.Second, the number of zero crossings and extrema is the same as that of the successive sifting step.
By using the above algorithm, the original time series data () can be decomposed into  modes and a residue as follows: where  is the number of IMFs,   () represents IMFs which are nearly orthogonal to each other and periodic, and   () is the final residue which is a constant or a trend.By the sifting process, each IMF is independent and specific for expressing the local characteristics of the original time series data.The set of IMFs is derived from high frequency to low frequency, while   () represents the central tendency of data series ().
In addition, EMD can also be taken as a filter of high pass, band pass, or low pass.

Phase Space Reconstruction (PSR).
Takens embedding theorem [26] provides theoretical foundation for the analysis of time series generated by nonlinear dynamical systems.Later, Sauer et al. [27] show that a phase space can be reconstructed from a univariate chaotic time series.Let a univariate time series {  }  =1 , where  is the length of the time series, generated from a -dimension chaotic attractor, then a phase space   of the attractor can be reconstructed by using the delay coordinate defined as where  is known as the embedding dimension of reconstructed phase space and  is the delay constant.
The selection of the embedding dimension  and the delay constant  is very important for prediction modeling [3].Therefore, for a given delay time , a time series {  }  =1 is represented in the so-called "phase space" by a set of delay vectors (DVs) x() = [ − , . . .,  − ] of a given embedding dimension .
(1) Determination of the Delay Time  and the Minimum Embedding Dimension .The entropy ratio (ER) method [28] is a novel method for determining the set of parameters for a phase space representation of a time series.Based upon the differential entropy, both the optimal embedding dimension and time lag are simultaneously determined.
Based upon the probability density function of data, the differential entropy is defined as Particularly convenient is the Kozachenko-Leonenko (K-L) estimate of the differential entropy where  is the number of samples in the dataset,   is the Euclidean distance of the th delay vector to its nearest neighbour, and   is the Euler constant.For a given embedding dimension  and time lag , let (, , ) denote the differential entropies estimated for time delay embedded versions of a time series .
The K-L estimates for the time delay embedded versions of the original time series (, , ) and ( , , , ) are computed using (7) for increasing  and  (index  refers to the th surrogate).To determine the optimal embedding parameters, the ratio  (, ) =  (, , ) ⟨ ( , , , )⟩  (8) needs to be minimized, where ⟨⋅⟩ denotes the average over .
To penalize for higher embedding dimensions, the minimum description length (MDL) method is superimposed, yielding the "entropy ratio" (ER) where  is the number of delay vectors, which is kept constant for all values of  and  under consideration.The minimum of the plot of the entropy ratio yields the optimal set of embedding parameters.
(2) Identification of Nonlinearity.The delay vector variance (DVV) method [29] is a novel analysis of a time series which examines a signal's unpredictability by observing the variability of the targets belonging to sets of similar delay vectors (DVs).The DVV method can be summarized as follows for a given embedding dimension.
(i) The mean   and standard deviation   are computed over all pairwise distances between DVs.
(ii) The sets Ω  are generated, which consist of all DVs that lie closer to () than a certain distance.The distances are taken from the interval [  −     ,   +     ], for example, uniformly spaced, where   is a parameter controlling the span over which to perform the DVV analysis.
(iii) For every set Ω  , the variance of the corresponding targets,  2  , is computed.The average over all sets, divided by the variance of the time series  2  , yields the measure of unpredictability  * 2 : In the following study, the linear or nonlinear nature of the time series is examined by performing DVV analyses on both the original and a number of surrogate time series, using the optimal embedding dimension of the original time series.
(3) Identification of Chaotic Characteristic.When the largest Lyapunov exponent of the system is larger than zero, it indicates that there is a chaotic attractor which can be used to measure the chaotic degree [30].The largest Lyapunov exponent is computed by the Wolf method [31].

Extreme Learning Machine (ELM).
Suppose there are  distinct samples {(  ,   )}  =1 , where   ∈   and   ∈   .The SLFN with  hidden neurons can be described as where Thus, there also exist parameters   ,   , and   such that The above equations can be compactly described as where Unlike the traditional function approximation theories which require the adjustment of input weights and hidden layer biases, the input weights and hidden biases are randomly generated.
The smallest norm least-squares solution of the above linear system is where  + is the Moore-Penrose generalized inverse of matrix .Owing to the Moore-Penrose generalized inverse, the learning speed is dramatically increased for the single hidden-layer feedforward neural network [12].

Uranium Resource Price Forecasting Method Based on EMD-PSR-ELM
The proposed hybrid approach for international uranium resource price forecasting, namely, EMD-PSR-ELM, combines EMD, PSR, and ELM, and it is composed of four main stages.These four stages are described as follows.
Stage 2 (combination of the decomposition components).In this stage, owing to each IMF with different time scales, high frequency and low frequency components can be obtained by combining IMFs according to the frequency from high to low, and the residue is treated separately.The process can be performed as follows.
Step 1.The mean of each IMF is evaluated in order.
Step 2. Determining the first IMF   () whose mean significant deviation from zero.
Stage 3 (ELM modeling).This stage can be subdivided into three steps as follows.
Step 1 (phase space reconstruction).Firstly, some parameters such as the delay time and the embedding dimension should be determined.
Step 2 (data normalization).The data needs to be represented in a normalized scale for ELM training and prediction.Thus, in this study, the dataset of three phase space domain are linearly scaled in the range of [0, 1] using the following expression: where   represents the th dimension of th sample after normalization;   is th dimension of th sample before normalization;  min values of th dimension, respectively.Then the data was divided into two sets, training set and testing data.
Step 3 (ELM prediction).Regression forecast model is set up for the high frequency component, the low frequency component, and the residue by using ELM, respectively.

Stage 4 (result composition). The final prediction results
were obtained by compositing the prediction values after denormalization.
The proposed EMD-PSR-ELM method is schematically depicted in Figure 1.

Data Description and Evaluation Criteria.
In this paper, for evaluating the performance of proposed EMD-PSR-ELM prediction model, the real time series about international uranium resource prices data were chosen as experimental samples.The data used in this study are monthly data which are freely available from the IndexMundi website (http://www.indexmundi.com)and cover the period from October 1982 to September 2012 with a total of 360 values.Firstly, the time series are analyzed by chaos theory.Delay time  = 4 and embedding dimension  = 3 can be simultaneously determined (Figure 2).Secondly, nonlinearity is analyzed by using DVV method.Due to the standardization of the distance axis, these plots can be conveniently combined in a scatter diagram, where the horizontal axis corresponds to the DVV-plot of the original time series and the vertical to that of the surrogate time series.If the surrogate time series yield DVV-plots similar to the original, the "DVV scatter diagram" coincides with the bisector line, and the original time series is probably linear.The deviation from the bisector line is, thus, a measure of nonlinearity; they can be seen from Figures 3 and 4. Additionally, Wolf et al. method [31] is employed to compute the largest Lyapunov exponent  = 0.0339.Therefore, we can conclude that the international uranium resource prices times series has chaotic characteristic owing to  > 0. When these parameters are gotten, the phase space can be reconstructed; that is to say, these optimal embedding dimensions and delay times are  used to construct the input matrix.It is easy to know that there are in total 351 data points in the phase space.The data was divided into two sets, training set and testing data, and the first 326 data points are used as the training samples while the remaining 25 data points are used as the testing samples.The prediction performance is evaluated using measures including mean absolute error (MAE), root mean squared error (RMSE), and directional prediction statistics (DS) [32].These measures are as follows: where   and x are the actual and prediction values, respectively,  is the sample size, and Obviously, the smaller RMSE, MAE, and larger DS mean better performance.

Forecasting Results.
According to previous steps shown in Section 3, we carried out the prediction experiments.First, using the EMD technique, the international uranium resource price series can be decomposed into seven independent IMFs and one residue.Figure 5 shows the decomposition results for uranium resource price series using EMD.
Through the analysis of the mean value of each IMF after the decomposition of international uranium resource prices time series data by using EMD, it can be seen from Figure 6 that the IMF whose mean significant deviation from zero is IMF3.Therefore, the partial reconstruction with IMF1 and IMF2 represents high frequency components of international uranium resource prices time series data, with the characteristics of small amplitudes, which contains the effects of markets' short-term fluctuations; and the partial reconstruction with IMF3, IMF4, IMF5, IMF6, and IMF7 represents low frequency components, which should be representative of the effect of these events [33].The residue is treated as the long-term trend during the evolution of uranium resource prices separately.Figure 7 shows the three components of uranium resource price monthly data series from Oct. 1982 to Sept. 2012.
After reconstructing phase space, high frequency, low frequency, and residue are individually used for building ELM prediction models, and the predicted values for the next 25 months uranium resource prices are then obtained by different models.In building EMD-PSR-ELM, the number of hidden nodes is set to 9 for ELM.The error results of the ELM model for the three components are shown in Table 1.As can be seen in Table 1, ELM prediction models perform well for all three components.
The forecasting results of the proposed EMD-PSR-ELM model are compared with the PSR-ELM model and PSR-RBFNN model, which uses non-EMD forecasting variables.The PSR-RBFNN forecasting model with three input nodes and one output node is built by using the same phase space reconstruction method as above.The neural network toolbox of MATLAB software is adapted in this study.Mean squared error is 0.1, the spread of radial basis functions is 1.6, and the default settings of neural network toolbox are used for   the remaining parameters.In PSR-ELM, the number of hidden nodes is set the same as the above.
Figure 8 depicts the actual prices and the predicted values from the EMD-PSR-ELM, PSR-RBFNN, and PSR-ELM models.From this figure, it can be observed that there is a smaller deviation between the actual and predicted values using the proposed EMD-PSR-ELM model.
Table 2 compares the prediction results obtained with the EMD-PSR-ELM, PSR-RBFNN, and PSR-ELM models for the next 25 months uranium resource price.It can be observed from Table 2 that the proposed EMD-PSR-ELM model provides a better forecasting result than the PSR-RBFNN and PSR-ELM models in terms of RMSE, MAE, and DS.

Conclusions
This study has presented a forecasting model for uranium resource price by integrating EMD, PSR, and ELM.In terms of the experimental results presented in this study, we can draw the following conclusions.
(1) EMD can fully capture the local fluctuations of data and can be used as a preprocessor to decompose the complicated raw data into a finite set of IMFs and a residue, which have simpler frequency components and high correlations.
(2) On the one hand, although empirical mode composition is an important tool for multiscale modeling, it suffers from mode mixing and mode splitting.To rectify this issue, the Ensemble EMD [34] can be used in the future.In addition, It is also worth trying to employ Multivariate EMD and in particular the Noise Assisted MEMD in order to calculate EMD free of any artifacts [35][36][37].On the other hand, other multiscale models, such as the synchrosqueezed transform in the forecasting task, could also be considered [38].
(3) The network topology of the model has an important influence on prediction performance for ELM and RBF neural networks.It is more objective to identify the chaotic characteristic of uranium resource price series and determine the embedding dimension of the reconstructed phase space by quantitative calculation, and then the determined embedding dimension can be served as the numbers of nodes in input layer for the single hidden-layer feedforward neural network and RBF networks.
(4) RMSE, MAE, and DS are used to measure the forecasting accuracy of the forecasting models.The experimental results reveal that the proposed hybrid EMD-PSR-ELM approach outperforms the other models such as PSR-RBFNN and PSR-ELM.Therefore, the proposed method is very suitable for prediction with nonlinear, nonstationary, and highly complex data and is an efficient method for uranium resource price prediction.

Figure 2 :
Figure 2: Plots of the ER for uranium resource price monthly data series.

Figure 5 :
Figure 5: The IMFs and residue for uranium resource price monthly data series.

8 Figure 6 :Figure 7 :Figure 8 :
Figure 6: The mean of each IMF of uranium resource price monthly data series.
is the actual output vector, and   is the threshold of the th hidden neuron.(⋅)represents the activation function of hidden neuron, and the operation   ⋅   denotes the inner product of   and   .If the SLFN can approximate the  samples with a zero error then we have ]  is the weight vector connecting the th hidden neuron and the input neurons,   = [ 1 ,  2 , . . .,   ]  is the weight vector connecting the th hidden neuron and the output neurons,   = [ 1 ,  2 , . . .,   ]

Table 1 :
The error statistics of the ELM model for the three components.

Table 2 :
Performance of the three forecasting methods for international uranium resource prices.