River Water Salinity Prediction Using Hybrid Machine Learning Models

: Electrical conductivity (EC), one of the most widely used indices for water quality assessment, has been applied to predict the salinity of the Babol-Rood River, the greatest source of irrigation water in northern Iran. This study uses two individual—M5 Prime (M5P) and random forest (RF)—and eight novel hybrid algorithms—bagging-M5P, bagging-RF, random subspace (RS)-M5P, RS-RF, random committee (RC)-M5P, RC-RF, additive regression (AR)-M5P, and AR-RF—to predict EC. Thirty-six years of observations collected by the Mazandaran Regional Water Authority were randomly divided into two sets: 70% from the period 1980 to 2008 was used as model-training data and 30% from 2009 to 2016 was used as testing data to validate the models. Several water quality variables—pH, HCO 3 − , Cl − , SO 42 − , Na + , Mg 2 + , Ca 2 + , river discharge (Q), and total dissolved solids (TDS)—were modeling inputs. Using EC and the correlation coe ﬃ cients (CC) of the water quality variables, a set of nine input combinations were established. TDS, the most e ﬀ ective input variable, had the highest EC-CC (r = 0.91), and it was also determined to be the most important input variable among the input combinations. All models were trained and each model’s prediction power was evaluated with the testing data. Several quantitative criteria and visual comparisons were used to evaluate modeling capabilities. Results indicate that, in most cases, hybrid algorithms enhance individual algorithms’ predictive powers. The AR algorithm enhanced both M5P and RF predictions better than bagging, RS, and RC. M5P performed better than RF. Further, AR-M5P outperformed all other algorithms (R 2 = 0.995, RMSE = 8.90 µ s / cm, MAE = 6.20 µ s / cm, NSE = 0.994 and PBIAS = − 0.042). The hybridization of machine learning methods has signiﬁcantly improved model performance to capture maximum salinity values, which is essential in water resource management. two standalone and eight novel hybrid-algorithms to predict river salinity. The algorithms used include M5P, RF, bagging-M5P, bagging-RF, RS-M5P, RS-RF, RC-M5P, RC-RF, AR-M5P, and AR-RF. Monthly data collected over a 36-year period were used to construct several sets of variable combinations as inputs for model building and evaluation. Qualitative (visual) and quantitative criteria were applied to validate and compare models. The results reveal that TDS is the most important variable for predicting EC (r = 0.99). Among 10 input combinations, the ﬁrst, a standalone TDS, was signiﬁcantly important to the results. The combinations demonstrate that not


Introduction
Rivers are the principal sources of water for human consumption, irrigation, municipal, and industrial demands, and provide habitat for aquatic species in many regions of the world [1]. The deterioration of water quality of rivers causes irreparable damage to the environments and human health [2,3]. Monitoring and assessments of water quality parameters (e.g., salinity, dissolved oxygen (DO), algae, nitrogen (N), among others) are necessary to develop water management plans [4], as more than one billion people lack access to safe drinking water. Thus, robust and flexible models are urgently needed to accurately predict water quality and to estimate future supplies.
Numerous physical or mathematical models have been developed to predict and plan for the management of water quality (i.e., QUAL2K, MOUSE), but the models are complex, time-consuming (especially during the calibration phase), and data-intensive. These models are challenging for users in developing countries where data are insufficient or where background information is scant. Statistical models of water quality have been developed based on both linear and non-linear relationships between input and output variables. These models, however, often fail to adequately represent the complexity of these phenomena in environments with multivariate non-linear relationships [5,6]. There are non-linear, stochastic, and lagging relationships among several water quality parameters, and it is challenging to create a mathematical model to predict events in these circumstances with traditional approaches [7].
Artificial intelligence (AI) algorithms for water quality modeling have been explored in recent years. These algorithms explore hidden and complex relationships among input and output variables to craft models that best represent these relationships. Several advantages of AI models over traditional physical and statistical models include: the data that are needed for the AI models can be collected relatively easily, sometimes from remote sensing platforms; AI models are less sensitive than traditional approaches to missing data; the structures of AI models are flexible, non-linear, and robust; and AI models can handle huge amounts of data and data at different scales [8][9][10]. Research over the last couple of decades has explored water quality modeling using numerous machine learning (ML) algorithms [11][12][13][14][15][16][17][18][19][20][21][22][23][24].
ANN, the most widely used ML algorithm, has poor prediction power, especially when validation data are not in the range of the training data, and when using small datasets [38][39][40][41][42]. Thus, a hybrid of ANN and fuzzy logic (FL) was proposed in an adaptive neuro-fuzzy inference system (ANFIS) to overcome this weakness. Although ANFIS was successful and had higher prediction accuracy than both ANN and FL, it was still unable to properly determine the weight of the membership function. To resolve this, bioinspired (or metaheuristic) algorithms have been used to automatically determine these weights. Some, using metaheuristic algorithms with ANFIS, have found that hybrid algorithms have higher powers of prediction than any standalone ANFIS [43,44].
Decision-tree algorithms are more powerful predictors than those algorithms with hidden layers in their structure (i.e., ANN, ANFIS, and support vector machine (SVM)) [8]. Barzegar et al. [45] used extreme learning machine (ELM) and wavelet-extreme learning machine hybrid (WA-ELM) to predict EC in the Aji-Chay River watershed and compared the prediction powers of those models to ANFIS. They found that WA-ELM outperformed the other models. ELM is popular and well-known for fast training, but has the disadvantage of being unable to encode more than one layer of abstraction.
Currently, data-mining algorithms like RF, logistic model tree (LMT), and naïve Bayes trees (NBT) are widely used in flood-susceptibility mapping [46,47], groundwater vulnerability assessments using the bootstrap aggregating (BA) algorithm [48], and landslide-susceptibility mapping with Bayesian algorithms [49,50]. It is clear from that these algorithms have been widely used for spatial modeling, and their accuracies have been vetted with field data, but they are seldom used for prediction with time-series data. They have been used to predict suspended sediment load with a reduced error pruning tree (REPT), M5P, and instance-based learning (IBK), and random committee-REPT (RC-REPT) [9], to predict apparent shear stress with RF, M5 Prime (M5P), and random tree (RT) algorithms [51] to predict concentrations of dissolved oxygen in rivers with M5P and SVM [52] and to estimate solar radiation with RF, RT, M5P, and REPT [53].
This paper develops, explores, and compares the effectiveness of two standalone (M5P and RF) and eight novel hybrid algorithms (bagging-M5P, bagging-RF, random subspace (RS)-M5P, RS-RF, random committee (RC)-M5P, RC-RF, additive regression (AR)-M5P, and AR-RF) to predict salinity in the Babol-Rood River in northern Iran. The models' performances are evaluated using model efficiency indices and graphical comparison. These hybrid algorithms have not previously been used to predict surface water quality (salinity). In fact, there have been no prior attempts to apply them for any geoscientific purpose. Thus, the main contribution of this study is that it develops and evaluates these hybrid algorithms for the first time and offers new insights into the applicability of hybrid ML algorithms. Not only will their use enhance the capacity to predict salinity in rivers, but also in other studies of hydrology and hydraulics.

Study Area
The Babol-Rood River is in the Savadkouh district of Mazandaran Province, Iran ( Figure 1) in the Alborz Mountains near the Caspian Sea. The elevation in the catchment ranges from 55.36 to 3317.38 m. The average elevation is 923.55 m. Rangelands are located between 1800 and 3100 m. Ninety percent of the area is between 2200 and 2900 m. Forestlands are found between 150 and 2500 m. The slopes of the watershed have been categorized into five classes. The largest proportion of the slopes are from 30 to 60%, whereas slopes less than 5% are least common. A large part of the region's geology, mostly in the central part of the watershed near the North Alborz fault, is Mesozoic bedrock and sediments are from the Cenozoic. The soils are dictated by geomorphology, climate, and vegetation. The soils of steep slopes with sparse vegetation cover, particularly in the watershed's upperparts, are semi-deep with light to medium texture. In contrast, soils of the middle and lower parts, where slopes are gentler and where vegetation is moderately dense to dense, are moderately deep to deep and have medium to heavy textures. The soils of the area are classified as B and C in terms of permeability. The catchment is characterized by high precipitation of up to 2000 mm annually with fall and winter maxima. Summers are hot, and winters are mild. In the coastal plains and in the highest mountain areas, the climate is Mediterranean with dry summers. The Ghorantalar hydrometric station at the catchment outlet provides observation records. Most river salinity in northern Iran, and especially in the study area, is due to the seawater intrusion, saline water discharge of agricultural regions, and the carbonate formations [54].

Data Collection and Preparation
The Babol-Rood River's water salinity was calculated as EC (µs/cm) and was modeled using ten water quality variables: water discharge (Q), total dissolved solids (TDS), pH, HCO 3 − , Cl − , SO 4 2− , Ca 2+ , Mg 2+ , and Na + . A 36-year (1980 to 2016) monthly record of water quality was obtained from the Mazandaran Regional Water Authority (MRWA). The data were randomly divided into two sets-70% and 30%, a ratio widely used for both spatial [46,47], temporal, and time-series modeling [9,49]. The first 28 years of data  were used for training to build the models, and the remaining eight years (2009-2016) were used for testing that validated the models (Table 1). EC values ranged from 0.03 to 900 µs/cm. To remove the impacts of the different scales of the variables and to guarantee prediction flexibility [38], the data values were normalized (x i ) to a range from 0 to 1 [55][56][57][58][59]: where, x i , x min , and x max are measured data, the minimum value of variable x i , and the maximum value of the of variable x i , respectively. The chemistry of rivers in the study region is mainly affected by effluents and nutrient loads from agricultural runoff and seawater intrusion. Water temperature affects solution concentrations in shallow non-flowing waters where evaporation can enhance salinity and push it to high levels. Therefore, modeling of water bodies, like lakes, reservoirs, or marshes, ought to factor temperature. Flowing rivers, by contrast, are more affected by agricultural contaminants, sea level rise, and saltwater intrusion in coastal areas than by temperature, and thus it is not an important factor.

Selection of Input Combinations
There are two main steps to devise appropriate input combinations: determination of the most effective input variables and identification of the optimum values for each model's operator. First, the proposed algorithms were tested using several input combinations to determine the most effective input variable or combination of variables. In total, nine input combinations were devised according to the correlation coefficient (CC) of EC with each water quality variable ( Table 2). The most effective variable (i.e., variable with the highest CC), in this case TDS, was introduced into each model. The working hypothesis is that the variable with the highest CC suffices to accurately predict EC. The variable with the next highest CC was added to the first and each subsequent variable was added in a stepwise fashion to the growing string until the variable with lowest CC was added (i.e., a combination of all 9) ( Table 3). To identify the best input combination, the individual and hybrid models were performed using a fixed set of parameters (operator) [8]. All models were evaluated using the 10 input combinations, and the root mean square error (RMSE) determined the most effective combination in the test set. Several other model forms are possible, considering all possible input combinations (2 9 − 1 = 511). Hence, we limited our candidate set of models to the input combination in accordance with previous studies [53,60].

Identification of Optimum Values of Operators
After determining the best input variable combination, optimum values for each model's operators were determined by trial and error [8]. There are no optimum global values for each operator [61][62][63][64]. Therefore, different values were tested and the optimum value was selected according to the RMSE of the testing phase as in the previous stage. In the first iteration, the models were run using the default values. Then according to the result of each iteration, values above and below the default values were examined randomly. This continued until the optimum values for each model's operator were obtained.

M5P Algorithm
Quinlan [65] introduced the M5 algorithm as a type of decision-tree. Mosavi et al., [66] expanded the original M5 algorithm to create the M5P. A remarkable capability of the M5P algorithm is that it can handle large datasets efficiently, as well as many variables and many dimensions. M5P is also able to handle missing data [67]. An example depicts the M5P model predicting electrical conductivity ( Figure 2). To process data and build a decision tree, the M5P algorithm performs 4 steps: Step 1: Divide the initial input space into multiple subspaces. Separation criteria are used to minimize the intra-subspace overlap (from root to the node). Standard deviation (sd) is used to measure the variability of the intra-subspace. The standard deviation reduction (SDR) criterion is used to construct the tree: where D is the dataset that reached the node and D i denotes the sub-space (splitting node).
Step 2: Develop a linear regression model for each sub-space using the sub-dataset.
Step 3: Prune the nodes of the developed tree. If SDR is lower than the expected error while developing the linear model (for sub-tree), over-fitting can occur. Pruning is undertaken to overcome this.
Step 4: Smooth the discontinuities caused by the pruning process. A final tree is created by combining the leaf and root to filter the predicted value. Combining the linear model and the node from leaf to root uses: where P' indicates the predicted value passed to higher nodes; p passed denotes the predicted values passed to the current node from the previous nodes; p model is the predicted values of the model at the current nodes; n train is the number of observations in the training dataset that reach the previous nodes; σ is a constant. The M5P model has a simple structure and also a few operators to set, thus hybridization can improve the prediction power of the RF standalone model.

Random Forest (RF)
RF is a decision-tree algorithm introduced by Breiman [68] enhanced and used in numerous studies [69][70][71][72][73][74][75]. RF is flexible, able to deal with classification and continuity issues (i.e., regression problems). It was extended based on bagging [76] and competed with the boosting approach [71]. To make predictions, RF employs bootstrapping to split data into multiple sub-datasets. Subsequently, decision trees are developed for each sub-dataset. Finally, predictions from the sub-decision trees are ensembled, and the final prediction reflects the entire RF model [77,78]. RF has become a favorite ML algorithm for practical engineers because it is able to handle both regression and classification data and learns and predicts very quickly. The performance of the RF model depends on only one or two hyper-parameter(s); it is premised upon generalized error estimation and is able to process data with a large number of dimensions. A pseudo-code of the RF algorithm was developed (Figure 3, adapted from [79]) to use RF for prediction of the water salinity levels of the Babol-Rood River. RF does not give precise continuous predictions in the form of regressions but achieves high performance in classification. Hence, hybridization can also improve the prediction power of the standalone RF model. Small changes in the dataset will lead to instability in ML algorithms. Breiman [76] proposed a technique called "bagging" to overcome model instability and to improve performance. Many studies have successfully employed bagging to increase the effectiveness of predictions [80]. The general concept of the bagging algorithm is: Consider a dataset D = (y EC n , x EC n ), n = 1, 2, . . . , N , in which y EC n is the output (i.e., EC) and x EC n denotes the input variables (i.e., TDS, Na + , Hco 3 − , Cl − , SAR, Ca 2+ , Mg 2+ , So 4 2− PH, Q). For any given ML algorithm or artificial intelligence model, a λ(x EC , D) model (classifier or regressor) is generated. Five forms of ML models are created: M5P, RF, random subspace (RS), random committee (RC), and additive regression (AR). Bagging tries to combine them to improve the accuracies of single models λ x EC , D p from D p , where p = 1, 2, . . . , P is the number of learning sets generated from D. For each p, a small change is made in D.
Bagging can be applied to improve many base learners, like decision trees. For example, bagging was used to train the M5P and RF base learners to predict river salinity.

Random Subspace (RS)
Like bagging, RS uses the powers of bootstrapping and aggregation to build a forecasting model. RS, however, uses the bootstraps in the feature space instead of training samples, as bagging does [81]. Like RF and bagging, RS can handle both regression and classification data. It consists of 4 main components: the training dataset x, the number of subspaces L, the classifier or regressor w, and the number of features d s [82]. In RS, the number of subsets is generated randomly in d s features. Subsequently, they are saved in the subspace L. In the second step, each base regressor/classifier is trained/learned on each of the subsets to create a different regressor/classifier. Then, these are combined to build an ensemble regressor/classifier E. More details of RS can be found in Kuncheva and Plumpton [83] and Kuncheva et al. [84]. Herein, RS is used to construct hybrid models with the M5P and RF for EC prediction.

Random Committee (RC)
A meta-algorithm, RC is an ensemble of random tree classifiers that can improve the learning performance of classification. RC uses a number of random seeds on the same dataset to build a group of random classifications, upon which predictions are made. Finally, averages are computed as outcome predictions based on the predictions generated by each of those classifiers [85]. RC has been employed in many fields. It has been combined with ANN for electrical disturbance classification [86], with the random tree by voting for classifying anomalies [87], and with bagging and J48 algorithms for efficient intrusions classification [88]. In this study, RC was used to predict EC.

Additive Regression (AR)
AR is a nonparametric algorithm proposed by Breiman and Friedman [89]. It is an essential part of the alternating conditional-expectations algorithm. To build a restricted class of nonparametric regression models, AR uses a one-dimensional smoother. It is more flexible than linear models. It is also more easily interpreted than is a general regression model. However, overfitting, model selection, and multicollinearity are the disadvantages of AR. A brief description of AR algorithm is: Consider a dataset D = (y EC n , x EC n ), n = 1, 2, . . . , N , in which y EC n is the output (i.e., EC) and x EC n denotes the input variables (i.e., TDS, Na + , HCO 3 − , Cl − , SAR, Ca 2+ , Mg 2+ , So 4 2− PH, Q). The AR model will take a form: In other words, the AR model can take a form: where f j (x ij ) functions are unknown smooth functions, which are used to fit the dataset. The back-fitting algorithm can be used to fit the AR model in this case [90].

Model Evaluation Criteria
To evaluate the models' performances, five statistical criteria-root mean square error (RMSE), r squared (R 2 ), mean absolute error (MAE), Nash-Sutcliffe efficiency (NSE), and percent bias (PBIAS)-were calculated in their testing phases. These criteria were calculated as follows [91][92][93][94][95][96]: where EC measured , EC measured , EC predicted , and EC predicted denote the measured, mean of measured, predicted, and mean of predicted values of EC, respectively, and n is the number of observations. In addition, the Taylor diagram [97] and box plot [98] were used to provide for visual assessment of the models. The Taylor diagram provides an overview of correlations and standard deviations of each model. The closer the predicted value of CC and SD are to the observed value, the higher is the prediction capacity. The box plot provided an overview of the models' capabilities using extreme values, medians, and the first-and third-quartile predictions.

Best Determinant Combination of EC
Ten water quality variables were selected to be tested as potential predictors of EC ( Table 2). The analyses of correlations between EC and the water quality variables show a strong linear relationship (R = 0.91) between EC and TDS. Therefore, because TDS was the most effective input variable, it was included in all ten input combinations. In this step, the analysis was limited and compared only with the RMSE statistic calculated in the testing phase ( Table 4). The results reveal that the best input combination was constructed with TDS alone (i.e., combination 1), which is consistent with the high CC of TDS to EC. In other words, including other variables in the model only distorts the models' abilities. This result, however, is case specific and studies of other variables may find that other combinations are more predictive. For instance, Khosravi et al. [9] predicted suspended sediment load with data mining algorithms using best-input combinations with eight, nine, and ten variables. Barzegar et al. [32,41] used Ca 2+ , Mg 2+ , Na + , So 4 2− , and Cl − as best-input variables to predict water salinity in the Aji-Chay River, Iran. Combining other explanatory variables with TDS consistently decreased the models' accuracies, reflected in increasing RMSE values. Alongside the model's structure, data accuracy, data type, and data structure have strong effects on the result. The best input combination, therefore, should be examined for each new dataset and should not rely solely on the variables with the highest correlation coefficients. To sum up, the best input combination is input No 1, in which TDS is the only input to the models.

Models' Performances and Validation
After determining the most effective input variable and the optimum value of each model's operator, every model was trained with the training dataset and was evaluated with the testing dataset [99]. The standalone-and the hybrid-models' results were compared and validated using time-series of EC predictions and observations, scatter plots, box plots, and Taylor diagrams (Figures 4-6). A deeper look into these comparisons was also provided and validated with several quantitative statistical indices (Table 5). Using time-series and scatter plots (Figure 4), the first observation is that there is high agreement between the values predicted by the algorithms and the values measured by the M5P, bagging-M5P, RS-M5P, and AR-M5P algorithms. There was no apparent bias between the measured and calculated values, and the models captured the extreme EC values well. This performance is partially explained by the high correlation coefficient between EC and TDS (R = 0.99, Table 2) and the data-mining techniques' high prediction ability.
There is no visible difference between the predicted EC values of the four best models (M5P, bagging-M5P, RS-M5P, AR-M5P) and the measured values ( Figure 5). Nor can one see differences between the predicted median values, the predicted 25th and 75th percentiles, and even the extreme values. The more flexible algorithms-M5P, bagging-M5P, RS-M5P, AR-M5P-marginally overestimate maximums above the measured 900 µs/cm at 906, 901, 906, and 906 µs/cm, respectively. Furthermore, they slightly underestimate the minimums above the measured 165 µs/cm to 173, 172, 172 and 172 µs/cm, respectively. Our results reveal that the models with highest predictive powers based on RMSE can generally predict extreme values properly, but the results are inconsistent. As Khosravi et al. [100] showed, although hybrids of ANFIS-bioinspired models generally fit the data well, they were less robust in their estimates of the extreme values of reference evaporation. By contrast, random-tree algorithms which have generally poorer overall performance, predict extreme values accurately. Khozani et al. [51] found results similar to Khosravi et al. [100] when predicting apparent shear stress.
It is also clear from the Taylor diagram ( Figure 6) that the models are systematically comparable to the measured values. This is good proof of the models' performance and point to their usefulness for water quality modeling if implemented properly using suitable input variables. Although all of the models have good performance metrics (i.e., CC > 0.99), EC predicted by M5P, bagging-M5P, RS-M5P, and AR-M5P have similar standard deviations, and higher CC with EC, while the rest of the models have smaller standard deviations. To sum up, AR-M5P is the best predictor among the algorithms.     Although all models performed well (NSE > 0.75 and PBIAS < ±10%), the best accuracy was achieved by the hybrid AR-M5P algorithm, which had excellent accuracy (NSE = 0.994, RMSE = 8.9, MAE = 6.20), approximately equal to the success of RS-M5P. These were followed by the bagging-M5P (NSE = 0.992, RMSE = 8.96, MAE = 6.24) and the M5P (NSE = 0.991, RMSE = 9.04, MAE = 6.29). These had significantly better results than the other models. RC-M5P, and AR-RF had comparable results (NSE > 0.97, MAE < 10 µs/cm, RMSE < 10 µs/cm). RF, bagging-RF, RS-RF, and the RC-RF had lower accuracy (RMSE and MAE > 18 µs/cm) and thus are less suitable for modeling EC. The RMSE was the greatest for the RS-RF and RC-RF hybrid models (19.30 µs/cm). The model with the next largest was bagging-RF (19.01 µs/cm) and then the individual RF model (18.50 µs/cm). It might be argued that such large RMSE values calculated using the individual and hybrid RF models leads to clearer conclusions. But comparing the individual models, M5P is a better prediction model than RF; the M5P had at least 51.13% lower RMSE and 38.21% lower MAE than the RF model. The M5P algorithm is a better predictor than hybrid-RF algorithms, as well. M5P has advantages over RF: it can manage large datasets with large numbers of attributes and with numerous dimensions, and it can deal with missing data without any ambiguity [9].

Discussion
Now adays application of machine learning algorithms as a practical tool in a different field of geoscience increased rapidly such as flood forecasting [101][102][103][104][105], relationship between fish community and environmental factor [106], groundwater modeling [48], and many other field of study which declared in introduction section.
Hybrid algorithms might increase the prediction power of individual algorithms (e.g., compare the M5P results with the results of M5P hybrids as well as RF with its hybrids). There is an exception to this: RC reduces prediction quality below the individual M5P and RF models. Therefore, though some hybrid algorithms enhance performance of individual algorithms, this isn't always true. Based on PBIAS, M5P and its hybrids tend to overestimate (except in the case of RC-M5P) and RF and its hybrids (and also RC-M5P) underestimate their predictions.
As the M5P model has a higher prediction power than RF (Table 5), hybrid M5P algorithms are more suitable for predicting EC than either the individual or hybrid RF models. In another words, hybrid algorithms generate results that are dependent upon the results of the individual models themselves. Generally, additive regression (AR) algorithms may increase the prediction power of the individual M5P and RF algorithms more than other hybridizing algorithms. This is a result of their flexibility and structural consistency with the other models. RC causes a lowering of the performance of the individual algorithms.
Since the dataset used to predict EC was the same for all models, the difference in performance is the result of each model's structure (i.e., its flexibility, prediction capacity, and tendency to over-fit). The computational capabilities and the complexity of each ML model dictate its results [107]. Generally, decision tree-based algorithms (M5P and RF) require no assumptions about the data distribution, adapt to outliers by isolating them in small regions of the feature space [108], have no hidden layers in their structure, and use tree algorithms to estimation by pattern recognition, and therefore perform better; particularly in comparison to the data-intelligence algorithms which have hidden layers in their structures [8]. This finding conforms to previous studies that indicated the superiority of tree-based models in terms of estimation capacity [109]. Due to the non-linear nature of many environmental phenomena (e.g., EC), more flexible models with non-linear structures will yield better results [110]. According to the literature, hybrid algorithms in most of the cases are more flexible and can enhance the prediction power of standalone models [9,49,51,111].
Previous studies [112][113][114][115][116], have applied several kind of ML models to estimate river EC. Compared to these studies, our results seem to be more accurate. Rohmer and Brisset [112] linked EC to seawater levels and Q using a kernel-based support vector machine (SVM) technique and reported a CC around 0.90; this is significantly less accurate than the results of our study. Ravansalar and Rajaee [113] applied the multilayer perception neural network (MLPNN) and the MLPNN combined with wavelet decomposition (WMLPNN) to predict EC in the Asi River at the Demirköprü gauging station, Turkey, using Q and the EC measured in the previous period. They demonstrated that wavelet decomposition significantly improves the model's accuracy and the WMLPNN (R 2 = 0.94) significantly outperformed the standard MLPNN (R 2 = 0.38). However, the WMLPNN was less accurate than the AR-M5P (R 2 = 0.99) in our study. Azad et al. [114] compared adaptive neuro-fuzzy inference system (ANFIS), ANFIS optimized with particle swarm optimization (PSO) called ANFIS-PSO, and ANFIS with ant-colony optimization for continuous domains (ACOR) called ANFIS-ACOR while predicting EC using TH, Cl − , and Na + . The results revealed high performances of all three models: R 2 = 0.99 for ANFIS-PSO, R 2 = 0.98 for ANFIS-ACOR, and R 2 = 0.94 for ANFIS. The hybrid algorithms in this study performed just a bit better, but the individual M5P and RF algorithms compared to ANFIS had higher prediction powers and it is in accord with Choubin et al. [109] which reported that decision-tree algorithms of classification and regression trees (CART) outperformed ANFIS, SVM, and ANN in their prediction of suspended sediment load. Tutmez et al. [115] used ANFIS model to predict EC using the TDS and achieved R 2 = 0.937. Al-Mukhtar and Al-Yaseen [116] compared the ANFIS, MLPNN, and the MLR models while predicting EC using several water quality variables and achieved NSE = 0.98. Barzegar et al. [41] used several water quality parameters to forecast EC with ELM, WA-ELM, and ANFIS, and reported NSE ranging from 0.6 to 0.95. And Ghorbani et al. [81] predicted EC with MLPNN and MLR models. The MLPNN was more accurate (R 2 = 0.965, RMSE = 50.810 µs/cm, and MAE = 37.495 µs/cm). In sum, the algorithms used in this study outperformed all of the algorithms used in previous efforts to predict EC.
The results show that these models can be used for river salinity prediction and even forecasting with high accuracy. The outcome is undoubtedly useful for water quality protection and management. It is better to use these types of models for real-time river salinity forecasting which is beneficial and practical.

Conclusions
The Babol-Rood River is the main source of irrigation and drinking water in its watershed. Deterioration of its water supply can cause irreparable damage to the environment and the resident population. Monitoring to predict its future conditions need accurate models. These tasks are the foundations of planning to conserve and manage crucial resources. This study is a pioneering evaluation of the use of two standalone and eight novel hybrid-algorithms to predict river salinity. The algorithms used include M5P, RF, bagging-M5P, bagging-RF, RS-M5P, RS-RF, RC-M5P, RC-RF, AR-M5P, and AR-RF. Monthly data collected over a 36-year period were used to construct several sets of variable combinations as inputs for model building and evaluation. Qualitative (visual) and quantitative criteria were applied to validate and compare models. The results reveal that TDS is the most important variable for predicting EC (r = 0.99). Among 10 input combinations, the first, a standalone TDS, was significantly important to the results. The combinations demonstrate that not only does a model's structure have significant effects on prediction accuracy, but so too do the input variables. According to the validation criteria, M5P has greater prediction power than RF and those hybrid algorithms with M5P performed better than those hybridized with RF. The M5P algorithm and its hybrids (except RC-M5P) can accurately predict EC, even in extreme conditions. The results also reveal that hybrid algorithms enhanced the performances of standalone algorithms, except in the case of the RC model, which reduced prediction power when hybridized with M5P and RF. Although all algorithms showed very good prediction capacities, AR-M5P outperformed all of the others. The rest, in order of accuracy, were RS-M5P, bagging-M5P, M5P, RC-M5P, AR-RF, bagging-RF, RF, RS-RF, and RC-RF. These results cannot be generalized to other study areas or with other hydrological data, but AR-M5P would undoubtedly be one of the more accurate algorithms, if not the most accurate. Hybridized models outperformed the standalone models in this study. Overall, the river water quality can be accurately predicted through such enhanced machine learning algorithms and several readily available input variables.