Time-series analysis of monthly rainfall data for the Mahanadi River Basin, India

Janhabi Meher , Ramakar Jha

Department of Civil Engineering, National Institute of Technology (NIT Rourkela), Orissa, India

1 Introduction

Rainfall is a natural phenomenon resulting from atmospheric and oceanic circulation (local convection, frontal or orographic patterns). In the Mahanadi River Basin in India,which is a sub-tropical/semi-arid region, the prediction of rainfall is extremely important for proper mitigation and management of floods, droughts, environmental flows, water demand by different sectors, maintaining reservoir levels,and disasters.

Many attempts have been made in the recent past to model and forecast rainfall using various techniques, with the use of time series techniques proving to be the most common (Gorman and Toman, 1966; Salaset al., 1980;Galeati, 1990; Lall and Bosworth, 1993; Hsuet al., 1995;Davidsonet al., 2003). In time series analysis it is assumed that the data consist of a systematic pattern (usually a set of identifiable components) and random noise (error) which usually makes the pattern difficult to identify. Time series analysis techniques usually involve some method of filtering out noise in order to make the pattern more salient. The time series patterns can be described in terms of two basic classes of components: trend and seasonality. The trend represents a general systematic linear or (most often) nonlinear component that changes over time and does not repeat or at least does not repeat within the time range captured by the data.The seasonality may have a formally similar nature; however, it repeats itself in systematic intervals over time. Those two general classes of time series components may coexist in real-life data.

The present work exclusively deals with a time series forecasting model, in particular, the autoregressive inte-grated moving average (ARIMA) model. These models were described by Box and Jenkins (1970) and further discussed in other resources such as Montgomery and Johnson (1967), Vandaele (1983), and Chatfield (1996). In fact, the ARIMA model is an important forecasting tool,and is the basis of many fundamental ideas in time-series analysis. An autoregressive model of orderpis conventionally classified as AR(p) and a moving average model withqterms is known as MA(q). A combined model that containspautoregressive terms andqmoving average terms is called ARMA(p,q) (Gujarati, 1995). If the object series is differenceddtimes to achieve stationarity, the model is classified as ARIMA(p,d,q), where the letter "I"signifies "integrated". Thus, an ARIMA model is a combination of an autoregressive (AR) process and a moving average (MA) process applied to a non-stationary data series.

The ARIMA model possesses many appealing features.It allows a researcher who has data only on past years (e.g.,rainfall), to forecast future events without having to search for other related time series data (e.g., temperature). The Box-Jenkins approach also allows for the use of one time series, for example temperature, to explain the behavior of another series, for example rainfall, if the other time series data are correlated with a variable of interest and if there appears to be some cause for this correlation. ARIMA modeling has been successfully applied in various water and environmental management applications and some of the relevant literature is discussed below.

Chiewet al.(1993) conducted a comparison of six rainfall-runoff modeling approaches to simulate daily, monthly and annual flows in eight unregulated catchments. They concluded that the time-series approach can provide adequate estimates of monthly and annual yields in the water resources of the catchments.

Kuo and Sun (1993) employed an intervention model for average 10 days stream flow forecast and synthesis which was investigated to deal with the extraordinary phenomena caused by typhoons and other serious abnormalities of the weather of the Tamshui River Basin in Taiwan.

Time series analysis was used by Langu (1993) to detect changes in rainfall and runoff patterns and to search for significant changes in the components of a number of rainfall time series.

Brathet al.(2002) applied the ARIMA model for forecasting short-term future rainfall of the Sieve River Basin in central Italy for improving the real-time flood forecasts issued by deterministic lumped rainfall-runoff models.

Tsenget al.(2002) proposed a hybrid forecasting model,which combines seasonal time series ARIMA (SARIMA)and neural network back propagation (BP) models, known as SARIMABP. These models were used to forecast two seasonal time series datasets of total production value for the Taiwan machinery industries.

Kihoroet al.(2004) compared the performance of artificial neural networks (ANNs) and ARIMA models in forecasting of seasonal (monthly) rainfall time series. A rule of selecting input lags into the input set based on their relevance/contribution to the model was also proposed.

Yurekli and Kurung (2004) used linear stochastic models(ARIMA) to simulate drought periods in hydrologically homogeneous regions. Predicted data were compared to the observed data using the best models over a period of five years. The results showed that the predicted data represented the actual data very well for each hydrologically homogeneous region.

Mohammadiet al.(2005) used three different methods(artificial neural network (ANN), ARIMA time series, and a regression analysis between some hydro-climatological data and inflows) to predict spring inflows. The spring inflows account for almost 60% of the annual inflow to the Amir Kabir Reservoir in Iran. Twenty-five years of observed data were used to train or calibrate the models and five years were applied for testing.

Weesakul and Lowanichchai (2005) used an ARIMA model to fit the time series of annual rainfall during 1951–1990 of 31 rainfall stations distributed in all regions of Thailand. This study found that the ARIMA model was more suited to describe the inter-annual variation of annual rainfall in Thailand,i.e., most of the rainfall stations were better fitted with the ARIMA model, while only eight stations were better fitted with ARMA model.

Somvanshiet al.(2006) presented tools for modeling and predicting the behavioral patterns in rainfall phenomena based on past observations. This paper introduced two fundamentally different approaches for designing a model, the statistical method based on ARIMA and the emerging computationally powerful techniques based on ANN. In order to evaluate the prediction efficiency, 104 years of mean annual rainfall data from 1901 to 2003 of the Hyderabad region (in India) were analyzed. The models were trained with 93 years of mean annual rainfall data. The ANN and the ARIMA approaches were applied to the data to derive the weights and the regression coefficients respectively. The performance of the model was evaluated by using the remaining 10 years of data.

Naill and Momani (2009) used ARIMA to model monthly rainfall data for Amman, Jordan and developed an ARIMA model to help forecast monthly rainfall. This model was used for forecasting the monthly rainfall for the upcoming 10 years to help decision makers establish priorities in terms of water demand management.

Otok and Suhartono (2009) investigated the best method to analyze rainfall index data in Indonesia by comparing the forecast accuracy among ARIMA, ASTAR, and single-input transfer function, and multi-input transfer function models,using both non-seasonal ARIMA and seasonal ARIMA(SARIMA) models to forecast rainfall in Indonesia. This study found that the SARIMA model performed better than the non-seasonal ARIMA model.

Rabenjaet al.(2009) forecasted both monthly rainfall and discharge of the Namorona River in the Vohiparara Riv-er Basin of Madagascar using ARIMA and SARIMA models and also concluded that the SARIMA model was more suitable than the non-seasonal ARIMA model.

Abuduet al.(2010) presented the application of ARIMA,SARIMA, and Jordan-Elman artificial neural networks(ANN) models in forecasting the monthly streamflow of the Kizil River in Xinjiang, China. In this study, two different types of monthly streamflow data (original and deseasonalized data) were used to develop time series and Jordan-Elman ANN models using previous flow conditions as predictors.

Tularam and Ilahee (2010) examined a large data set involving more than 50 years of rainfall and temperature data using spectral analysis ARIMA methodology to analyze climatic trends and interactions. Fourier analysis, linear regression and ARIMA based time-series models were used to analyze the large data sets using MATLAB, SPSS and SAS software.

Chattopadhyay and Chattopadhyay (2010) forecasted the Indian Summer Monsoon Rainfall (ISMR) using a univariate ARIMA model.

Mauludiyantoet al.(2010) modeled tropical rain attenuation in Surabaya, Indonesia adopting the ARIMA model.

Helman (2011) analyzed of monthly average temperature and precipitation sum time series recorded at 44 measurement stations in the Czech Republic over the period of 1961–2008. The two objectives of this study were the construction of SARIMA models based on Box-Jenkins methodology and a comparison of different models constructed according to the given factors of particular measurement stations’ elevation, longitude and latitude.

Shamsniaet al.(2011) analyzed 20 years of statistics on relative humidity, monthly average temperature and precipitation of the Abadeh Station, Iran using ITSM time series analysis software. According to the ARIMA model, ACF,PACF and evaluation of all the eventual samples, an ARIMA precipitation model of monthly average temperature and relative humidity was obtained.

Babuet al.(2011) used rainfall flow data from a meteorological station located in Vellore in Tamil Nadu, India, to determine mean monthly flow by making use of an autoregressive approach. This approach can be used for regenerating the future sequence while preserving the inherited properties of the observed data. Their comparison of the observed rainfall flow and the synthetically generated data showed that the statistical characteristics were satisfactorily preserved.

Mahsinet al.(2012) used the Box-Jenkins methodology to build a SARIMA model for monthly rainfall data from the Dhaka, Bangladesh Station for the period of 1981–2010 with a total of 354 readings. In that paper, the SARIMA model was found adequate and was used to forecast the monthly rainfall for the upcoming two years to help decision makers establish priorities for water demand management.

All of these models discussed above used mean rainfall figures from across a basin, obtained using various geo-statistical methods, for simulating and forecasting rainfall data. In the present work, a generic ARIMA model has been developed for (a) simulating and forecasting mean rainfall, obtained using Theissen weights, over the Mahanadi River Basin in India, and (b) simulating and forecasting mean rainfall for 38 district towns’ rain-gauge stations, independently. Various error statistics are used to test the validity and applicability of the developed ARIMA model.

2 Study area and data collection

The Mahanadi River is one of the major inter-state east flowing rivers in peninsular India (Figure 1). It originates at an elevation of about 442 m a.s.l. from the Amarkantak Hills of the Bastar Plateau near Pharasiya village in Raipur district of Chhatisgarh, India. During the course of its traverse, it drains fairly large areas of the state of Chhatisgarh-Orissa and comparatively small areas in the state of Jharkhand-Maharashtra. The Mahanadi River Basin is located at 80°30'E–86°50'E and 19°20'N–23°35'N.

The total catchment area of the river is 141,134 km2with mean annual river flow of 59,155 million cubic meters. The average annual discharge is 1,895 m3/s, with a maximum of 6,352 m3/s during the summer monsoon.Minimum discharge is 759 m3/s and occurs during the months from October to next May. The river passes through the tropical zone and is subjected to cyclonic storms and seasonal precipitation. The monsoon season usually begins in the second half of June, and lasts until the first half of September. Normal annual precipitation of the basin is 1,360 mm of which about 86% (1,170 mm) occurs during the monsoon season. Monthly rainfall data from the 38 district town rain-gauge stations, as shown in Figure 2,were collected for the period of 1901–2002, with a total of 1,224 data sets at each station.

We made time series plots and computed basic statistics to understand the trends, seasonality and statistical variations at each district town rain-gauge station. Figure 3 shows the time-series plots of two stations, as representative stations among the 38 in the study area. Some basic statistics of all the stations are shown in Figure 4. The standard deviations in rainfall were high in the district towns away from coastal areas, and the maximum rainfall amount also increased inland away from the coastal areas. However, the mean annual rainfall and coefficient of variation values were almost constant except at the town of Bilaspur, where the rainfall was significantly lower.

Figure 4 shows that the mean monthly rainfall ranged from 61.37 mm in Bilaspur to 121.88 mm in Kendrapada with a coefficient of variation ranging from 107.10 in Puri to 143.70 in Jashpur. Maximum monthly rainfall ranged from 424.88 mm in Puri during July, 1951 to 763.45 mm in Jashpur during August, 1943. During the 102-year study period the highest monthly rainfalls in the month of July were observed in 1903, 1929, and 1956.

Figure 1 Mahanadi River Basin Map

Figure 2 Location of district towns having rain-gauge stations

Figure 3 Monthly rainfall data of Angul (a) and Bastar (b)

Figure 4 Some basic characteristics of rainfall in all the district towns of the Mahanadi River Basin: (a) maximum rainfall;(b) standard deviation in rainfall data; (c) mean rainfall; (d) coefficient of variation in rainfall data

3 Methodology

3.1 ARIMA model

The ARIMA model is an extension of the ARMA model in the sense that by including auto-regression and moving average it has an extra function for differencing the time series. If a dataset exhibits long-term variations such as trends, seasonality and cyclic components, differencing a dataset in ARIMA allows the model to deal with them. Two common processes of ARIMA for identifying patterns in time-series data and forecasting are auto-regression and moving average.

3.1.1 Autoregressive process

Most time series consist of elements that are serially dependent in the sense that one can estimate a coefficient or a set of coefficients that describe consecutive elements of the series from specific, time-lagged (previous) elements. Each observation of the time series is made up of random error components (random shock, ε) and a linear combination of prior observations.

3.1.2 Moving average process

Independent from the autoregressive process, each element in the series can also be affected by the past errors (or random shock) that cannot be accounted for by the autoregressive component. Each observation of the time series is made up of a random error component (random shock, ε)and a linear combination of prior random shocks.

3.1.3 General form of non-seasonal and seasonal ARIMA

ARIMA models are sometimes called Box-Jenkins models.An ARIMA model is a combination of an autoregressive (AR) process and a moving average (MA) process applied to a non-stationary data series. As such, in the general non-seasonal ARIMA-(p,d,q) model, AR-(p)refers to order of the autoregressive part, I-(d)refers to degree of differencing involved and MA-(q)refers to order of the moving average part. The equation for the simplest ARIMA-(p,d,q)model is:

or in a general form:

where,φirefers toith term autoregressive parameter,θirefers toith term moving average parameter,cmeans constant,emeans error at timet,Bprefers topth order backward shift operator, andXtrefers to time series value at timet.

Seasonal ARIMA (SARIMA) is a generalization and extension of the ARIMA method in which a pattern repeats seasonally over time. In addition to the non-seasonal parameters, seasonal parameters for a specified lag (established in the identification phase) need to be estimated. Analogous to the simple ARIMA parameters, these are: seasonal autoregressive (P), seasonal differencing (D), and seasonal moving average parameters (Q). The seasonal lag used for the seasonal parameters is usually determined during the identification phase and must be explicitly specified. In addition to the non-seasonal ARIMA-(p,d,q) model introduced above, we could identify SARIMA-(P,D,Q) parameters for our data. The general form of the SARIMA-(p,d,q)-(P,D,Q)Smodel using backshift notation is given by:

wheresrefers to number of periods per season,ΦA(chǔ)Rrefers to non-seasonal autoregressive parameter,ΦSARrefers to seasonal autoregressive parameter,θMArefers to non-seasonal moving average parameter, andθSMArefers to seasonal moving average parameter.

Four phases are involved in identifying patterns of time series data using non-seasonal and seasonal ARIMA. These are: model identification, parameter estimation, diagnostic checking and forecasting. The first step is to determine if the time series is stationary and if there is any significant seasonality that needs to be modeled. Mahanadi River Basin mean rainfall data obtained using the Theissen method and monthly rainfall data of all 38 district rain-gauge stations were examined to check for the most appropriate class of ARIMA processes through selecting the order of the consecutive and seasonal differencing required to make the series stationary.

We identified the stationary component of a data set by performing the Ljung and Box test. We tested this hypothesis by choosing a level of significance for the model adequacy and compared the computed Chi-square (Χ2) values with the Χ2values obtained from the table. If the calculated value is less than the actual Χ2value, then the model is adequate, otherwise not. TheQ(r) statistic is calculated by the following formula:

wherenis the number of observations in the series andr(j)is the estimated correlation at lagj.

Furthermore, we tested the data to specify the order of the regular and seasonal autoregressive and moving average polynomials necessary to adequately represent the time series model. For this purpose, model parameters were estimated using a maximum likelihood algorithm that minimized the sums of squared residuals and maximized the likelihood (probability) of the observed series. The maximum likelihood estimation is generally the preferred least square technique.

The major tools used in the identification phase are plots of the series, correlograms (plots of autocorrelation and partial autocorrelation verses lag) of the autocorrelation function (ACF) and the partial autocorrelation function (PACF).The ACF and the PACF are the most important elements of time series analysis and forecasting. The ACF measures the amount of linear dependence between observations in a time series that are separated by a lagk. The PACF plot helps to determine how many autoregressive terms are necessary to reveal one or more of the following characteristics: time lags where high correlations appear, seasonality of the series, and trend either in the mean level or in the variance of the series.

In diagnostic checking, the residuals from the fitted model were examined against their adequacy. This is usually done by correlation analysis through the residual ACF plots and by goodness-of-fit test using means of Chi-square statistics.

At the forecasting stage, the estimated parameters were used to calculate new values of the time series with their confidence intervals for the predicted values.

3.2 Performance evaluation

To choose the best model among the class of plausible models, we used the Akaike Information Criterion (AIC)proposed by Akaike (1974). The model which had the minimum AIC value was considered as the best model for the Mahanadi River Basin rainfall analysis of cases (a) and (b)above.

wherek=1 ifc≠0 and 0 otherwise, andLis the maximized likelihood of the model fitting to the differenced data(1-Bs)D(1-B)dXt.

At the forecasting stage, the estimated parameters were tested for their validity using error statistics such as coefficient of determination (R2), mean square error (MSE), and mean absolute error (MAE) criteria.

4 Results and discussion

The time series analyses of the mean rainfall (obtained using Theissen weights), over the Mahanadi River Basin and the mean rainfall for each of the 38 district towns, were done separately. A total of 1,224 data sets obtained during the period of 1901–2002 were used for the analysis. The results indicate that the data sets were non-stationary in nature (in terms of both mean and variance) and they reflected seasonal cycles. This was confirmed: when the ACF and PACF plots of the original data were prior to any transformation and differencing, they were obtained (Figure 5).

Figure 5 ACF (a, c) and PACF (b, d) obtained from observed data

In order to fit an ARIMA model, a stationary series (in terms of both in mean and variance) is needed. To establish the stationarity of the variance of the time series, a Box-Cox power transformation (α=0.5) was applied to all the data sets of case (a) and case (b) above. Stationarity of the mean could be attained by differencing the series. Differencing for non-seasonal ARIMA was not done due to no trends being detected in the data sets. However, for SARIMA, first difference (D=1) of the original data was done in order to establish stationarity in the series with no seasonal impact. The ACF and PACF plots for the differenced series were obtained again to check the stationary (Figure 6). The figure confirms that the ACF and PACF plots for the differenced and de-seasonalized rainfall data were nearly stable and the ARIMA model (p,0,q)(P,1,Q)12could be identified for further analysis.

In the next step, model parametersp,q,PandQwere identified. The ACF and PACF plots of the ARIMA model(p,0,q)(P,1,Q)12with first order seasonal differencing (Figure 5) suggested that at the initial stage the tentative model should be (1,0,1)(1,1,1)12because there was one autoregressive and one moving average parameter in the plots.We found similar condition for all the data sets of cases (a)and (b).

In ARIMA modeling it is necessary to minimize the sum of squared residuals (SSR) between the actual and estimated values to represent the data most appropriately. It is important to note that: the best ARIMA model should have the least number of parameters to acquire the minimum AIC along with the minimum SSR. Therefore, in the stage of identifying the number of autoregressive and moving average parameters, an ARIMA model (p,0,q)(P,1,Q)12with the least number of parameters was attempted. We evaluated nine different ARIMA models (p,0,q)(P,1,Q)12indicating very low SSR values for all the cases to obtain the best model among them (Table 1).

To choose the best model among the nine models selected above, the AIC shown in Equation(5)was used. The model which had the minimum AIC value was considered as the best model for Mahanadi River Basin rainfall analysis of cases (a) and (b).

As mentioned earlier, a total of 1,224 observations were used, of which 1,080 observations obtained for the years 1901-1990 were used for model calibration and 144 observations obtained for the years 1991-2002 were used for forecasting. The AIC values for all nine models were estimated using Equation(5); the ARIMA model (1,0,0)(0,1,1)12provided the best results (with minimum AIC values) for case (a) and case (b) (Figure 7). The results obtained were highly satisfactory. At some stations, the ARIMA model(0,0,1)(0,1,1)12also showed minimum AIC values similar to the ARIMA model (1,0,0)(0,1,1)12. However, ARIMA model (0,0,1)(0,1,1)12was not applicable for all the cases and was not a generic model.

Figure 6 ACF (a, c) and PACF (b, d) plots after transformation and differencing

Table 1 ARIMA models with different numbers of parameters

As discussed earlier, the ARIMA (1,0,0)(0,1,1)12model could be written in the following form:

After estimating the parametersp,q,PandQ, 1,080 observations were obtained for the years 1901-1990 and were used for model calibration. Figure 8 illustrates the comparison of results of two representative rain gauge stations Angul and Bastar obtained using the ARIMA model(1,0,0)(0,1,1)12.

The goodness-of-fit of the ARIMA model(1,0,0)(0,1,1)12was tested using the Ljung-Box statistic as shown in Equation(4). The goodness of fit values for the autocorrelations of residuals from the (1,0,0)(0,1,1)12model up to lag 36 was ≥0.05 for all 38 district towns and for the mean rainfall obtained by the Theissen method. These results substantiate the acceptance of the null hypothesis of model adequacy at the 5% significance level and the set of autocorrelations of residuals was considered white noise.

Time-series plots and normal probability plots of residuals resulting from the ARIMA model (1,0,0)(0,1,1)12were essential to find the existence of any correlation between the residuals. Time-series plots of the residuals showed a few outliers (Figure 9), and the normal probability plots of the residuals showed a few departures from the expected normal value due to the outliers (Figure 10). All the tests and verifications resulted in successful calibration of the ARIMA model (1,0,0)(0,1,1)12.

Figure 7 Minimum AIC values obtained for all 38 district town rain gauge stations in the Mahanadi Basin

Figure 8 Calibration results of ARIMA model (1,0,0)(0,1,1)12

Figure 9 Time series plot of residuals in two representative district towns

The ARIMA model (1,0,0)(0,1,1)12was also tested for its validity to forecast 144 observations obtained for the years 1991-2002 for cases (a) and (b). The results obtained using the model described above are shown in Figure 11. The observed mean rainfall was found to be closely aligned to the forecasted values of mean rainfall that were calculated for each district town, as well as to the mean rainfall for the Mahanadi River Basin. From the results presented in this study, it is apparent that the chosen model should be sufficiently accurate to forecast rainfall in this region.

Figure 10 Normal probability plot of residuals in two representative district towns

Various error statistics such as coefficient of determination (R2), mean square error and mean absolute error were used to test the validity and applicability of the ARIMA model (1,0,0)(0,1,1)12. Figure 12 illustrates all of the error statistics. TheR2value ranged from 0.682 for Bilaspur to 0.88 for Sambalpur. The lowestR2value of 0.68 for Bilaspur can be justified with its lowest mean value of 61.36 mm for monthly rainfall data. MSE ranged from 2,055.44 to 4,569.39 and MAE ranged from 27.79 to 42.95. The highest MSE and second highest MAE were in Jashpur, which make sense because Jashpur had the highest coefficient of variation, 143.70, for monthly rainfall data.

Figure 11 Forecasting of mean rainfall using developed ARIMA model (1,0,0)(0,1,1)12

Figure 12 Error statistics of ARIMA model (1,0,0)(0,1,1)12 (a) coefficient of determination,(b) mean square error, (c) mean absolute error

5 Conclusions

Time series analysis is an important tool in modeling and forecasting rainfall data. In this study we used the ARIMA model to simulate and forecast (a) mean rainfall, obtained using Theissen weights, over the Mahanadi River Basin and,(b) mean rainfall for 38 district towns. The ARIMA model(1,0,0)(0,1,1)12was developed considering step-wise analysis, non-seasonal and seasonal parameters, and various diagnostic checks. Interestingly, one unique ARIMA model fits both cases (a) and (b). The forecasting results for the upcoming 12 years are considered to be excellent and accurate. This will certainly assist policy makers and decision makers in planning for any kind of disaster or extreme condition in every district town of the Mahanadi River Basin by generating scenarios for the next few years.

The authors would like to thank National Institute of Technology, Rourkela for providing all the facilities needed for preparing this manuscript.

Abudu S, Cui CL, King JP, Abudukadeer K, 2010. Comparison of performance of statistical models in forecasting monthly stream flow of Kizil River, China. Water Science and Engineering, 3(3): 269–281.

Akaike H, 1974. A new look at the statistical model identification. IEEE Transactions on Automatic Control, 19(6): 716–723.

Babu SKK, Karthikeyan K, Ramanaiah MV, Ramanah D, 2011. Prediction of rain-fall flow time series using Auto-Regressive Models. Advances in Applied Science Research, 2(2): 128–133.

Box GEP, Jenkins GM, 1970. Time Series Analysis, Forecasting and Control,Holden-Day, San Francisco, CA.

Brath A, Montanari A, Toth E, 2002. Neural networks and non-parametric methods for improving real time flood forecasting through conceptual hydrological models. Hydrology and Earth System Sciences, 6(4):627–640.

Chatfield C, 1996. The Analysis of Time Series: An Introduction, 5th Ed.Chapman and Hall, London.

Chattopadhyay S, Chattopadhyay G, 2010. Univariate modelling of summer-monsoon rainfall time series: Comparison between ARIMA and ARNN. Comptes Rendus Geoscience, 342(2): 100–107.

Chiew FHS, Stewardson MJ, Mcmahon TA, 1993. Comparison of six rainfall-runoff modeling approaches. Journal of Hydrology, 147(1–4): 1–36.

Davidson JW, Savic DA, Walters GA, 2003. Symbolic and numerical regression: experiments and applications. Journal of Information Science,150(1–2): 95–117.

Galeati G, 1990. A comparison of parametric and non-parametric methods for runoff forecasting. Journal of Hydrological Sciences, 35(1–2): 79–94.

Gorman JW, Toman RJ, 1966. Selection of variables for fitting equation to data. Technometrics, 8(1): 27–51.

Gujarati DN, 1995. Basic Econometrics, 5th Ed. McGraw-Hill Book Co.,New York.

Helman K, 2011. SARIMA models for temperature and precipitation time series in the Czech Republic for the period 1961–2008. Journal of Applied Mathematics, 4(3): 281–290.

Hsu K, Gupta, HV, Sorooshian S, 1995. Artificial neural network modeling of the rainfall-runoff process. Water Resources Research, 31(10):2517–2530.

Kihoro JM, Otieno RO, Wafula C, 2004. Seasonal time series forecasting: A comparative study of ARIMA and ANN models. African Journal of Science and Technology, 5(2): 41–50.

Kuo JT, Sun YH, 1993. An intervention model for average 10 day stream flow forecast and synthesis. Journal of Hydrology, 151: 35–56.

Lall U, Bosworth K, 1993. Multivariate kernel estimation of functions of space and time hydrologic data. In: Hipel KW (ed.). Stochastic and Statistical Methods in Hydrology and Environmental Engineering. Springer Velag, Kluwer, New York.

Langu EM, 1993. Detection of changes in rainfall and runoff patterns.Journal of Hydrology, 147: 153–167.

Mahsin M, Akhter Y, Begum M, 2012. Modeling rainfall in Dhaka Division of Bangladesh using time series analysis. Journal of Mathematical Modelling and Application, 1(5): 67–73.

Mauludiyanto A, HendrantoroG, Purnomo MH,Ramadhany T, Matsushima A, 2010. ARIMA modeling of tropical rain attenuation on a short 28-GHz terrestrial link. IEEE Antennas and Wireless Propagation Letters,9: 223–227.

Mohammadi K, Eslami HR, Dardashti SD, 2005. Comparison of Regression,ARIMA and ANN models for reservoir inflow forecasting using snowmelt equivalent (a case study of Karaj). Journal of Agricultural Science and Technology, 7: 17–30.

Montgomery DC, Johnson LA, 1967. Forecasting and Time Series Analysis.McGraw-Hill, New York.

Naill PE, Momani M, 2009. Time series analysis model for rainfall data in Jordan: Case study for using time series analysis. American Journal of Environmental Sciences, 5(5): 599–604.

Otok BW, Suhartono, 2009. Development of rainfall forecasting model in Indonesia by using ASTAR, transfer function, and ARIMA methods.European Journal of Scientific Research, 38(3): 386–395.

Rabenja AT, Ratiarison A, Rabeharisoa JM, 2009. Forecasting of the rainfall and the discharge of the Namorona River in Vohiparara and FFT analyses of these data. Proceedings, 4th International Conference in High-Energy Physics, Antananarivo, Madagascar, pp. 1–12.

Salas JD, Deulleur JW, Yevjevich V, Lane WL, 1980. Applied Modeling of Hydrologic Time Series. Water Resources Publications, Littleton, CO.

Shamsnia SA, Shahidi N, Liaghat A, Sarraf A, Vahdat SF, 2011. Modeling of weather parameters using stochastic methods (ARIMA Model) (Case study: Abadeh Region, Iran). International Conference on Environment and Industrial Innovation, 12: 282–285.

Somvanshi VK, Pandey OP, Agrawal PK, Kalanker NV, Prakash Ravi M,Chand R, 2006. Modelling and prediction of rainfall using artificial neural network and ARIMA techniques. Journal of the Indian Geophysical Union, 10(2): 141–151.

Tseng FM, Yub HC, Tzeng GH, 2002. Combining neural network model with seasonal time series ARIMA model. Technological Forecasting &Social Change, 69: 71–87.

Tularam GA, Ilahee M, 2010. Time series analysis of rainfall and temperature interactions in coastal catchments. Journal of Mathematics and Statistics, 6(3): 372–380.

Vandaele W, 1983. Applied Time Series and Box-Jenkins Models. Academic Press, Orlando, FL.

Weesakul U, Lowanichchai S, 2005. Rainfall forecast for agricultural water allocation planning in Thailand. Thammasat International Journal of Science and Technology, 10(3): 18–27.

Yurekli K, Kurunc A, 2004. Simulation of drought periods using stochastic models. Turkish Journal of Engineering and Environmental Science, 28:181–190.