## 3 Results
### Experimental Analyses  
To assess the performance of the neural networks, we adopted the following calculations for errors and correlation coefficients: mean absolute error (MAE), root mean square error (RMSE), index of agreement (IA), determination coefficient (R\({}^{2}\)) and mean absolute percentage error (MAPE), which are defined as follows:  
\[\text{MAE}=\frac{1}{\text{N}}\ \sum_{\text{i=1}}^{\text{N}}\lvert\text{O}_{ \text{i}}-\text{P}_{\text{i}}\rvert\]  
\[\text{RMSE}=\sqrt{\frac{1}{\text{N}}\ \sum_{\text{i=1}}^{\text{N}}\left(\text{O}_{ \text{i}}-\text{P}_{\text{i}}\right)^{2}}\]  
\[\text{IA}=1-\frac{\sum_{\text{i=1}}^{\text{N}}\left(\text{O}_{\text{i}}-\text {P}_{\text{i}}\right)^{2}}{\sum_{\text{i=1}}^{\text{N}}\left(\lvert\text{O}_{ \text{i}}-\overline{\text{O}}\rvert+\lvert\text{P}_{\text{i}}-\overline{ \text{O}}\rvert\right)^{2}}\]  
\[\text{R}^{2}=(\frac{1}{\text{N}-1}\sum_{\text{i=1}}^{\text{N}}(\frac{\text{O}_ {\text{i}}-\overline{\text{O}}}{\sigma_{\text{O}}})(\frac{\text{P}_{\text{i}} -\overline{\text{P}}}{\sigma_{\text{P}}}))^{2}\]  
\[\text{MAPE}=\text{median}(\lvert\text{O}_{\text{i}}-\text{P}_{\text{i}}\rvert /\text{P}_{\text{i}})\]  
where N is the number of samples. O and P represent the observed and predicted values. Overbar and \(\sigma\) represent the average and the standard deviation of corresponding variables. In all experiments, the PM\({}_{2.5}\) sequence between 19 May 2014 and 31 December 2016 was used for training (90%) and testing (10%). The sequence from 1 January 2017 to 31 December 2017 was used for prediction. Figure 4 shows the wavelet spectrum of the original PM\({}_{2.5}\) sequence, the differences of wavelet spectrum between the original sequence and the ones trained from Exp. 1 by the BP algorithm and Exp. 2 by the EEMD+BP algorithm. Comparing Figure 4b and Figure 4c, we see that, when BP is applied directly to the original sequence (Exp. 1, Figure 4b), the errors are reduced all over the spectrum, with maximum errors between 64 to 128 days. However, if applying an EEMD decomposition before BP fitting (Exp. 2, Figure 4c), the errors can be further reduced, especially in the time period of 8 to 64 days. It indicates that individual BP training can produce an optimized fit for each IMF, especially for medium to low-frequency sequences. However, if the BP training is applied to the original sequence, although it still produces the best fit in principal, the fitting errors span the entire frequency spectrum. Figure 5 shows the comparison of observed and predicted values for high-frequency modes (IMFs 1 and 2) and low-frequency modes (IMFs 3 to 8). For high-frequency modes, the MAE is reduced from 27.6 in Exp. 1 to 24.8 in Exp. 2. In contrast, the MAE is reduced from 21.0 to 4.8 for the low-frequency modes. Note that the R\({}^{2}\) of the Exp. 2 reaches 0.99 (Figure 5d), indicating that the individually trained BP network works perfectly for the low-frequency IMFs. However, we also see that the BP algorithm faces the same challenge in dealing with high-frequency IMFs for both Exp. 1 and 2. This is understandable, as the IMF-1 includes stochastic processes that cannot be trained by a neural network based on a variable-variable relationship. As shown in the EOF maps in Figure 3, the PM\({}_{2.5}\) variability in Beijing can be affected by the surrounding area, implying that spatial information would be helpful for improving the Beijing PM\({}_{2.5}\) prediction. In Exp. 3, a LSTM network was applied directly to the original PM\({}_{2.5}\) sequence. In contrast to the BP algorithm, LSTM is able to train the PM\({}_{2.5}\) sequence without additional meteorological variables, as it tests the autocorrelation of the sequence itself. In Exp. 4, a CNN algorithm was used to gather the spatial information, that was then transferred to the LSTM for prediction. This is to test how spatial information can help the prediction of PM\({}_{2.5}\) in Beijing. Figure 6 compares the MAEs for all experiments in five pollution levels. Comparing Exp.