Predicting River Flow Using an AI-Based Sequential Adaptive Neuro-Fuzzy Inference System

: Artiﬁcial intelligence (AI) techniques have been successfully adopted in predictive modeling to capture the nonlinearity of natural systems. The high seasonal variability of rivers in cold weather regions poses a challenge to river ﬂow forecasting, which tends to be complex and data demanding. This study proposes a novel technique to forecast ﬂows that use a single-input sequential adaptive neuro-fuzzy inference system (ANFIS) along the Athabasca River in Alberta, Canada. After estimating the optimal lead time between four hydrometric stations, gauging data measured near the source were used to predict river ﬂow near the mouth, over approximately 1000 km. The performance of this technique was compared to nonsequential and multi-input ANFISs, which use gauging data measured at each of the four hydrometric stations. The results show that a sequential ANFIS can accurately predict river ﬂow ( r 2 = 0.99, Nash–Sutcli ﬀ e = 0.98) with a longer lead time (6 days) by using a single input, compared to nonsequential and multi-input ANFIS (2 days). This method provides accurate predictions over large distances, allowing for ﬂow forecasts over longer periods of time. Therefore, governmental agencies and community planners could utilize this technique to improve ﬂood prevention and planning, operations, maintenance, and the administration of water resource systems. This study adopts the adaptive neuro-fuzzy inference system (ANFIS) with di ﬀ erent methods and compares the results with the existing literature to understand if it is possible to employ a data-limited modeling approach that can accurately forecast daily water ﬂow at Fort McMurray. The ANFIS has been largely used in the literature for streamﬂow forecasting worldwide, primarily in mild and temperate areas. Anusree and Varghese [16] compared the performance of the ANFIS, artiﬁcial neural networks (ANN) and multiple nonlinear regression (MNLR) for predicting daily ﬂow at the outlet of Karuvannur River Basin, India. The results showed that the ANFIS model predicts daily ﬂow more accurately compared to the ANN and MNLR models. Sabzi et al. [17] investigated how data preprocessing and data mining techniques can improve the accuracy of streamﬂow predictive models, such as autoregressive integrated moving average (ARIMA), ANN, a hybrid model of ANN and ARIMA (ANN–ARIMA), and the ANFIS. The authors concluded that the ANFIS model achieved a superior streamﬂow prediction performance overall. Dariane and Azimi [18] successfully combined two ANFIS methods: subtractive (sub)-ANFIS and fuzzy C-means (FCM)-ANFIS to forecast streamﬂow in two sub-basins of the Urmia Lake Basin, which is located within two Azerbaijan provinces in northwest Iran. Poul et al. [19] adopted multi-linear regression (MLR), ANN, the ANFIS, and k-nearest neighbors (KNN) to predict the monthly ﬂow in the St. Clair River between the US and Canada. The authors demonstrated that the performances of three nonlinear models of ANN, the ANFIS, and KNN were highly satisfying and that among them, the ANFIS model was superior. Ehteram et al. [20] used the ANFIS to predict the Aidoughmoush monthly streamﬂow in Iran. Their results demonstrated the high capability of the ANFIS in capturing the variability in streamﬂow based on di ﬀ erent climatic indices inputs. The literature demonstrated that the ANFIS generally would perform


Introduction
The modeling of large watersheds is challenging because of the complexity of hydroclimatic regimes due to intra-and inter-basin variations in topography, climatic patterns, land cover, basin drainage density, soil drainage capacity, and other associated factors. For example, simulated flows along the mainstream of rivers located in cold weather regions are usually more sensitive to climate data inputs, while in other cases, where the runoff cycle is interflow-dominated, the hydrologic response is more sensitive to the regional topography [1].
The Athabasca River Basin (ARB) has been subject to several hydrological studies over the past decade because of the increasing population and industrial/agricultural activities that this region has been experiencing over the past 40 years. There is particular interest in understanding the variability in the Athabasca River flow, because it represents an important resource for oil and gas extraction and operational processes, as well as agricultural irrigation. Changes in the magnitude of river flow and seasonality may lead to decreases in water supply, which will impact natural ecosystems, including This study adopts the adaptive neuro-fuzzy inference system (ANFIS) with different methods and compares the results with the existing literature to understand if it is possible to employ a data-limited modeling approach that can accurately forecast daily water flow at Fort McMurray. The ANFIS has been largely used in the literature for streamflow forecasting worldwide, primarily in mild and temperate areas. Anusree and Varghese [16] compared the performance of the ANFIS, artificial neural networks (ANN) and multiple nonlinear regression (MNLR) for predicting daily flow at the outlet of Karuvannur River Basin, India. The results showed that the ANFIS model predicts daily flow more accurately compared to the ANN and MNLR models. Sabzi et al. [17] investigated how data preprocessing and data mining techniques can improve the accuracy of streamflow predictive models, such as autoregressive integrated moving average (ARIMA), ANN, a hybrid model of ANN and ARIMA (ANN-ARIMA), and the ANFIS. The authors concluded that the ANFIS model achieved a superior streamflow prediction performance overall. Dariane and Azimi [18] successfully combined two ANFIS methods: subtractive (sub)-ANFIS and fuzzy C-means (FCM)-ANFIS to forecast streamflow in two sub-basins of the Urmia Lake Basin, which is located within two Azerbaijan provinces in northwest Iran. Poul et al. [19] adopted multi-linear regression (MLR), ANN, the ANFIS, and k-nearest neighbors (KNN) to predict the monthly flow in the St. Clair River between the US and Canada. The authors demonstrated that the performances of three nonlinear models of ANN, the ANFIS, and KNN were highly satisfying and that among them, the ANFIS model was superior. Ehteram et al. [20] used the ANFIS to predict the Aidoughmoush monthly streamflow in Iran. Their results demonstrated the high capability of the ANFIS in capturing the variability in streamflow based on different climatic indices inputs. The literature demonstrated that the ANFIS generally would perform more accurately than ANN for river flow forecasting. In fact, the ANFIS can overcome the disadvantages of ANN models, such as the disregard for data-related uncertainty, which leads ANN models to correlate inputs to outputs using a strict if-then set of rules. At the same time, ANN models are very efficient in adapting and learning. By using the learning capability of ANN and introducing ambiguity in the data inputs by fuzzification, the ANFIS can automatically generate fuzzy if-then rules and optimize its parameters from mathematical algorithms. More details regarding the ANFIS is provided in the Materials and Method section.
Three different methods were adopted in this study using the ANFIS: "Nonsequential", "Sequential", and "Multi-input". The "Nonsequential ANFIS" uses flow data inputs from one station upstream to predict river flow at the station of interest located downstream. The "Sequential ANFIS" uses gauging data collected near the source to sequentially predict flow at different stations downstream. The "Multi-input ANFIS" simultaneously uses multiple gauged flow data located upstream to predict flow at the downstream station of interest.
Existing hydrological models for the ARB usually require a large amount of data in the form of explanatory variables for calibration. Moreover, physical hydrological models are often expensive and necessitate expert personnel in order to properly function. This study proposes a novel application of the ANFIS for streamflow forecasting in cold weather regions using a data-limited modeling approach that can accurately forecast daily water flow over an extended area.

Study Area and Data Source
The Athabasca River Basin is approximately 159,000 km 2 and it represents about 24% of Alberta's landmass. The Athabasca River is the second largest river in Alberta and its average flows are 2.79 × 10 9 m 3 at Jasper, 1.36 × 10 10 m 3 at Athabasca, and 2.09 × 10 10 m 3 at Fort McMurray, per year. The river originates at the Columbia Glacier in Jasper National Park, flowing northeast across Alberta for over 1300 km into Lake Athabasca (Figure 1). The upper reaches of the Athabasca River are located within a mountainous topography characterized by alpine, sub-alpine, and montane ecoregions. This area is historically significant as a waterway for First Nations and the fur trade, as well as the mapping of western Canada. For this reason, the portion of the Athabasca River located within Jasper National Park has been designated a Canadian Heritage River. Industrial developments such as forestry, open pit coal mines, limestone quarries, and growing agricultural areas are located in the middle portion of the Athabasca River Basin. The lower reaches of the Athabasca River begin at Fort McMurray and finish with the confluence of the Peace and Athabasca rivers with Lake Athabasca, forming a vast wetland called the Peace-Athabasca delta. This is known as one of the world's most ecologically significant wetlands and has been designated as a Ramsar Convention wetland and a United Nations Education, Scientific and Cultural Organization (UNESCO) World Heritage Site [21,22]. The lower portion of the Athabasca river basin has undergone an extensive urban and industrial development over the past 40 years due to the extraction of energy resources, primarily oil and gas. Here, surface water assessment is crucial to understand what impact this development is having on the area, because the oil and gas industry relies on the water uptake from the Athabasca River for operational purposes. In addition, the growing energy sector results in specific land uses that influence surface water quality and, subsequently, affect settlements and a variety of people who live along the river.
Water 2020, 12, x FOR PEER REVIEW 4 of 18 wetland called the Peace-Athabasca delta. This is known as one of the world's most ecologically significant wetlands and has been designated as a Ramsar Convention wetland and a United Nations Education, Scientific and Cultural Organization (UNESCO) World Heritage Site [21,22]. The lower portion of the Athabasca river basin has undergone an extensive urban and industrial development over the past 40 years due to the extraction of energy resources, primarily oil and gas. Here, surface water assessment is crucial to understand what impact this development is having on the area, because the oil and gas industry relies on the water uptake from the Athabasca River for operational purposes. In addition, the growing energy sector results in specific land uses that influence surface water quality and, subsequently, affect settlements and a variety of people who live along the river. Generally, the Athabasca River flow is influenced by the large variations in climatic conditions over the year, with long, cold winters and short, warm summers. Near the source of the river in Jasper, the months with the lowest average high temperature are December and January (−6 °C) while the warmest month is July (21 °C ). The average precipitation is highest in July (69 mm) and lowest in April (29 mm). In Fort McMurray, closer to the mouth of the river, the month with the lowest average high temperature is January (−12.2 °C ) while the warmest month is July (23.7 °C ). The average precipitation is highest in July (80.7 mm) and lowest in January (0.4 mm). In cold regions, climatic conditions dictate a river's water sources: there is no contribution of precipitation and snowmelt during the winter, while an abundant rainfall-runoff and snowmelt occur during spring and summer [23]. The large annual variability of water systems in cold weather regions represents a challenge in hydrological modeling. Thus, a data-driven modeling technique that can capture such variability, and bypasses the need to model the complex underlying hydrologic processes governing the flow at Fort McMurray, is selected.
Fort McMurray is the largest urbanized center in the Regional Municipality of Wood Buffalo. This area draws attention from around the world as the residential and commercial focal point of Canada's oil sands industry. The Regional Municipality of Wood Buffalo counts 111,687 people over 66,361 km 2 , where approximately 82,724 people live in Fort McMurray [24,25]. This area is of global significance, as it represents the third largest oil deposit in the world. Although its significance has Generally, the Athabasca River flow is influenced by the large variations in climatic conditions over the year, with long, cold winters and short, warm summers. Near the source of the river in Jasper, the months with the lowest average high temperature are December and January (−6 • C) while the warmest month is July (21 • C). The average precipitation is highest in July (69 mm) and lowest in April (29 mm). In Fort McMurray, closer to the mouth of the river, the month with the lowest average high temperature is January (−12.2 • C) while the warmest month is July (23.7 • C). The average precipitation is highest in July (80.7 mm) and lowest in January (0.4 mm). In cold regions, climatic conditions dictate a river's water sources: there is no contribution of precipitation and snowmelt during the winter, while an abundant rainfall-runoff and snowmelt occur during spring and summer [23]. The large annual variability of water systems in cold weather regions represents a challenge in hydrological modeling. Thus, a data-driven modeling technique that can capture such variability, and bypasses the need to model the complex underlying hydrologic processes governing the flow at Fort McMurray, is selected.
Fort McMurray is the largest urbanized center in the Regional Municipality of Wood Buffalo. This area draws attention from around the world as the residential and commercial focal point of Canada's oil sands industry. The Regional Municipality of Wood Buffalo counts 111,687 people over 66,361 km 2 , where approximately 82,724 people live in Fort McMurray [24,25]. This area is of global significance, as it represents the third largest oil deposit in the world. Although its significance has been recognized for decades, the economic and technological conditions necessary for commercial production have only been recently developed. Moreover, the strong demand for oil and gas, the population and economic growth around the community of Fort McMurray, and the tension between industrial development and environmental protection have attracted attention. Oil sands development requires large amounts of water and energy; the current surface water intake is two to five barrels of water to produce one barrel of oil by mining. The industrial processes used, and the large scale of oil sands development, can result in negative impacts on the aquatic environment if deliberate action is not taken to protect these ecosystems [26][27][28].
The historical average daily flow data from 1971 to 2014 were downloaded from the Water Survey of Canada (WSC) at four stations: Jasper (07AA002), Hinton (07AD002), Athabasca (07BE001), and below Fort McMurray (07DA001) [29]. These locations were selected based on the data consistency and completeness. Table 1 provides general information regarding the four gauging stations used in this study to calibrate and validate the models. To forecast flows at Fort McMurray, antecedent flows at Jasper, Hinton, and Athabasca, which are located along the Athabasca River, were used as independent variables. Two sets of calibration-validation data were selected to forecast river flow at Fort McMurray: (1) data between 1971 and 2000 were used to calibrate the models, while data from 2001 to 2014 were used for model validation and (2) data from odd years (i.e., 1971,1973,1975, . . . , 2013) were used to calibrate the models, while data from even years (i.e., 1972, 1974, 1976, . . . , 2014) were used for model validation. The results of these two different approaches should help to detect possible bias in the calibration data.

Methods
This study focuses on developing a hydrological model to forecast river flow using different methods to the ANFIS, a method that has been successfully used in hydrologic modeling because of its high capability in representing nonlinear natural systems [16][17][18][19][20][30][31][32][33]. The disadvantages of the ANFIS often include a large amount of input data, a long computational time and memory, and mathematical complexity. This study aims to simplify the use of the ANFIS for hydrological modeling, while maintaining a similar level of accuracy. Figure 2 provides a conceptual diagram to show the methods considered in this study. After gathering data and pre-processing into calibrating and validating sets, the flow between stations was compared and correlated to identify the optimal lead time (n), which indicates the amount of time (in days) that is necessary for water to pass from a station upstream to another station downstream. Once the optimal lead time between stations was estimated, three different modeling methods were developed using the ANFIS: "Nonsequential", "Sequential", and "Multi-input". Finally, the models were validated and their performance compared. Among the three types of models, the Sequential ANFIS represents a novel approach. The details and assumptions of these models are described in the following subsections.

Identification of Optimal Lead Time
A correlation analysis was carried out between the flow at Fort McMurray at time t and the flow at other gauging stations (i.e., Jasper, Hinton, and Athabasca) in order to determine the optimal lead time for each pair of stations. Time lags between 1 and 10 days (i.e., t-1, t-2, . . . , t-10) were considered in this analysis. The selection of 10 days was based on the different regimes within the catchment that primarily depend on regional rainfall, topography, and land use. As result, the runoff generated downstream would be influenced by upstream catchments. Considering more than 10 days would not provide much contribution on the catchment of interest. The highest correlation coefficient was considered as the optimal lead time between Fort McMurray and the other stations. Historical daily flow records at Jasper, Hinton, Athabasca, and Fort McMurray between 1971 and 2000 were used to determine the optimal lead time for the first approach, and between 1971 and 2014, only considering odd years (1971,1973,1975, . . . , 2013), for the second approach. By estimating the optimal lead time, it is possible to understand how far into the future a model can predict.
Water 2020, 12, x FOR PEER REVIEW 6 of 18 primarily depend on regional rainfall, topography, and land use. As result, the runoff generated downstream would be influenced by upstream catchments. Considering more than 10 days would not provide much contribution on the catchment of interest. The highest correlation coefficient was considered as the optimal lead time between Fort McMurray and the other stations. Historical daily flow records at Jasper, Hinton, Athabasca, and Fort McMurray between 1971 and 2000 were used to determine the optimal lead time for the first approach, and between 1971 and 2014, only considering odd years (1971, 1973, 1975, …, 2013), for the second approach. By estimating the optimal lead time, it is possible to understand how far into the future a model can predict.

Figure 2.
Conceptual diagram of the approaches and methods adopted in this study.

Models Calibration Using the ANFIS
Suparta and Alhasa [34] and Jang [35] thoroughly described the mechanisms and the mathematics underlying the ANFIS and how this technique is well-suited for highly nonlinear systems. This study considers the ANFIS using grid partitioning (supervised learning algorithm) and adopting the Takagi-Sugeno type inference system. A hybrid algorithm, which is a combination of a least squares estimator and the gradient descent method, is adopted. This means that, during the model training process, a forward (from Layer 1 to Layer 5) and backward (from Layer 5 to Layer 1) propagation algorithm ( Figure 3) adjusts the parameters of the membership functions. The gradient descent method is used to find the nonlinear function minimum, resulting from the weights generated by the fuzzy rules.
Layer 1: For each input variable, there is a set of membership functions that contain function parameters. Each node generates an output that is a degree of membership value given by the input of the membership functions. In this study, membership functions are set as Gaussian distributions because it requires the least number of parameters for calibration compared to other membership function types and the smoothness of the curve allows for a more homogeneous trend in the validation phase.

Models Calibration Using the ANFIS
Suparta and Alhasa [34] and Jang [35] thoroughly described the mechanisms and the mathematics underlying the ANFIS and how this technique is well-suited for highly nonlinear systems. This study considers the ANFIS using grid partitioning (supervised learning algorithm) and adopting the Takagi-Sugeno type inference system. A hybrid algorithm, which is a combination of a least squares estimator and the gradient descent method, is adopted. This means that, during the model training process, a forward (from Layer 1 to Layer 5) and backward (from Layer 5 to Layer 1) propagation algorithm ( Figure 3) adjusts the parameters of the membership functions. The gradient descent method is used to find the nonlinear function minimum, resulting from the weights generated by the fuzzy rules.
Layer 1: For each input variable, there is a set of membership functions that contain function parameters. Each node generates an output that is a degree of membership value given by the input of the membership functions. In this study, membership functions are set as Gaussian distributions because it requires the least number of parameters for calibration compared to other membership function types and the smoothness of the curve allows for a more homogeneous trend in the validation phase.
where µ is the degree of membership functions for the given fuzzy set, x is one of the input variables, and a and b are the parameters of a membership function. Layer 2: Every node is fixed (non-adaptive), and the circle node is labeled as Π. The output node results from the multiplication of incoming signals and is delivered to the next node. The T-norm operator with general performance (AND) is applied to obtain the output, because all the explanatory variables occur simultaneously. w where w j is the output that represents the firing strength of each rule, j represents each node in this layer, and f V i indicates the various forms of membership functions. Layer 3: Every node is fixed (non-adaptive), and the circle node is labeled as N. Each node represents the ratio between the j-th rule firing strength and the sum of all firing strengths. It is also called the normalized firing strength.
Layer 4: Every node is an adaptive node to an output, with a node function defined as follows: where w j is the normalized firing strength from the previous layer and p j x + q j y + r j is a parameter in the node. The parameters in this layer are referred to as consequent parameters. Layer 5: The single node is a fixed (non-adaptive) node that computes the overall output as the sum of all the incoming signals from the previous node. This circle node is labeled as .
The first and the fourth layer contain the parameters that can be modified over time until the gradient descent converges to a minimum error. The first layer contains a nonlinear set of premise parameters, while the fourth layer includes linear consequent parameters. To update both parameter types, a learning algorithm is necessary so that they can adapt to the model's environment. A hybrid algorithm is used in this study. The hybrid learning algorithm consists of two parts: the forward and the backward propagation. The premise parameters in the Gaussian function (a and b) must be steady in Layer 1. A recursive least square estimator (RLSE) method is applied to repair the consequent parameters in the fourth layer. After the consequent parameters are computed, the backward propagation allows for comparison between the generated output and the observed output through the adaptive network input of initial data. The error obtained during the comparison between the generated and actual output is propagated back to the first layer. At the same time, the premise parameters in Layer 1 are updated. One level of hybrid learning (one forward and one backward propagation) is called epoch. Using a hybrid learning algorithm, which combines RLSE and the gradient descent methods, the convergence can be reached faster than using the backward propagation algorithm only, because the dimensional search space of the error is reduced. More details regarding the ANFIS can be found in [36][37][38][39]. Figure 3 shows an example of the fuzzy reasoning mechanisms for this study using the "Multi-input" approach, which employs three input variables simultaneously for forecasting river flow. The "Nonsequential" and "Sequential" models only use one input variable, which simplifies the software's computation exercise.
Three different methods were adopted in this study using the ANFIS: "Nonsequential", "Sequential", and "Multi-input". The "Nonsequential ANFIS" predicts river flow at Fort McMurray using flow inputs from one single station upstream (i.e., Fort McMurray-Jasper, Fort McMurray-Hinton, and Fort McMurray-Athabasca, in three different sets). The "Sequential ANFIS" uses gauging data at Jasper, near the source, to predict flow at Hinton. Subsequently, the forecasted flow at Hinton is automatically entered to predict flow at Athabasca, which in turn is used to predict flow at Fort McMurray. The "Multi-input ANFIS" uses measured data at Jasper, Hinton, and Athabasca simultaneously to predict flow at Fort McMurray. Table 2 provides an overview of the ANFIS settings used per each model. Although the ANFIS codes could be replicated with any programming language, MATLAB ® was adopted in this study with running times between 2 and 10 s using a regular laptop.
Water 2020, 12, x FOR PEER REVIEW 8 of 18 at Fort McMurray. The "Multi-input ANFIS" uses measured data at Jasper, Hinton, and Athabasca simultaneously to predict flow at Fort McMurray. Table 2 provides an overview of the ANFIS settings used per each model. Although the ANFIS codes could be replicated with any programming language, MATLAB ® was adopted in this study with running times between 2 and 10 s using a regular laptop. Figure 3. Example of the fuzzy reasoning mechanisms for this study using the "Multi-input" approach. V1, V2, and V3 are input variables i.e., river flow at Jasper, Hinton, and Athabasca). A, B, and C are the membership functions for each input variable. Π represents the firing strength of the fuzzy logic rules and N is the ratio between the i-th rule firing strength and the sum of all firing strengths. ∑ is the sum of all the incoming signals from the previous node.

Model Validation
The developed models were compared in terms of performance by using quantitative statistical metrics, including the coefficient of determination (r 2 ), the root mean square error (RMSE), and the Nash-Sutcliffe coefficient of efficiency (ENS). The estimated statistics were also used to compare Figure 3. Example of the fuzzy reasoning mechanisms for this study using the "Multi-input" approach. V 1 , V 2 , and V 3 are input variables i.e., river flow at Jasper, Hinton, and Athabasca). A, B, and C are the membership functions for each input variable. Π represents the firing strength of the fuzzy logic rules and N is the ratio between the i-th rule firing strength and the sum of all firing strengths. is the sum of all the incoming signals from the previous node.

Model Validation
The developed models were compared in terms of performance by using quantitative statistical metrics, including the coefficient of determination (r 2 ), the root mean square error (RMSE), and the Nash-Sutcliffe coefficient of efficiency (E NS ). The estimated statistics were also used to compare model performance to the existing literature. The r 2 indicates the goodness-of-fit between measured and predicted flow, while the RMSE is the normalized error represented by the distance between the predicted and the measured flow at Fort McMurray. E NS is a widely used statistic for assessing specifically the goodness of fit of hydrologic models. The quantitative statistical metrics are calculated as follows: where, Y is the predicted flow; Y is the mean of the predicted flows; X is the observed antecedent flow; X is the mean of the observed antecedent flows; n is the number of observations. Note that the r 2 ranges between 0 and 1, where 1 indicates a perfect fit between the observed and predicted values. The E NS can range between −infinity (−∞) and 1, where 1 corresponds to a perfect fit. Finally, RMSE should be close to zero to indicate good model performance and its magnitude can vary between +infinity (+∞) and −infinity (−∞).

Identification of the Optimal Lead Time
A correlation analysis was performed to estimate the time (in days) required for a mass of water to flow from one station to another. This quantity is also called optimal lead time. The coefficients of determination were calculated by comparing the flow at Fort McMurray with the flow at station i (i = Jasper, Hinton, and Athabasca) at different time lags such as t (same day), t-1, t-2, …, t-10. A similar analysis was carried out for Approach 2. Table 3 shows the correlation parameters between each station from t to t-10 for Approach 1 and 2.
Between Jasper and Fort McMurray, r 2 shows a poor correlation, where the highest coefficient of determination is observed at 5 days. Similar results indicate that between Hinton and Fort McMurray, the optimal lead time is 4 days. A strong correlation was found between Athabasca and Fort McMurray, where the highest r 2 was 0.92 with two-day lead time. A similar correlation analysis was carried out between Jasper and Hinton and Hinton and Athabasca to estimate the optimal lead time between each station along the river. The optimal lead time between Jasper and Hinton is 1 day (r 2 = 0.96), while that between Hinton and Athabasca is 3 days (r 2 = 0.63). Figure 5 schematically summarizes the main findings in terms of how far in advance each developed model can predict flow at Fort McMurray. It should be noted that between Jasper and Fort McMurray, the best correlation was found at 5 days, while when summing the optimal lead time between each station (i.e., Jasper-Hinton = 1 day, Hinton-Athabasca = 3 days, and Athabasca-Fort McMurray = 2 days) the total optimal lead time is 6 days. This might be due to the actual optimal lead time between Jasper and Fort McMurray being in between 5 and 6 days (for Jasper-Fort McMurray, r 2 = 0.504 at t-5 and r 2 = 0.502 at t-6).

Identification of the Optimal Lead Time
A correlation analysis was performed to estimate the time (in days) required for a mass of water to flow from one station to another. This quantity is also called optimal lead time. The coefficients of determination were calculated by comparing the flow at Fort McMurray with the flow at station i (i = Jasper, Hinton, and Athabasca) at different time lags such as t (same day), t-1, t-2, . . . , t-10. A similar analysis was carried out for Approach 2. Table 3 shows the correlation parameters between each station from t to t-10 for Approach 1 and 2.
Between Jasper and Fort McMurray, r 2 shows a poor correlation, where the highest coefficient of determination is observed at 5 days. Similar results indicate that between Hinton and Fort McMurray, the optimal lead time is 4 days. A strong correlation was found between Athabasca and Fort McMurray, where the highest r 2 was 0.92 with two-day lead time. A similar correlation analysis was carried out between Jasper and Hinton and Hinton and Athabasca to estimate the optimal lead time between each station along the river. The optimal lead time between Jasper and Hinton is 1 day (r 2 = 0.96), while that between Hinton and Athabasca is 3 days (r 2 = 0.63). Figure 5 schematically summarizes the main findings in terms of how far in advance each developed model can predict flow at Fort McMurray. It should be noted that between Jasper and Fort McMurray, the best correlation was found at 5 days, while when summing the optimal lead time between each station (i.e., Jasper-Hinton = 1 day, Hinton-Athabasca = 3 days, and Athabasca-Fort McMurray = 2 days) the total optimal lead time is 6 days. This might be due to the actual optimal lead time between Jasper and Fort McMurray being in between 5 and 6 days (for Jasper-Fort McMurray, r 2 = 0.504 at t-5 and r 2 = 0.502 at t-6). The correlation analysis performed on the calibration-validation datasets for Approach 2 showed results identical to those for Approach 1. This indicates that the optimal lead time is independent of the calibration dataset.
Water 2020, 12, x; doi: FOR PEER REVIEW www.mdpi.com/journal/water The correlation analysis performed on the calibration-validation datasets for Approach 2 showed results identical to those for Approach 1. This indicates that the optimal lead time is independent of the calibration dataset.

Evaluation of the Results
The performance of the developed models was compared and the results are summarized in Table 4. Overall, "Nonsequential Jasper-Fort McMurray" and "Nonsequential Hinton-Fort McMurray" can poorly predict river flow at Fort McMurray because, although their r 2 values indicate a good correlation between the measured and predicted flow, the RMSE is considerably large in both approaches. A good performance is observed in the "Nonsequential Athabasca-Fort McMurray" model (ENS = 0.99 and RMSE = 49 m 3 /s for Approach 1 and ENS = 0.99, and RMSE = 46 m 3 /s for Approach 2), which can predict river flow at Fort McMurray over 2 days. Although this model is accurate, the limited predictive capability represents a disadvantage. The "Sequential" model helps to cope with this limitation, providing accurate predictions over 6 days (ENS = 0.98 and RMSE = 66 m 3 /s for Approach 1, and ENS = 0.99 and RMSE = 43 m 3 /s for Approach 2). Among the models proposed in this study, the most accurate is "Multi-input", with ENS = 0.98 and RMSE = 39 m 3 /s using Approach 2. However, similarly to "Nonsequential Athabasca-Fort McMurray", this model allows predictions over 2 days. Table 4. Summary of the results to show the models' performance in terms of coefficient of determination (r 2 ), Nash-Sutcliffe efficiency coefficient (ENS), root mean square error (RMSE), and predictive capability.

Evaluation of the Results
The performance of the developed models was compared and the results are summarized in Table 4. Overall, "Nonsequential Jasper-Fort McMurray" and "Nonsequential Hinton-Fort McMurray" can poorly predict river flow at Fort McMurray because, although their r 2 values indicate a good correlation between the measured and predicted flow, the RMSE is considerably large in both approaches. A good performance is observed in the "Nonsequential Athabasca-Fort McMurray" model (E NS = 0.99 and RMSE = 49 m 3 /s for Approach 1 and E NS = 0.99, and RMSE = 46 m 3 /s for Approach 2), which can predict river flow at Fort McMurray over 2 days. Although this model is accurate, the limited predictive capability represents a disadvantage. The "Sequential" model helps to cope with this limitation, providing accurate predictions over 6 days (E NS = 0.98 and RMSE = 66 m 3 /s for Approach 1, and E NS = 0.99 and RMSE = 43 m 3 /s for Approach 2). Among the models proposed in this study, the most accurate is "Multi-input", with E NS = 0.98 and RMSE = 39 m 3 /s using Approach 2. However, similarly to "Nonsequential Athabasca-Fort McMurray", this model allows predictions over 2 days. Table 4. Summary of the results to show the models' performance in terms of coefficient of determination (r 2 ), Nash-Sutcliffe efficiency coefficient (E NS ), root mean square error (RMSE), and predictive capability.  Figure 6 shows the hydrographs comparing the measured to the predicted flow at Fort McMurray using Approach 1. The scatter plots are shown in Figure 7. In general, the "Sequential" and "Multi-input" models perform more accurately during the springtime increase and late summer decrease than the "Nonsequential" models. The advantage of using the "Sequential" model is the higher predictive capability, as indicated in Figure 8. The "Nonsequential Athabasca-Fort McMurray" model is more accurate in predicting the base flow in the colder months, when there is no contribution of rainfall and snowmelt, while the "Multi-input" model could not perform as accurately. The scatter plots in Figure 7 show the low accuracy of the "Nonsequential" method using the Jasper and Hinton stations to predict river flow at Fort McMurray.
Water 2020, 12, x FOR PEER REVIEW 13 of 18 Figure 6 shows the hydrographs comparing the measured to the predicted flow at Fort McMurray using Approach 1. The scatter plots are shown in Figure 7. In general, the "Sequential" and "Multi-input" models perform more accurately during the springtime increase and late summer decrease than the "Nonsequential" models. The advantage of using the "Sequential" model is the higher predictive capability, as indicated in Figure 8. The "Nonsequential Athabasca-Fort McMurray" model is more accurate in predicting the base flow in the colder months, when there is no contribution of rainfall and snowmelt, while the "Multi-input" model could not perform as accurately. The scatter plots in Figure 7 show the low accuracy of the "Nonsequential" method using the Jasper and Hinton stations to predict river flow at Fort McMurray.  The results for Approach 2 are shown in Figure 8 in the form of annual hydrographs, and in Figure 9 as scatter plots. The predictive performance of the ANFIS improved in all methods (i.e., "Nonsequential", "Sequential", and "Multi-input") adopted in this study when using a calibration dataset that is time independent. The "Multi-input" method was able to perform accurately during Water 2020, 12, x FOR PEER REVIEW 13 of 18 Figure 6 shows the hydrographs comparing the measured to the predicted flow at Fort McMurray using Approach 1. The scatter plots are shown in Figure 7. In general, the "Sequential" and "Multi-input" models perform more accurately during the springtime increase and late summer decrease than the "Nonsequential" models. The advantage of using the "Sequential" model is the higher predictive capability, as indicated in Figure 8. The "Nonsequential Athabasca-Fort McMurray" model is more accurate in predicting the base flow in the colder months, when there is no contribution of rainfall and snowmelt, while the "Multi-input" model could not perform as accurately. The scatter plots in Figure 7 show the low accuracy of the "Nonsequential" method using the Jasper and Hinton stations to predict river flow at Fort McMurray.  The results for Approach 2 are shown in Figure 8 in the form of annual hydrographs, and in Figure 9 as scatter plots. The predictive performance of the ANFIS improved in all methods (i.e., "Nonsequential", "Sequential", and "Multi-input") adopted in this study when using a calibration dataset that is time independent. The "Multi-input" method was able to perform accurately during The results for Approach 2 are shown in Figure 8 in the form of annual hydrographs, and in Figure 9 as scatter plots. The predictive performance of the ANFIS improved in all methods (i.e., "Nonsequential", "Sequential", and "Multi-input") adopted in this study when using a calibration dataset that is time independent. The "Multi-input" method was able to perform accurately during the colder months, which represents an improvement from Approach 1, shown in Figure 6. In addition, the "Nonsequential" models were capable of better predicting the late summer decrease when compared to Approach 1. Similar to Approach 1, the scatter plots in Figure 9 show the lower accuracy of the "Nonsequential" method when the Jasper and Hinton stations were used to predict river flow at Fort McMurray. Table 5 provides the results to show the inter-annual variations in terms of r 2 , E NS , and RMSE.
Water 2020, 12, x FOR PEER REVIEW 14 of 18 the colder months, which represents an improvement from Approach 1, shown in Figure 6. In addition, the "Nonsequential" models were capable of better predicting the late summer decrease when compared to Approach 1. Similar to Approach 1, the scatter plots in Figure 9 show the lower accuracy of the "Nonsequential" method when the Jasper and Hinton stations were used to predict river flow at Fort McMurray. Table 5 provides the results to show the inter-annual variations in terms of r 2 , ENS, and RMSE.  Generally, the ANFIS is superior to other modeling techniques reported in the literature to predict the Athabasca River flow at Fort McMurray. All the three methods used in this study, namely the "Nonsequential Athabasca-Fort McMurray", "Sequential", and "Multi-Input" ANFISs, performed better than other modeling techniques previously explored by other authors on the ARB,  the colder months, which represents an improvement from Approach 1, shown in Figure 6. In addition, the "Nonsequential" models were capable of better predicting the late summer decrease when compared to Approach 1. Similar to Approach 1, the scatter plots in Figure 9 show the lower accuracy of the "Nonsequential" method when the Jasper and Hinton stations were used to predict river flow at Fort McMurray. Table 5 provides the results to show the inter-annual variations in terms of r 2 , ENS, and RMSE.  Generally, the ANFIS is superior to other modeling techniques reported in the literature to predict the Athabasca River flow at Fort McMurray. All the three methods used in this study, namely the "Nonsequential Athabasca-Fort McMurray", "Sequential", and "Multi-Input" ANFISs, performed better than other modeling techniques previously explored by other authors on the ARB, Generally, the ANFIS is superior to other modeling techniques reported in the literature to predict the Athabasca River flow at Fort McMurray. All the three methods used in this study, namely the "Nonsequential Athabasca-Fort McMurray", "Sequential", and "Multi-Input" ANFISs, performed better than other modeling techniques previously explored by other authors on the ARB, using both Approach 1 and 2. This indicates that the ANFIS is highly capable of capturing the nonlinearity of the natural river cycles over the year in cold weather regions, while bypassing the physical explanation of the input-output variables' dependence. Not only can the ANFIS predict more accurately, it also uses a simpler set of input-output variables compared to the more complex dataset used by the VIC or SWAT, which employ a large amount of data for climate and runoff information to calibrate the model. Previous studies carried out in the ARB, which were discussed in the introduction section, show a lower performance in terms of Nash-Sutcliffe coefficient, compared to the method proposed in this paper. By using SWAT, Shresta et al. [14] achieved the highest accuracy level among past attempts found in the literature for the ARB, despite their highly data-demanding set of explanatory variables. Eum et al. [1] also showed a highly valuable modeling performance using a combination of VIC and GCM, although this approach could lead to longer computational times and more expensive budgets. The method proposed in this study shows that it is possible to achieve a higher accuracy when a limited number of inputs are employed, and a more simplistic input-output relationship is outlined. Interestingly, the use of the two calibration-validation dataset pairs (Approach 1 and 2) led to a difference in performance using the ANFIS for the three methods adopted in this study. This was also indicated by Zheng et al. [40], who investigated the statistical behavior of data splitting methods to achieve representative evaluation performance for flow forecasting [40][41][42]. Other modeling techniques should be investigated in a similar fashion to better understand the contribution of calibration-validation datasets on the accuracy of the model output. Finally, Approach 2 provided more accurate results for the three ANFIS methods, possibly because the model's outcome is not influenced by time-dependent variables. For example, the variation in rainfall and snowmelt contributions, earlier springtime increases or late summer decreases, and the growing water uptake from the oil and gas industry from more recent years should be investigated and correlated to the river flow variations overtime.

Conclusions
This study used the adaptive neuro-fuzzy inference system (ANFIS), which is an artificial intelligence (AI) technique for machine learning, to forecast river flow at Fort McMurray, located on the lower reaches of the Athabasca River in Alberta, Canada.
Different techniques using the ANFIS were developed and compared to the existing literature. Initially, a correlation analysis was carried out between the flow at Fort McMurray and the flow at other gauging stations at different times. The highest correlation coefficient indicated the optimal lead time between Fort McMurray and the stations upstream. Three distinct techniques were then adopted: "Nonsequential", "Sequential", and "Multi-input". Although the "Nonsequential" and "Multi-input" models were capable of accurately predicting river flow at Fort McMurray (r 2 = 0.99, E NS > 0.98), they only allowed predictions with a two-day notice, while the "Sequential ANFIS" could forecast accurate flow regimes and allowed modeling with a six-day notice. Subsequently, a different set of calibration and validation data were adopted to perform the same analyses and compare the accuracy of the results. The latter approach provided more accurate results for the three ANFIS methods, possibly because the model's outcome was not influenced by time-dependent variables (i.e., variation in rainfall and snowmelt contributions, earlier springtime increases or late summer decreases, and the growing water uptake from the oil and gas industry from more recent years).
In conclusion, the "Sequential" ANFIS modeling technique is recommended to forecast daily river flow at Fort McMurray because of its capability in capturing the nonlinearity of the natural river cycles over the year in cold weather regions, while bypassing the physical relationship of the input-output variables. This study thus demonstrates the successful application of the ANFIS for sequential river flow forecasting in cold weather over an extended geographical area. The simplistic approach and the lower computational resources and time required for this exercise could find a use of this model in assisting governmental agencies and communities to improve flood prevention and the planning of water resource systems, operations, maintenance, and administration.