Water management policy is a crucial global issue intended to govern water acquisition and determine proper ways to minimize the misapplication and spoil of water resources. River discharge is very important in the hydrological cycle and hydro-environmental management, and plays a basic designation in designing and operating irrigation schemes at a watershed scale. One of the most important problems in analysis and management of river discharge is the modeling/simulation. For this purpose, the complexity and uncertainty of the process has been widely analyzed using artificial intelligences (AI) by hydrologists and scientists.

AI approaches such as artificial neural networks (ANN) are black box modeling approaches that have been applied in a variety of areas of hydrology, including rainfall–runoff modeling (e.g.,Hsuet al. 1995;Coulibalyet al. 1999;Anmalaet al. 2000;Tokar & Markus 2000;Nilssonet al. 2006;Aqilet al. 2007). Likewise, machine learning approaches as AI approaches, are remarkable prediction approaches which in recent years have been implemented in the river modeling field (e.g.,Azamathullaet al. 2009;Chenet al. 2010;Shiri & Kisi 2010;Sunet al. 2012;Roushangaret al. 2014a,2014b,2016a,2016b;Roushangar & Alizadeh 2015;Roushangar & Koosheh 2015;Nouraniet al. 2016). Also, some researchers have carried out reviews on good state-of-the-art application of AI and their advantages (Maier & Dandy 2000;Labat 2005;Abrahartet al. 2012;Sang 2013;Nouraniet al. 2014). The extreme learning machines (ELMs) as a new AI approach is a fast learning technique with high generalization performance that basically uses single-hidden layer feed forward neural networks (SLFNNs) (Huanget al. 2006). This technique has been successfully applied in various fields of research (Abdullahet. al. 2015;Limaet al. 2016).

Likewise, wavelet transform (WT) provides remarkable vision into the physical form of the data by presenting information in both time and frequency domains of the time series (Farajzadeh & Alizadeh 2017;Roushangar & Alizadeh 2017). The WT is an appropriate temporal pre-processing method that can be employed to extract the diversity of features from the time series, such as short-term and long-term fluctuations, by decomposing the time series into different sub-components. A time series in the WT breaks down into a series of linearly independent detail signals and one approximation signal by using discrete WT with a specific wavelet function (Foufoula-Georgiou & Kumar 1995).Mallat (1998) presented a complete theory for wavelet multi-resolution signal decomposition (also mentioned as a pyramid decomposition algorithm). Researches confirmed that proper data pre-processing by applying the WT can lead the models to adequately illustrate the real specifications of the basic system. WT decomposes a non-stationary signal into a given quantity of stationary sub-signals. Then, AI approaches can be combined with WT to improve preciseness of the prediction. Precise hybrid models have been employed in recent years to forecast hydrological and hydrogeological processes (e.g.,Partal 2009;Adamowski & Sun 2010;Tiwari & Chatterjee 2010).

In most former AI-based hydrological modeling techniques, temporal variables (as historical time series) of the process were used as inputs of the model to predict runoff (Anmalaet al. 2000). Furthermore, in some studies, the AI-based and conceptual models have been combined to estimate river flow (Tokar & Markus 2000;Ashu & Sanaga 2006;Corzoet al. 2009;Huoet al. 2012). These models which combine conceptual hydrological models and AI-based models usually include two parts: one is a conceptual model for runoff of every sub-catchment and the other is the AI-based model to ensemble the runoff from each catchment.Nouraniet al. (2013) developed a neuro-fuzzy-based geomorphological approach to handle the Eel River watershed runoff. The quality of the results supported the beneficial application of the approach in spatio-temporal modeling of hydrological processes. Multi-station modeling is a conceptual approach whose capability has been confirmed and, in some cases, its imposing geomorphological features caused an improvement in results (Zhang & Govindaraju 2003;Sarangiet al. 2005).Nouraniet al. (2016) proposed a spatio-temporal framework in which SVM used conjugate geostatistics to monitor Ajichay watershed runoff-suspended sediment load process.

The objectives of present study are to: (1) use the classic modeling strategy and wavelet-based ELM (W-ELM) to predict one-step-ahead runoff for all hydrometric stations; (2) propose a WT-based spatio-temporal model to forecast the Ajichay watershed runoff with an ELM approach to improve the model ability; and (3) propose an integrated ELMs model conjugated to WT to forecast runoff and catch peak values of the process. Results of the proposed models were compared by using cross-validation technique for station 3 to find the suitability of the models. Finally, results of ELM networks and ANNs were compared to analyze the ELM's capability.

ELM is a recently developed learning scheme for SLFFN, proposed byHuanget al. (2004). Almost all learning algorithms for SLFFNs require adjustment of parameters that results in dependence between different layers of parameters like weights and biases. Therefore, many iterative tuning steps are required by traditional learning algorithms (Huanget al. 2006,2012) whereas the ELM algorithm avoids slow iterative learning procedure and only requires a one-pass operation to learn SLFFNs. This is mainly due to the fact that there is no bias for the output layer neuron. In brief, to initiate one-pass learning operation, the hidden node network parameters (weights and biases) are randomly generated without any prior knowledge or training procedure. Consequently, the ELM turns into a system of linear equations and the unknown weights between the hidden and output layer nodes can be obtained analytically by only applying Moore–Penrose generalized inverse procedure. The learning scheme of ELM can be summarized in three steps (Algorithm 1).

Algorithm 1. Learning scheme of an ELM

Assume

-n inputs,m outputs, hidden nodes (k =1…⁠)

Require

- An activation functiong(x)

• 1:

  Randomly assign parameters of hidden nodes, i.e., weights and bias (w_k,bias_k).

• 2:

  Obtain the hidden layer output matrixH.

• 3:

  Find the output weight matrixβ:β=H^†T, whereH^† represents the Moore–Penrose generalized inverse solution for the hidden layer output matrixH.

With respect to the reviewed literature, most issues about runoff prediction could be related to the AI application. AI-based models which are fed only by temporal features cannot predict/interpolate the variable of interest in an ideal point across the watershed accurately. Such an AI-based model may not be able to cover the natural uncertainty of the hydrological process. Furthermore, the rainfall–runoff data include a broad domain of values and fluctuations; to resolve the mentioned deficiencies, in this research two strategies in rainfall–runoff modeling were proposed. Geomorphology integrated ELM (G-ELM) and integrated ELMs (I-ELMs) conjugated to WT are the strategies of modeling which are going to be discussed. In order to develop the models, m files were coded in the MATLAB environment.

In order to find the nonlinearity and uncertainty of the rainfall–runoff process over the watershed, the G-ELM was employed in this study. For this end, the G-ELM model was established according to the observed rainfall, runoff time series along with spatially varying physical and geomorphological parameters such as upstream drainage area and mean slope of sub-basins. In other words, the model inputs were determined using rainfall–runoff temporal features and geomorphological characteristics of the sub-basins. These physical parameters were selected as inputs of the model since most conceptual models take advantage of such physical/geomorphological parameters to convert independent variables to the variable of interest (e.g.,Melesseet al. 2003). Likewise, these parameters are the preliminary features used by several experimental equations (e.g., SCS and rational equations). Area (upstream drainage area) and mean slope are two geomorphology-based parameters utilized in most physically based rainfall–runoff models. In order to examine the capability of the mentioned parameters in the modeling process, upstream drainage area (A_i) and mean slope (S_i) as spatial variables were selected to determine the input layer of ELM. These features were used as the physical characteristics of the watershed in the proposed model, in addition to hydro-climatologic parameters (i.e., rainfall and runoff temporal features).

In Equation (7),i indicates the sub-basin number.α andβ are the lag times for rainfall (I) and runoff (Q) data, respectively.A_i andS_i are upstream drainage area and slope for sub-basini, respectively. As shown in Equation (7), in addition to temporal variables, spatial variables were considered as inputs of the model. The proposed G-ELM model could be used for monthly rainfall–runoff modeling of the entire watershed (Figure 1). In the G-ELM model, the input variables consisted of different sets of antecedent and current rainfall and runoff values and physical characteristics of the sub-basins relevant to all stations to estimate the runoff time series of the stations (Q_i(t),i = 1, 2, 3, …,n, whereQ is the runoff values at station i,n is the number of sub-basins). As shown in 1, the G-ELM model exhibits a singular model which is capable enough to be used instead of several models (such as AI) for stations inside the watershed. Furthermore, the proposed G-ELM model is able to predict runoff values in stations located in any desired point within the watershed (e.g., station 3 inFigure 1).

In most previous studies, a simple AI approach was used to predict the river runoff. In other words, some antecedent-based or wavelet-based inputs were determined first in which rainfall–runoff data or other relevant features were provided as input data, then by using the input data, river runoff could be predicted. Due to their simple structure, they could not be applied to predict runoff at watershed outlets (Huoet al. 2012). To overcome this weak point, this study planned a conceptual cascade reservoir model based on integrated form of the ELM models. The I-ELMs strategy used herein utilizes the concept of representing watershed as a cascade of reservoirs with different runoff capacity (Nash & Sutcliffe 1970). Similarly, forn reservoirs placed as a cascade, representing each of the sub-basins, the outflow fromnth reservoir could be determined based on the channel drainage network (Figure 2). The integrated form of ELMs (I-ELMs) is developed to take advantage of the semi-distributed models considering the heterogeneity and spatial variation of the rainfall in the watershed. Flow charts of the I-ELMs are illustrated inFigure 2. According toFigure 2, the I-ELMs were designed to employ a meta-learning system; meta-learning is a privileged technique to unify the results of multiple learners (Liao & Feng 2014). I-ELMs model were achieved by forming several base ELMs and one meta-learner. The meta-learner learns from the outputs of the base ELMs, as shown inFigure 2. Each base ELM generates a base-approximator and the meta-learner generates a meta-approximator. In other words, each base-learner is a trained and tested form of ELM for each station located at locationi. According toFigure 2, the approximated results from runoff data were aggregated and were imposed into meta-level training along with decomposed rainfall data. Finally the meta-learner forecasts the runoff for station 13 (approximately located at downstream of the watershed) by enjoying all temporal features which were imposed into the system. As indicated inFigure 2, the entire watershed was divided into sub-basins based on the characteristics of the watershed so that several ELMs could be trained to estimate runoff from the sub-basins, respectively. The developed I-ELMs were formed in a way that the runoff outputs of all base ELMs in the first stage could be aggregated and used as the inputs of meta-learner for forecasting runoff downstream according to the conceptual reservoirs cascade model (Figure 2). The meta-learner is a two-fold training system. In the first stage (base-learners) the ELM was used to train the 13 stations located at upstream of station 13. Next, in the meta-learning stage, the trained models are aggregated to form a unit time series of runoff and then train the meta-learner.

This research was extended to investigate the nonlinear relationship between rainfall–runoff processes in monthly scale for Ajichay watershed which is located in a semi-arid area of Iran. The monthly rainfall and runoff dataset used in this study (Table 1) were provided by East Azerbaijan regional water company. Ajichay watershed's hydrometric stations and rain gauges are located in northwest Iran in Azerbaijan province (between 47° 45′ and 45° 30′ east longitude and 38° 30′ and 37° 45′ north latitude).Figure 3 demonstrates the geographical location and digital elevation map (DEM) of Ajichay watershed. The watershed drainage area is about 10,853 km² and covers about 25% of Urmia Lake. Watershed elevation varies between 1,228 m and 3,755 meters above sea level. River runoff is shed into Urmia Lake, and is momentous for survival of the lake. Mean daily temperature varies from 2.5°C in winter to 20°C in summer. Relative humidity varies from about 55 to 60% over the year. Wind speed in the western area is higher than in the eastern and is 17.5 and 11.5 km/h, respectively. Relative information about rainfall–runoff and watershed spatial information, such as statistical parameters of the data, are given inTable 1. For the purpose of forecasting, which includes the calibration and verification steps, the dataset was divided into two parts. The first 75% of total data was used in the training set and the remaining 25% of data was used for testing purposes. Runoff of station 3, which covers 12 years of rainfall–runoff data, was predicted by using cross-validation technique. Data duration is from 1993 to 2014 (except for station 3 which is from 2002 to 2014). Monthly streamflow time series from 1993 to 2014 for station 13 and summation of other stations were used in this study.Figure 4 illustrates the sum of runoff values for all stations except station 13 and respective linear trends. It is deduced from the figure that all stations have a behavior similar to station 13.

Input and output variables are normalized by scaling between 0 and 1 to eliminate their dimensions and to ensure that all variables receive equal attention during calibration of the model. Therefore,RMSE values are to be transformed to between 0 and 1 (Nouraniet al. 2016).

First, a SLFNN model was trained to predict one-step-ahead runoff with normalized rainfall and runoff features for all stations. ELM is a fast learning emerging technique with high generalization performance. It basically uses single-hidden layer (feature mapping) feed forward neural networks (SLFNNs), with number of hidden nodes that can be randomly generated; in the traditional FFNNs, all the parameters are required to be tuned, which generates a case of dependency between different layers of parameters (weights and biases) (Huanget al. 2006). The significance of this kind of ELM is that the hidden layer does not need to be tuned; learning speed is faster than FFNNs, in addition to better generalization performance.

In order to determine the best classic ELM structure, temporal features were used in different combinations. For this end, lag-1 autocorrelation coefficients (ACF) of rainfall–runoff was performed for all stations. It was observed that generally all stations had satisfactory correlation of rainfall–runoff values with the prior 1, 2, and 12 months (t− 1,t − 2, andt− 12). Hence, the following input combinations (1–5) which included different numbers of input variablesQ andR were considered in determining the input matrix to predict the watershed runoff.

• Comb. 1.Q_(t)

• Comb. 2.Q_(t),R_(t)

• Comb. 3.Q_(t),Q_(t−1),R_(t),R_(t−1)

• Comb. 4.Q_(t),Q_(t−12),R_(t),R_(t−12)

• Comb. 5.db(4,3)Q_(t),db(4,3)R_(t)

In all cases,t demonstrates previous time step,Q represents runoff features, andR represents rainfall features. The output consists of a sole variable, i.e.,Q at future time steps (Q_(t+1)).

In order to determine the ELM outcome, the performance evaluation results are presented for different combinations inFigure 5 (DC and R). According to the results of one-step-ahead fromFigure 5, although Combs. 3 and 4 might have shown adequate efficiency, Comb. 5 showed the best performance. Also, Combs. 1 and 2 could not predict the runoff accurately. From this point of view, it could be concluded that Comb. 5 was the best choice as input combination. Incorporating WT to ELM caused a great increase in prediction accuracy compared to the outcome of Comb. 1. As an illustration,Figure 6 shows the scatterplot of predictions against observed values for stations 2 and 10. It could be inferred that Comb. 5 (W-ELM) had the best outcome.

Performance evaluation was performed to determine the efficient proposed G-ELM architecture. For this purpose, the input combinations (1–6) which include different numbers of input variables were considered in the input matrix (Table 2).Table 2 shows the G-ELM performance for different combinations of input matrix. As shown inTable 2, in order to find out the effect of physical and geomorphological parameters, Comb. 1 was defined with only temporal variables (i.e., only rainfall and runoff data). According toTable 2, the impact of imposing dimensionless physical parameters along with the temporal variables was eligible. In general, Combs. 2–5 showed higher efficiency with respect to Comb. 1. Furthermore, substituting spatially varying parameters (i.e.,A andS) instead of rainfall or runoff variables (i.e.,I(t − 1)/Q(t − 2)) in the input combination led to better performance. By comparing the obtained results for Combs. 1–6 (Table 2), it could be concluded that Comb. 6 proved to be the best structure of G-ELM. The input structure of this combination consisted of two temporal variables (db(4,3)I(t) anddb(4,3)Q(t)) and two spatial variables (i.e.,A andS). The drainage area (A) is an important physical characteristic of the watershed. Drainage area reflects the volume of water that can be generated from rainfall. Mean slope of watershed (S) is a significant factor for runoff yield from rainfall. By considering these facts, the volume of available water for runoff would be the product of rainfall depth and the drainage area and move along the watershed affected by the slope. As well, flood magnitudes reflect the momentum of the runoff and it is obvious that slope is an important factor in the momentum. Therefore, the importance of the roles of drainage area and slope in the rainfall–runoff process could be reconfirmed by the outcome of the current modeling. Conjugating WT to temporal features in Comb. 6 gave an incredible insight into temporal features. Comb. 6 had a supreme performance in the G-ELM.

Figure 7 presents decomposed times series of runoff values bydb(4,3) for all stations formed by the input matrix. InFigure 7, S refers to the main signal which was decomposed to approximation and detail signals.a3 represents approximation sub-series at level 3;d1,d2, andd3 represents detail sub-series at level 1, 2, and 3 respectively.

Figure 8 shows the computed time series versus the observed time series for stations 1–15 (except station 3) in the calibration and validation steps. InFigure 8, St represents an abbreviation for the word station. Each time series was separated into calibration (189 months) and verification (65 months) periods. The consecutive form of multiple-stations time series of Ajichay watershed and location of stations is demonstrated inFigure 8. The figure proves the ability of the proposed G-ELM model for predicting runoff at different parts of the watershed. The reason for this ability could be due to the detection of interaction among stations when their data were imposed onto a unique framework. Additionally, the spatial relationship among the stations could be detected as well as their temporal pattern. As a result, the model should be classified as a time–space model. This model can operate as a nonlinear time–space regression approach rather than a fully lumped model. However,Figure 8 also shows that the proposed model could not catch the peak values. Therefore, the potential risk of flooding could be very high when using the model to predict flow. According to the presented results inTable 2, the G-ELM led to sufficient performance. The main reason for such a suitable efficiency was probably related to the input structure of the G-ELM. Since the runoff variables (i.e.,db(4,3)Q(t)) were considered in the input structure, the autoregressive and Markovian property and seasonality of runoff time series could be caught accurately. It is worth noting that implication of geomorphological parameters, in addition to the temporal variables, caused an improvement in the model performance (Table 2).

The temporal data used in I-ELMs were first pre-processed by WT. Such a pre-processing can cause an improvement in the performance of I-ELMs. Next, each base-ELM was used to forecast the relative target runoff (seeFigure 3). Outcome results of these models were aggregated to form a single time series. In the next stage, the decomposed time series of runoff (aggregated time series) and decomposed rainfall features of station 13 were imposed into meta-learner to forecast station 13's runoff time series. Such a modeling strategy has the privilege of enjoying all datasets so that the meta-learner can handle station 13's time series' specification. One of the great differences between the G-ELM and I-ELMs is the use of spatial characteristics of the watershed as input data. I-ELMs did not employ spatial data, whereas G-ELM was designed based on employing spatial data. In other words, I-ELMs did not apply spatial data of the watershed; in exchange, it was designed in a way that it can explore spatial variation in the rainfall–runoff process and did not require any physical features of the watershed. In order to train the I-ELM, each base-learner was trained usingdb(4,3)Q_i(t),db(4,3)R_i(t) as input data. In other words, decomposed rainfall and runoff values viadb(4,3) were used to calibrate the base models. Next,db(4,3)∑Q_i(t),db(4,3)R₁₃(t) was used as input data to train the meta-learner. Since the meta-learner is trained by using the base-learners, the outcome of the base-learners were aggregated to create a united time series of runoff (∑Q_i(t)). The runoff time series was then decomposed viadb(4,3) and was used to train the I-ELM along withR₁₃(t) (rainfall time series of station 13). Results demonstrated that in the calibration periodDC = 0.95,R = 0.95, andRMSE = 0.028 and for the verification periodDC = 0.88,R = 0.91, andRMSE = 0.42 was obtained. According to the obtained results, I-ELMs performed well and in comparison to G-ELM, outperformed it relatively. However, both these strategies are unique and applicable. It means that these strategies could be set for different spatial and temporal conditions.

As previously mentioned, station 3 was selected to be predicted through cross-validation technique. Cross-validation is a way of qualifying the performance of proposed strategies. For cross-validating the conceptual models, temporal features of station 3 were omitted from the calibration process and it was predicted in the verification process. In order to apply this technique in this study, the time series of station 3 was imposed into the trained system of the models. The G-ELM, which utilized spatio-temporal features in its structure, was set in a way to interpolate station 3's runoff accurately. For this purpose, calibration dataset was arranged according to the G-ELM structure; the only difference was the input dataset, in which station 6's spatio-temporal features were entered into the input matrix because of its approximately similar topographical position. The relation and interaction among the stations’ data, learnt through the calibration phase and collaborated with the model to have appropriate predictions for the cross-validation procedure (Nouraniet al. 2013). Moreover, for verifying the capability of I-ELMs in the cross-validation procedure for station 3, I-ELMs were set as proposed. The only diversity was at meta-learning formation, which was set for 12 years of monthly time steps (calibration and verification, from December 2002 to October 2014). The results of predictions were compared with the available observed data to evaluate performance of the developed models for spatial predictions. In order to evaluate the proposed G-ELM model and I-ELMs model, the predicted values are shown inFigure 9. According to the results, a reliable performance was obtained (DC = 0.8,R = 0.83, andRMSE = 0.1 for G-ELM andDC = 0.78,R = 0.82, andRMSE = 0.12 for I-ELMs) by considering limitations in the modeling procedure. It is obvious that if the G-ELM and I-ELMs were trained using more stations, the model accuracy would be increased in the validation and cross-validation steps. According to the obtained results, G-ELM performed better than I-ELMs. It can be stated that G-ELM, which employed spatial data and similar temporal features, led to better performance. The results of cross-validation in terms of computed versus observed time series and scattered plot are shown inFigure 9, which demonstrates the ability of the G-ELM model for spatiotemporal prediction of the rainfall–runoff process of the Ajichay watershed.

Generally, the ELM algorithm shows faster learning speed over traditional algorithms for SLFNs. For example, it learns faster than the support vector machine by a factor of up to thousands. ELM does not suffer from problems like imprecise learning rate, presence of local minima, etc. This underlines the suitability of ELM for real complex applications that need fast prediction and response capabilities. In addition, the ELM algorithm requires less human involvement, because it does not have any control parameters to be manually tuned, which shows its better applicability for real world applications (Bhatet al. 2008). Finally, recent works confirm the advantages of ELM over earlier approaches for ANN (Malathiet al. 2011). Nevertheless, three key aspects should be pointed out:

• Random initialization of hidden node parameters may affect the performance of ELM (Bhatet al. 2008) that also require high complexity of SLFFNs for improved performance. This may lead to unfavorable conditions, which means that an ELM may not be robust enough to capture variations in data.

• Considering the network complexity issue, it is essential to carefully choose hidden neuron activation functions that show good ability to handle complexity, improve convergence of algorithm, and result in a compact structure of network (Huang & Chen 2008). These aspects are taken into account in the proposed ELM-based model.

• ELMs have been widely applied to hydrology to predict runoff (Abdullahet al. 2015;Limaet al. 2016). One limitation of ELMs is that the prediction can only be made based on historical rainfall–runoff relations. Therefore ELMs will show weakness when trying to predict some high peaks that have not yet been covered in the training data. In this study, pre-processing by WT recovered some weaknesses. In proposed models (i.e., W-ELM and I-ELM) peak values were predicted accurately. However, in G-ELM in which a wide range of fluctuation in data was observed, the model could not catch the peak values.

In order to compare the performance of the presented models, the obtained results via classic models (W-ELM), G-ELM and I-ELMs for station 13 are shown inFigure 10 (as mutual outcome). A zoom window is presented for the first 75 months because of intense fluctuations in these months in comparison with the whole. Also, calculated data were plotted against observed data in order to judge the obtained results more carefully. As indicated fromFigure 10, the G-ELM showed less efficiency with respect to the I-ELMs and W-ELM. W-ELM outperformed proposed strategies in this research by catching maximum and minimum values with better precision. As a general conclusion, it must be stated that ELM is a powerful and fast forecasting approach. For example, by considering W-ELM in which ELM benefits from WT, the modeling procedure possesses a strong insight which makes it supreme in predicting shorter times series (here 252 months for one station). Therefore, such a compound of approaches for short-term modeling performed better than other proposed models. The G-ELM model, which took advantage of spatial characteristics of watershed, was designed in a way to forecast imposed data for all stations (14 stations). Adding WT also decreased error in modeling, but the very wide range of fluctuation caused some defection in some levels of forecasting. The I-ELMs strategy, which benefited from temporal features by taking advantage of WT also could not perform as well as W-ELM. This is because this model did not use the target data in the modeling procedure but took advantage of other station features in two stages including rainfall and runoff data.

Generally, the results of different wavelet-based models were satisfactory. This is due to the capability of WT and nonlinear ability of ELM and mostly the proposed strategies which were designed in a way to catch the peak and least values in modeling procedure.

FFNN is a powerful model in prediction of hydrological variables and was employed as a benchmark to analyze the performance of ELM for station 13. Since Comb. 5 in classic modeling led to the best outcome among all input data, it was used as input for FFNN. According toTable 3, FFNN with eight input data, one hidden layer and one output (8-1-1) used 200 epochs to achieve the best structure. Both wavelet-based FFNN (W-NN) and W-ELM forecast the runoff values of station 13 accurately (seeTable 3). However, results indicated that W-ELM outperformed W-NN by means ofDC andRMSE. The FFNN's structure was optimized via trial and error process, which means FFNN was optimized by several times running in MATLAB environment.

Numerous recently developed rainfall–runoff studies have focused on AI-based modeling in which AI approaches were used along with WT. The present study took advantage of different rainfall–runoff modeling strategies to improve their performances. Calibrating and designating sufficient inputs and the method of employing them are the key to taking advantage of them in a way that they can operate as semi-distributed models or physical-based models but with precision. Among different AI-based models, the ELM model has been used by researchers very recently because of its unique susceptibility. To this end, in this study different strategies of modeling were employed and examined (i.e., geomorphology integrated model and integrated form of ELMs) which were designed according to spatial and temporal features. At first, the lumped form of ELM was developed by using input combinations for one-step-ahead predictions. Comb. 5, which was conjugated to WT and utilized rainfall–runoff features (W-ELM viadb4), led to the best outcome for all hydrometric stations. Subsequently, two approaches were proposed, namely, G-ELM and I-ELMs. G-ELM presented a multi-station modeling of the Ajichy watershed runoff in continuous form which used the geomorphological properties of sub-basins along with decomposed rainfall and run off time series. The presented results for the G-ELM model indicated that the proposed integrated geomorphological approach is an efficient way to predict runoff in multi-station form for the watershed. A conceptual cascade of the reservoir model was developed according to the integrated form of ELM. This strategy did not employ any physical features of the sub-basins. However, the structure of I-ELMs was formed in a way to be able to explore spatial variation of rainfall–runoff process. I-ELMs took advantage of all correlated temporal information and performed better than G-ELM in terms ofDC,R, andRMSE. W-ELM, I-ELM, and G-ELM performance was ranked, respectively, based on statistical evaluations. The cross-validation technique was employed to predict the runoff values for station 3. Based on the results from cross-validation, G-ELM predicted the runoff values of station 3 better than I-ELMs. This could be due to the structure of the proposed models. The I-ELMs were designed to forecast runoff values at the outlet of the watershed; the spatial pattern and relationship among the stations could be detected better through the G-ELM model.

Finally, the W-ELM and conventional W-NN-based modeling was performed for station 13. Results indicate that W-ELM was more suitable than the W-NN model, which could not cope as well with the nonlinear characteristics of the rainfall–runoff process as W-ELM.

The presented models were tuned according to the physical properties and time series of the case study. The presented methodologies can be adopted for other watersheds in order to take advantage of the watershed geomorphological information for AI-based modeling. It is also suggested for future studies to focus on utilizing the spatial and temporal clustering approach. Evidently, long-term data for the watershed will be needed to achieve this goal. Based on the results, W-ELM outperformed I-ELMs and G-ELM; however, it must be considered that each strategy has weak points and advantages, therefore the unique structure of these strategies must be utilized for appropriate situations.

The authors would like to thank East Azerbaijan regional water company for data preparation. The authors also would like to express their gratitude to the anonymous reviewers for their constructive criticism.