ARTICLE IN PRESS
Atmospheric Environment 40 (2006) 1774–1780 www.elsevier.com/locate/atmosenv
Air quality forecasting using a hybrid autoregressive and nonlinear model Asha B. Chelani, S. Devotta National Environmental Engineering Research Institute, Nagpur 440 020, India Received 10 October 2005; received in revised form 27 October 2005; accepted 9 November 2005
Abstract The usual practices of air quality time-series forecasting are based on applying the models that deal with either the linear or nonlinear patterns. As the linear or nonlinear behavior of the time series is not known in advance, one applies the number of models and finally selects the one, which provides the most accurate results. The air pollutant concentration time series contain patterns that are not purely linear or nonlinear and applying either technique may give inadequate results. This study aims to develop a hybrid methodology that can deal with both the linear and nonlinear structure of the time series. The hybrid model is developed using the combination of autoregressive integrated moving average model, which deals with linear patterns and nonlinear dynamical model. To demonstrate the utility of the proposed technique, nitrogen dioxide concentration observed at a site in Delhi during 1999 to 2003 was utilized. The individual linear and nonlinear models were also applied in order to examine the performance of the hybrid model. The performance is compared for one-step and multi-step ahead forecasts using the error statistics such as mean absolute percentage error and relative error. It is observed that hybrid model outperforms the individual linear and nonlinear models. The exploitation of unique features of linear and nonlinear models makes it a powerful technique to predict the air pollutant concentrations. r 2005 Elsevier Ltd. All rights reserved. Keywords: Time-series forecasting; ARIMA; Nonlinear dynamics; Hybrid model
1. Introduction In the air quality literature, time-series analysis is generally carried out to understand the cause and effect relationships, which in turn helps in forecasting the future concentrations. In this direction, a class of techniques including autoregressive integrated moving average (ARIMA) or Box–Jenkins models (Shi and Harrison, 1997; Milionis and Davies, 1994; Zennetti, 1990) and structural models Corresponding author.
E-mail address:
[email protected] (A.B. Chelani). 1352-2310/$ - see front matter r 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.atmosenv.2005.11.019
(Schlink et al., 1997) have been applied to analyze air pollutant concentrations. These approaches are widely applied in the air-quality literature due to the lack of data on emissions of air pollutants. Although these models are quite flexible as they can represent several different types of time series, their major limitation is the pre-assumed linear form of the model. The approximation of linear models to real-world problems is not always satisfactory. For example, the air pollutant concentrations are influenced by several factors in the atmosphere and prediction using linear models may not always give reasonable results (Benarie, 1987). As an
ARTICLE IN PRESS A.B. Chelani, S. Devotta / Atmospheric Environment 40 (2006) 1774–1780
alternative, nonlinear models have been proposed in the literature. Artificial neural networks are one of the potential examples of nonlinear models that are applied to model and predict air pollutant concentrations (Gardner and Dorling, 1998). These models are generally developed using the external inputs such as meteorology and emissions and output is the air pollutant concentration at a site (Gardner and Dorling, 1999; Chelani et al., 2002). The application of these models is, however, restricted to the some particular cases where the data on emissions and meteorological parameters are available. Recently the time-series forecasting based on the nonlinear dynamical theory or chaos theory has been extensively studied and used in the forecasting of air pollutant concentrations (Raga and LeMoyne, 1996; Li et al., 1994; Chen et al., 1998; Kocak et al., 2000). The basic assumption involved in the application of these techniques is that the single air pollutant concentration time series contains the effect of all the influencing factors. The major advantage of these techniques over ARIMA is their ability to take into account the nonlinear dynamics involved in the time series. The information about the linearity or nonlinearity of the time series is however, not available in advance. So one applies the number of linear and nonlinear models and finally selects the one, which provides the most accurate results. Also, the real time series contains patterns that are not purely linear or nonlinear and applying either of the techniques may give inadequate results. Hence the models need to be developed that consider both linearity and nonlinearity involved in the time series. This helps in improving the forecasting ability of the model. Also, combining different models can increase the chance to capture different patterns in the data and improve forecasting performance (Clemen, 1989; Newbold and Granger, 1974). In this study, a hybrid methodology is therefore proposed to tackle the problem of modeling the air pollutant time series with linear and nonlinear patterns. For this, the concepts from ARIMA model and nonlinear dynamical systems theory are utilized. The proposed technique is applied to the time series of nitrogen dioxide (NO2) concentrations in ambient air measured at a site in Delhi during 1999–2003. In order to compare the forecasting efficiency of the proposed hybrid model, ARIMA and nonlinear models are also developed individually and the results of these models are then compared with the hybrid model.
1775
2. Box–Jenkins ARIMA models ARIMA linear models have dominated many areas of time series forecasting. As the application of these models is very common, it is described here briefly. In general, a nonseasonal time series, xt;t¼1...n (n being the number of observations) of air pollutant concentrations measured at an equal time intervals, can be modeled as a combination of past values and past errors as xt ¼ a1 xt1 þ a2 xt2 þ þ ap xtp þ et b1 et1 b2 et2 bq etq , ð1Þ where a and b are the coefficients, p and q are the order of the autoregressive and moving average polynomials, respectively. The further details to estimate the parameters and order of the model are given in Box and Jenkins (1970). 3. Nonlinear dynamical modeling The nonlinear dynamical modeling involves the reconstruction of phase space of the time series to describe the behavior of a nonlinear system. A phase space is an abstract construct whose coordinates are the components of the state (Cambel, 1993). In general, phase space is nothing but the collection of all possible variables underlying the system. The phase space portrait can be analyzed mathematically to demonstrate the presence of an attractor and its dimension. An attractor characterizes the long-term behavior of the system in the phase space (Martelli, 1999). If an attractor exists, then the minimum number of independent variables describing the system can be estimated by computing the dimension of the attractor. The general nonlinear prediction method is to reconstruct the phase space from the set of data in a minimum embedding space and then predict the future using a local approximation function computed from the set of given data (Farmer and Sidorowich, 1987; Abarbanel et al., 1993; Takens, 1981). According to Takens’ embedding theorem, the predictions can be obtained from the set of previous data points using the functional relationship X nþT ¼ f ðX n Þ,
(2)
where X n is a vector of data points defined by X n ¼ ðxn ; xnt ; xn2t ; . . . ; xnðm1Þt Þ,
(3)
ARTICLE IN PRESS A.B. Chelani, S. Devotta / Atmospheric Environment 40 (2006) 1774–1780
1776
where n is the number of data points in the time series, T is the prediction lead-time, t is the time delay and m is the embedding dimension. The time lag t can be obtained by using the mutual information IðtÞ (Fraser and Swinney, 1986) defined as IðtÞ ¼
X
pðxt ; xtþt Þ ln
t;tþt
pðxt ; xtþt Þ , pðxt Þpðxtþt Þ
(4)
where pðxt Þ and pðxtþt Þ are the probabilities to find a time-series value in the tth and ðt þ tÞth interval, respectively, and pðxt ; xtþt Þ is the joint probability that an observation falls into the tth and ðt þ tÞth interval. These probabilities can be obtained by plotting the histogram of the data. This function needs to be computed for various t and that value can be considered as optimum t, where this function exhibits a first minimum. An embedding dimension m of a dynamical system is an integer that gives the necessary number of coordinates to unfold its dynamics. It can be obtained by using false nearest-neighbor method (Kennel et al., 1992). According to this, for every vector Xt, its nearest neighbors X 0t can be obtained by computing the distance jjX t X 0t jj between the two vectors. The false nearest neighbor can be obtained if, Ft ¼
jjX t ðm þ 1Þ X 0 t ðm þ 1Þjj , jjX t ðmÞ X 0t ðmÞjj
(5)
exceeds a threshold value say ‘r’. The details of this method are given in Kennel et al. (1992). After computing the appropriate values of t and m, the next step is to establish the functional relationship in Eq. (2) using local approximations. Given the vector Xn from which the predictions are to be made, one selects its nearest neighbors X 0n by using the Euclidean distance between the two vectors. With this, the local functions can then be built, which take each point in the neighborhood to the next neighborhood, i.e., X 0n to X 0nþT . The approximation of function f then can be obtained in terms of local polynomial maps, which can be expressed as X 0nþT
¼AþB
X 0n
þ
CðX 0n Þ2 ,
(6)
where A, B and C are the coefficients that are to be determined from the learning sets by using leastsquares estimation. The predicted point is then set as the new starting vector and the above process can be repeated to predict the other values.
4. The hybrid methodology The ARIMA and nonlinear models are useful for modeling linear and nonlinear time series, respectively. For modeling the nonlinear time series, ARIMA models, however, do not provide accurate results and applying the nonlinear models to linear problems is not a reasonable step. Hence hybridizing the linear and nonlinear models would give better performance as compared to applying the individual models (Gooijer and Kumar, 1992). The hybrid methodology is based on the combination of linear autocorrelation structure and nonlinear part, which can be given as yt ¼ l t þ nt ,
(7)
where lt and nt denote the linear and nonlinear component of time series xt and these are to be estimated from the data. To compute these two components, the first step is to apply ARIMA model using the procedure described above to the time-series data. The next step is to obtain the residuals of the ARIMA model. The residuals now represent the nonlinearity part of the data. Let the residuals be denoted as rt at time t as (8) rt ¼ yt l^t , where l^t is the forecast value for time t from Eq. (1). These residuals can be used to evaluate the model performance. Accepting the model implies that the linear correlations are not significant in the residuals. However, there may be some nonlinear dependence among the residuals and residual analysis would not be able to capture this nonlinearity. Hence modeling residuals using nonlinear techniques can provide an insight into the nonlinear relationships in the data. The next step in the hybrid methodology is to model residuals using the nonlinear modeling technique described in earlier sections. Let the forecast from the nonlinear model be denoted as n^ t . The time-series forecast can therefore be obtained as y^ t ¼ l^t þ n^ t .
(9)
Describing briefly, the hybrid technique consists of two steps; first an ARIMA model is used to model the linearity in the data and then a nonlinear model is developed to model the residuals from the ARIMA model. The results from the nonlinear model can be used as predictions of the residual terms of the ARIMA model. To evaluate the performance of the proposed model, the error
ARTICLE IN PRESS
where ‘observed’ denotes the observed and ‘predicted’ denotes the estimated values by the proposed model. 5. Study area and data used To examine the effectivity of proposed hybrid model over real data, the time series of nitrogen dioxide concentration observed in ambient air during 1999–2003 at a site in Delhi was utilized. Delhi, the capital city of India (Latitude 281350 N, Longitude 771150 E) with a population of around 13.8 million, is among the most polluted cities in the world. The region has a semi-arid climate that is often described as tropical plain, with extremely hot summers, heavy rainfalls in the monsoon months (approximately 73 cm) and cold winters. The selected site is termed as ‘commercial’ based on the relative activities in and around the site. The sampling frequency of NO2 was 4 hourly monitored round the clock in a day and twice a week. The daily averages were taken; hence in each year, 104 measurements were available (NEERI Report 2001,2002,2003,2004). 6. Results and discussions The time series of NO2 concentration observed during 1999–2003 is plotted in Fig. 1. The concentration of NO2 varies between 33 and 128 mg m3 with an average and standard deviation of 61 and 17.8 mg m3, respectively. Maximum concentration is generally observed in winter months from November to February with minimum concentration in monsoon season during July to October. During winter months, the levels have exceeded the regulatory limits stipulated by Central Pollution Control Board (CPCB). It can be observed from
1777
180 160 140 120 100 80 60 40 20 0 1999
2000
2001 Year
2002
2003
Fig. 1. NO2 time series observed at a site in Delhi during 1999–2003. 200 Learning set NO2 (µg m-3)
statistics such as correlation between observed and predicted concentrations, root mean square error (RMSE), mean absolute percentage error (MAPE), relative error (RE) are utilized. These test statistics can be obtained as n 1X observedt predictedt , MAPE ¼ 100x n t¼1 observedt rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 Xn ðobservedt predictedt Þ2 , RMSE ¼ t¼1 n sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Pn 2 t¼1 ðobservedt predictedt Þ , ð10Þ RE ¼ Pn1 2 t¼1 ðobservedt observedtþ1 Þ
NO2 (µg m-3)
A.B. Chelani, S. Devotta / Atmospheric Environment 40 (2006) 1774–1780
Prediction set
150 100 50 0 1999
2000 Observed
2001 Year
2002
2003
Predicted using AR(2) model
Fig. 2. Prediction performance of AR(2) model for NO2 concentration observed at a site in Delhi.
Fig. 1 that NO2 time series is stationary with oscillating characteristics. For modeling purpose, the time series is divided into two parts; learning and prediction sets. The data observed during 1999–2002 is considered as learning set and the data observed during 2003 is considered as prediction set. The codes in MATLAB (Beale, 1997) were written for all the computations. For ARIMA modeling, the order of the model is selected by plotting the autocorrelation and partial autocorrelation functions. An autoregressive model of order 2, i.e. AR(2) is found to be appropriate. With this model order, autoregressive model was fitted to the learning data. The model parameters were obtained by adopting the Box–Jenkins methodology. The one-step ahead predictions were then obtained and plotted in Fig. 2. The residual analysis was carried out to examine the significance of autocorrelations in the residuals. The insignificant autocorrelations were observed in the residuals of AR(2) model fitted to NO2 time series. RMSE, MAPE, RE and correlation between observed and predicted NO2 concentration is given in Table 1. For nonlinear modeling, phase space was reconstructed using Eq. (3). Here choice of m and t is crucial for proper unfolding of the dynamics
ARTICLE IN PRESS
Mean absolute Root mean Corr percentage error square error
ARIMA Nonlinear prediction Hybrid model
0.24 0.23
17.3 11.6
0.19
5.37
58.78 55.37
0.90 0.91
13.93
0.93
Corr indicates the correlation between observed and predicted time series.
involved in the time series. The time delay t is computed using the average mutual information function, which is plotted against various lags in Fig. 3. That lag value can be selected as optimum t, where average mutual information function exhibits the first minimum. Mutual information function I(t) shows a clear first minimum at t ¼ 13. Using this value of t, false nearest-neighbor search was carried out to compute embedding dimension m. The threshold ‘r’ is selected as equal to 15 following Abarbanel (1996), who found that for many nonlinear systems ‘r’ approaches 15. Fig. 4 shows the percentage of false nearest neighbors as a function of dimension. Minimum embedding dimension is that value for which the false nearest neighbors goes to zero. For the case of NO2 time series, the function decreases rapidly at m ¼ 6 and remains approximately constant at zero. Hence an embedding dimension of 6 should be sufficient to represent the system. After selecting the appropriate values for t and m, nonlinear modeling was performed using the above described method. The learning set (1999–2002) was used for reconstructing the phase space with an embedding dimension of 6 and time lag of 13. For example, the phase space vector Xn for n ¼ 66 can be formed as ½x66 ; x53 ; x40 ; x27 ; x14 ; x1 . The nearest neighbors X 0n of X n (X 066 is denoted as the nearest neighbor of X 66 ) were searched out in the learning set. The one-step predictions (T ¼ 1 in Eqs. (2) and (6)) were then obtained by least-squares fitting. In this way, to predict each point in the phase space, function f is approximated using local polynomial equation and the coefficients A, B, and C are learned from the learning set. The coefficients A, B and C are not the constants rather depend on the learning samples. For each point, its nearest neighbors are obtained and then projected to obtain
0.8 0.6 0.4 0.2 0.0 0
5
10
15 Time Lag
20
25
30
Fig. 3. Average mutual information function for NO2 time series.
false nearest neighbours (%)
Relative error
1.0
45 40 35 30 25 20 15 10 5 0 0
1
2
3
4
5
6
7 8 9 Time Lag
10 11 12 13 14 15
Fig. 4. Percentage of false nearest neighbours for NO2 time series.
Observed
NO2 (µg m-3)
Prediction technique
1.2
180 160 140 120 100 80 60 40 20 0
Predicted using Nonlinear model
1999
2000
Observed 180 160 140 120 100 80 60 40 20 0
2001 Year
2002
2000
2001
2003
Predicted using Hybrid model
Learining set
1999
(b)
Prediction set
Learning set
(a)
NO2 (µg m-3)
Table 1 Forecasting performance of three models for the prediction of NO2 concentration
Average Mutual Information (bits)
A.B. Chelani, S. Devotta / Atmospheric Environment 40 (2006) 1774–1780
1778
2002 Year
Prediction set
2003
2004
2005
Fig. 5. Prediction performance of (a) nonlinear and (b) hybrid model for NO2 time series observed at a site in Delhi.
the next point. The prediction results are plotted in Fig. 5a. The performance statistics were also computed for this model and presented in Table 1.
ARTICLE IN PRESS A.B. Chelani, S. Devotta / Atmospheric Environment 40 (2006) 1774–1780
For hybrid modeling, first an autoregressive model of order 2 was applied to the NO2 time series and the residuals were obtained. The local approximation model was then used to model these residuals. For this purpose, the time delay t ¼ 1 and an embedding dimension m ¼ 3 were found appropriate for the residuals of the AR(2) model. Using these measures, nonlinear modeling technique given in Eqs. (2) and (6) was applied to the reconstructed phase space of residuals and the predictions were then combined to obtain the predictions of NO2 concentrations. The observed and predicted NO2 concentration is plotted in Fig. 5b. Evaluating the prediction performance quantitatively, the correlation between observed and predicted time series of AR(2) model is lower than the other two models. Considering the RMSE, MAPE and RE, it is again evident that AR(2) model does not perform well as compared to other two models for NO2 concentration prediction. Nonlinear model or local approximation model, on the other hand, provides the forecasts with less RE and MAPE than AR (2) model. However, comparing the performance of nonlinear and hybrid model, the error statistics is lower for hybrid model than the nonlinear model. Hybrid model results indicate that it outperforms the AR(2) and nonlinear model. In order to assess the model performance in case of the high concentrations, peaks in December 2000 and January 2001 were selected. The NO2 concentrations exceeded the CPCB guideline of 80 mg m3 for these peaks. The observed concentrations for these peaks were 150 and 155 mg m3, where as the fitted concentrations were 121 and 127 mg m3 by AR(2) model, 129 and 145 mg m3 by nonlinear model and 137 and 148 mg m3 by hybrid model. Considering the high concentration of 136 mg m3 in prediction set during December 2003, the predicted concentration by AR(2), nonlinear and hybrid model was 115, 122 and 129 mg m3, respectively. It can be noted that AR(2) model does not perform well for high concentration as it underestimate the levels, whereas nonlinear and hybrid models perform fairly well for high concentration levels. As the hybrid model gives sensible results with less error than individual AR(2) and nonlinear models for the peak NO2 concentration, this type of modeling approach can be useful for the cases where lower forecasting error is desired. It is also tried to obtain the predictions for 2004 and 2005 to check the capability of the model to produce the next year’s forecast. Although this
1779
exercise was performed for all the three models, only the results of hybrid model are shown. Also, as the observed data was not available to validate the results, only the forecasts are plotted in Fig. 5b. It can be observed that the forecasts follow nearly the same variations as previous years. As and when the latest data are available, the model predictions can further be validated. 6.1. Effect of prediction horizon on prediction performance In order to examine the model sensitivity towards multi-step forecasting, two forecast horizons of periods 8 (one month ahead forecast) and 24 (3 months ahead forecast) were used. The error statistics MAPE and RE were used for the comparison of prediction performance. It can be observed form Table 2 that applying nonlinear model alone can improve the forecasting accuracy over AR(2) model for the 8-period prediction horizon, whereas the performance of the nonlinear model is getting worse as the time horizon extends to 24 periods. This may suggest that neither the nonlinear model nor the AR(2) model captures all of the patterns in the data. On the other hand, the MAPE and RE for the hybrid model are lower as compared to individual models. The results of the hybrid model show that by combining two models together, the overall forecasting errors can be significantly reduced. In terms of RE, the hybrid model improved over the AR(2) and nonlinear model with a reduction of 26.31% and 21.05%, respectively for one step forecast. For 8-period forecast, the improvement of hybrid model over AR(2) and nonlinear model is 22.22% and 7.4%, respectively and for 24-period forecasts, the reduction is 19.35% and 9.6%, respectively. This shows that hybrid model outperforms the individual models for multi-step forecasting. Table 2 Effect of prediction horizon on forecasting performance of three models Prediction technique
8-period ahead
24-period ahead
RE
MAPE
RE
MAPE
22.7 19.3
0.37 0.34
26.3 29.8
0.31
13.69
ARIMA 0.33 Nonlinear 0.29 prediction Hybrid model 0.27
9.58
ARTICLE IN PRESS 1780
A.B. Chelani, S. Devotta / Atmospheric Environment 40 (2006) 1774–1780
7. Conclusions In this paper, a hybrid approach is proposed to forecast the air pollutant concentrations. This is achieved by extracting the unique features of linear autoregressive model and nonlinear model based on the chaos theory. The hybrid model is applied to the time series of NO2 concentration observed at a site in Delhi. The autoregressive and nonlinear models were also applied in order to compare the results of the hybrid model. The one-step ahead, 8-period and 24-period ahead predictions were obtained to evaluate the model’s capabilities in forecasting for different time horizons. The univariate approaches applied above assess the characteristics of the time series and provide the predictions without having an understanding of the mechanisms that govern the system. The prediction performance results show that the hybrid modeling can be an effective tool to forecast the air pollutant concentrations instead of applying individual models. The study is useful for the cases where the data on other explanatory variables that influence the air pollutant concentrations is not available. The developed model can also be applied to predict other pollutants like ozone concentrations. References Abarbanel, H.D.I., Brown, R., Sidorowich, J.J., Tsimring, L.Sh., 1993. The analysis of observed chaotic data in physical systems. Review of Modern Physics 65, 1331–1392. Beale, M., 1997. MATLAB Users Guide (Version 5.1). The Math Works Inc, Natick, MA. Benarie, M., 1987. The limits of air pollution modeling. Atmospheric Environment 21 (1), 1–5. Box, G.E.P., Jenkins, G.M., 1970. Time Series Analysis, Forecasting and Control. Holden-Day, San Francisco, CA. Cambel, A.B., 1993. Applied Chaos Theory—a paradigm for complexity. Academic Press Inc., San Diego, CA. Chelani, A.B., Gajghate, D.G., Hasan, M.Z., 2002. Prediction of ambient PM10 and toxic metals using artificial neural networks. Journal of the Air and Waste Management Association 52, 805–810. Chen, J.L., Islam, S., Biswas, P., 1998. Nonlinear dynamics of hourly ozone concentrations: nonparametric short-term prediction. Atmospheric Environment 32 (11), 1839–1848. Clemen, R., 1989. Combining forecasts: a review and annotated bibliography with discussion. International Journal of Forecasting 5, 559–608. Farmer, D.J., Sidorowich, J.J., 1987. Predicting chaotic time series. Physics Review Letters 59, 85–848.
Fraser, A.M., Swinney, H.L., 1986. Independent coordinates for strange attractors from mutual information. Physics Review A 33, 1134–1140. Gardner, M.W., Dorling, S.R., 1998. Artificial neural networks (the multilayer perceptron)—a review of applications in the atmospheric sciences. Atmospheric Environment 32 (14/15), 2627–2636. Gardner, M.W., Dorling, S.R., 1999. Neural network modeling and prediction of hourly NOx and NO2 concentrations in urban air in London. Atmospheric Environment 33, 709–719. Gooijer, J.G., Kumar, K., 1992. Some recent developments in non-linear time series modeling, testing, and forecasting. International Journal of Forecasting 8, 135–156. Kennel, M.B., Brown, R., Abarbanel, H.D.I., 1992. Determining embedding dimension for phase space reconstruction using a geometric method. Physics Review A 45, 3403–3411. Kocak, K., Saylan, L., Sen, O., 2000. Nonlinear time series prediction of O3 concentration in Istanbul. Atmospheric Environment 34, 1267–1271. Li, I.F., Biswas, P., Islam, S., 1994. Estimation of dominant degrees of freedom for air pollutant concentration data: applications to ozone measurement. Atmospheric Environment 28, 1707–1714. Martelli, M., 1999. Introduction to Discrete Dynamical Systems and Chaos. Wiley Inter Science Series in Discrete Mathematics and Optimization. Wiley, New York. Milionis, A.E., Davies, T.D., 1994. Regression and stochastic models for air pollution—I, review, comments and suggestions. Atmospheric Environment 28 (17), 2801–2810. NEERI Report, 2001. Air Quality Status for ten cities of India. NEERI Nagpur, India, No. 9–10. NEERI Report, 2002. Air Quality Status for ten cities of India. NEERI Nagpur, India, No. 11. NEERI Report, 2003. Air Quality Status for ten cities of India. NEERI Nagpur, India, No. 12. NEERI Report, 2004. Air Quality Status for ten cities of India. NEERI Nagpur, India, No. 13. Newbold, P., Granger, C.W.J., 1974. Experience with forecasting univariate time series and the combination of forecasts (with discussion). Journal of the Royal Statistical Society Series A 137, 131–164. Raga, G.B., LeMoyne, L., 1996. On the nature of air pollution dynamics in Mexico city—I, non linear analysis. Atmospheric Environment 30 (23), 3987–3993. Schlink, U., Herbarth, O., Tetzlaff, G., 1997. A component timeseries model for SO2 data: forecasting, interpretation and modification. Atmospheric Environment 31, 1285–1295. Shi, J.P., Harrison, R.M., 1997. Regression modeling of hourly NOx and NO2 concentrations in urban air in London. Atmospheric Environment 31 (24), 4081–4094. Takens, F., 1981. Detecting strange attractors in turbulence. In: Rand, D.A., Young, L.S. (Eds.), Lecture Notes in Mathematics, vol. 898. Springer, New York. Zennetti, P., 1990. Air Pollution Modeling—Theories, Computational Methods and Available Software. Computational Mechanics Publications.