Groundwater quality analyses from across the entire state of Punjab were retrieved from various published reports and research articles. We collected a total of 17,317 F⁻ observations: 1,705 data points from Central Groundwater Board (CGWB), 433 observations from Central University of Punjab (CUPB), 745 observations from Duggal and Sharma (41), 11,226 from Khattak et al. (8), 59 observations from Sharma et al. (42), 38 observations from British Geological Survey (BGS) (43), and 3,111 from Department of water supply and sanitation, Punjab government. Besides F⁻, the availability of groundwater physicochemical parameters such as pH, Electrical Conductivity (EC), Total Dissolved Solids (TDS), Chloride (Cl⁻), Nitrate (NO₃ ⁻), Sulphate (SO₄ ²⁻), Phosphate (PO₄ ³⁻), Bicarbonate (HCO₃ ⁻), Sodium (Na⁺), Potassium (K⁺), Calcium (Ca²⁺), and Magnesium (Mg²⁺) is essential for model development. In conjunction with the groundwater quality determinants, the depth of the collected samples was also considered an important parameter for this study. These variables were selected based on their established or suspected association with the discharge and accumulation of F⁻ in groundwater and were further used to screen the data.

Although these attributes are often measured during groundwater monitoring assessments, some were missing from the datasets collected from different sources. Data collected from research papers (8, 41, 42) did not contain all of this information and, hence, was not included in the final database. In comparison, almost all variables were present in data collected from BGS, CUPB, and CGWB. Although CGWB data was collected from 2013–2015 and 2018–2020 (44–49), only the recent data of the year 2020 was used for prediction modeling. Also, the observations from locations monitored in the previous years but not in 2020 were considered for prediction modeling. This resulted in a total of 298 observations from CGWB. Furthermore, CUPB data consisted of some samples collected from canals and other surface water sources, which were excluded, resulting in a total of 420 data points. Besides the water quality variables and sampling depth, information on the geographical coordinates of the sampling locations was also considered as an essential screening criterion. For sampling points lacking the georeferenced location, Google Earth Pro was used to determine the same by using the name of the sampling location. Therefore, the final database that was finally used for creating a suitable prediction model for groundwater F⁻ concentration contained a total of 756 data points from CGWB, CUPB, and BGS ( Figure 1A , Table 1 ). The final data was also classified into two classes as per the depth at which the samples were collected (optimum depth was considered to be 60 m) (52, 53).

The values of some of the abovementioned variables were missing in the final dataset, which was estimated with the help of standard formulas. One such missing parameter was TDS in CGWB data that was determined using the following formula (Equations 1, 2) (54):

Similarly, a few other parameters were reported below their respective detection limit (BDL). These values were then replaced by dividing the BDL value by two. Furthermore, the values of all the parameters were converted to their respective similar unit to ensure uniformity in the dataset. All the anions, cations and TDS were represented in mg/L, EC in µS/cm, and depth in meters, while pH is unitless. Of the final 756 F⁻ measurements, 609 (~81%) were under the permissible limit of 1.5 mg/L, 100 (~13%) ranged from 1.5–3 mg/L, and the remaining 47 (~6%) were greater than 3 mg/L.

To understand the nature of the distribution of groundwater F⁻ levels and their corresponding physicochemical parameters at different depths, graphical and statistical inference methods were adopted. The compiled dataset was characterized by enumerating its descriptive statistics (minimum, maximum, mean, median, coefficient of variation, first and third quantiles, and percentage of samples exceeding the respective permissible limits). The normality for all the variables was tested using Kolmogorov–Smirnov test. Testing whether the data is normally distributed is necessary, especially for geochemical and other environmental data, because they are generally skewed, consisting of outliers and originating from varied sources (55). Normality testing further aided in selecting the appropriate statistical treatments for the data. Since most parameters are not normally distributed, Spearman’s rank correlation coefficient was enumerated to identify the potential associations of the F⁻ concentrations with the concurrently evaluated physicochemical attributes and sampling depth. All the statistical analyses and graphical plotting were performed in the R software version 4.3.2.

Machine learning algorithms were further applied to uncover the hidden patterns between the compiled F⁻ concentrations and the physicochemical variables and well depth and develop an optimized model for predicting F⁻ concentration in the study domain. In this study, F⁻ is the output or target (y) variable that will be determined using the input or predictor (x) variables, i.e., the abovementioned physicochemical attributes and the sampling well depth. Five different machine learning models, i.e., Extreme Gradient Boosting (Xgboost), Random Forest (RF), Support Vector Machine (SVM), Extreme Learning Machine (ELM), and Multilayer Perceptron (MLP), were implemented and tested on the final dataset. All of these models have been frequently used in the literature for groundwater-based investigations and, hence, considered for groundwater F⁻ prediction modeling. R software version 4.3.2 was used to develop these proposed models. Before the implementation of these models, a pre-processing step was involved in which data standardization was performed using the Z-score method with the following formula (Equation 3) (56):

where x i ´ denotes the standardized i ^th variable, x i   is the i ^th variable, σ   denotes standard deviation, and μ is the mean. Following standardization, the entire dataset was randomly shuffled, and a cross-validation technique was employed to further split the data for training and testing the model. 80% of the data was used for training the model, and the remaining 20% was used for validation.

Knowledge of the hydrochemical conditions of groundwater is indispensable for identifying potential contaminants to safeguard human health. The descriptive statistics for summarizing the hydrochemical characteristics of the groundwater samples compiled from different sources for this study are presented in Table 2 . The pH value for the entire dataset ranged from 6.0 to 9.1 ( Table 2 ), with median values of 8.05 and 7.33 at depths within and exceeding 60 meters ( Figure 3 ), respectively, indicating the predominance of alkaline conditions in the aquifers of this region. Likewise, EC and TDS ranged from 41–16,760 μS cm⁻¹ and 29–13,408 mg/L ( Table 2 ), respectively, with median values higher in shallow (809 μS cm⁻¹ and 531.6 mg L⁻¹, respectively) than in deeper (601 μS cm⁻¹ and 419.2 mg L⁻¹, respectively) waters ( Figure 3 ). According to Freeze and Cherry’s groundwater classification (78), shallow groundwaters in Punjab can be majorly considered as brackish (1000 < TDS < 10000 mg L⁻¹), while deeper waters are classified as freshwater (TDS <1000 mg L⁻¹). Furthermore, the results also show the occurrence of both cations and anions in excess, particularly in shallow depths. Dominant cations in shallow groundwater include Ca²⁺, Mg²⁺, Na⁺, and K⁺, and anions such as Cl⁻, NO₃ ⁻, SO₄ ²⁻, F⁻, and HCO₃ ⁻, whereas, in deeper waters were Mg²⁺, Na⁺, Cl⁻, and HCO₃ ⁻. Based on overall median concentrations, the cations and anions were arranged in the following order: Na⁺ > Mg²⁺ > Ca²⁺ > K⁺ and HCO₃ ⁻ > Cl⁻ > SO₄ ²⁻ > NO₃ ⁻ > F⁻ > PO₄ ³⁻, respectively. Furthermore, it was also observed that the median concentrations of Cl⁻, NO₃ ⁻, SO₄ ²⁻, Mg²⁺, Na⁺, Ca²⁺, and K⁺ were also elevated in shallow waters compared to deeper waters. The median values of Mg²⁺ and F⁻ were slightly higher in deeper waters than in shallow aquifers (Mg²⁺: 35 and 34.39 mg L⁻¹, respectively; F⁻: 0.71 and 0.57 mg L⁻¹, respectively). Similarly, the median concentration of HCO₃ ⁻ was much higher in deeper waters than in shallow waters. Therefore, it is evident that the majority of the ions, along with other water quality parameters, are in excess in the shallow aquifers than the deeper groundwater ( Figure 3 ).

All the hydrochemical parameters except pH displayed a greater degree of coefficient of variation [CV (%)] ( Table 2 ), as well. This clearly implies a wide range of variability within each of the water quality parameters in the study region, arising from various natural sub-surface and surface phenomena and anthropogenic influences. In addition to this, the physicochemical attributes of most of the groundwater samples, particularly sampled from the shallow aquifers, exceeded the recommended safe limit by the World Health Organization (1) ( Table 2 ). The overall percentage of samples exceeding their respective permissible limit for each parameter is as follows (total % exceeding/% exceeding in shallow groundwater samples): EC: 32%/29%; TDS: 48%/40%; Cl⁻: 13%/12%; NO₃ ⁻: 20%/18%; SO₄ ²⁻: 26%/24%; HCO₃ ⁻: 38%/24%; Na⁺: 20%/18%; K⁺: 19%/17%; Ca²⁺: 3%/2.6%; Mg²⁺: 54%/34%; F⁻: 19%/15%.

Furthermore, the Kolmogorov–Smirnov test verified that all the groundwater quality variables did not follow a normal distribution. In shallow waters, F⁻ had a weak to moderate positive correlation with almost all the variables except Ca²⁺. The Spearman’s rank correlation coefficients of F⁻ with all the variables are: Depth = 0.14 (p<0.01); pH = -0.23 (p<0.01); EC = 0.41 (p<0.01); TDS = 0.44 (p<0.01); Cl⁻ = 0.31 (p<0.01); NO₃ ⁻ = 0.14 (p<0.01); SO₄ ²⁻ = 0.35 (p<0.01); HCO₃ ⁻ = 0.34 (p<0.01); Na⁺ = 0.30 (p<0.01); K⁺ = 0.22 (p<0.01); Mg²⁺ = 0.31 (p<0.01); Ca²⁺ = 0.05 (p>0.01). EC and TDS had the highest influence on F⁻ concentration, thus indicating an increase in its concentration with an increase of these parameters.

The utilization of groundwater physicochemical parameters as predictor variables for forecasting F⁻ contamination levels through ML approaches has been well established (31). In this study, five different models with diverse architectures, such as RF, SVM, Xgboost, ELM, and MLP, were employed for predicting the groundwater fluoride concentration in the aquifers of Punjab. The model performance was evaluated based on the R² , RMSE, and MAE values. These are some of the commonly used metrics for determining the predictive ability of ML models. In the case of the Xgboost model, the Adaptive random search function among the other functions had the highest R² value and the lowest RMSE and MAE values in the testing stage ( Table 3 ). This implies that the Adaptive random search function is the best activation function for Xgboost in the current study, which was further considered for comparing the prediction performance with all the proposed models. For SVM, the radial basis kernel function (RBF) performed better than polynomial, sigmoid, and linear kernel functions as it can handle non-linear datasets (31) and, therefore, selected for prediction purposes in our study. This superior performance of RBF over other kernel functions was also confirmed by Rajasekaran et al. (79), Wu and Wang (80), Amirmojahedi et al. (81). Also, among the various activation functions in the ELM model, the ‘tansig’ function had the most satisfactory output and, therefore, was selected for further prediction-performance comparison between the selected models.

The overall statistical evaluation criteria for all the models yielded poor to satisfactory results, implying that a few models outperformed others in predicting the F⁻ levels. Based on the 80% of the total dataset used for training purpose, the R² achieved for different models are 0.42 (RF), 0.52 (SVM), 0.34 (Xgboost), 0.85 (ELM), and 0.21 (MLP) ( Table 3 ). Ideally, R² close to unity display greater proximity between the observed and simulated values. Although R² provides an indication of how well the model fits the data, with values close to 1 implying a better fit, it does not provide information about the magnitude of the errors between the actual and predicted values. Hence, RMSE and MAE values were computed along with R² to assess the performance of the different ML models. The RMSE values were 0.77 (RF), 0.70 (SVM), 0.84 (Xgboost), 0.46 (ELM), and 0.10 (ML), and MAE was 0.47 (RF), 0.41 (SVM), 0.49 (Xgboost), 0.36 (ELM), and 0.06 (MLP). Both RMSE and MAE values closer to 0 suggest little error between the actual and predicted values. Based on these values in the testing stage, MLP had the least amount of error followed by ELM, SVM, RF, and Xgboost. Despite the lowest RMSE and MAE values, MLP had the lowest R² , suggesting unreliable performance for F⁻ determination in this study. Furthermore, RF also trained very poorly, which is evident from the low R² and significantly greater RMSE and MAE. In addition to MLP and RF, SVM and Xgboost were also trained unsatisfactorily as per the R² , RMSE, and MAE values. On the contrary, relatively lower MAE and RMSE values and greater R² value of the ELM model indicates superior training ability relative to the other four models.

After model training, the remaining 20% of the dataset was utilized for testing the model, and the same evaluation metrics were applied to analyze each model’s predictability. The trend of the performance evaluation criteria for all the models was almost similar. The order of the proposed models in terms of the R² values was ELM (0.95) > Xgboost (0.70) > SVM (0.66) > RF (0.44) > MLP (0.33) ( Table 3 ). Satisfactory R² values were observed for ELM, Xgboost, and SVM. Furthermore, error metrics RMSE and MAE for MLP were 0.10 and 0.06, respectively, which were the least among all the models. However, the lowest R² value for MLP obtained from the trained data implies poor model performance and proved unreliable for groundwater fluoride prediction in this region. After MLP, the RMSE (0.31) and MAE (0.33) values for ELM in the testing phase were minimal among the remaining models, emphasizing good prediction ability. From comparing the statistical performance metrics of the training and testing stages of different models, it is evident that only ELM had the optimum values and can be considered for modeling F⁻ concentrations in Punjab. It is noteworthy that MLP and ELM have relatively less complex topology and training algorithms than the remaining three models (30). Nevertheless, their performance varied greatly in predicting the groundwater fluoride concentration in the study domain.

In order to better comprehend the accuracy of model prediction, the observed F⁻ concentrations and their corresponding predicted values after model training were plotted in a scatter diagram ( Figure 4 ). From Figure 4 , it is quite evident that the distribution of predicted F⁻ values in relation to the observed F⁻ concentrations is quite closely placed to the best fitting line as opposed to other models, which validates the robustness of the ELM model. Besides ELM, the predicted values of SVM, Xgboost, RF, and MLP did not closely match the actual values, which is substantiated by poor R² values ( Figure 4 ).

Spatial distribution maps were prepared to better visualize the actual and predicted F⁻ concentration values for all five models ( Figure 5 ). The predicted values of ELM, Xgboost, SVM, RF, and MLP were compared with the original F⁻ concentrations, and a significant difference between the model outcomes was noted. The south and southwestern regions of Punjab have elevated F⁻ levels (> 1.5 mg/L) in its groundwater system ( Figure 5A ), whereas the remaining areas exhibited relatively lower concentrations. A significantly similar F⁻ distribution pattern of the ELM predicted values ( Figure 5B ) with the original F⁻ concentrations was observed, implying a substantial prediction accuracy. Central Punjab had concentrations ranging between 0.5–1.0 mg/L, with parts of northwestern districts having groundwater fluoride surpassing 1.5 mg/L, identical to the original F⁻ distribution and ELM predicted map. Furthermore, based on the concentration values predicted by the remaining four models, excess F⁻ levels were evident across the entire Punjab state. Although central and northern regions had relatively safe groundwater fluoride levels (< 1.5 mg/L) ( Figure 5A ), contradictory F⁻ distribution as per the predicted values of Xgboost, SVM, RF, and MLP was observed ( Figures 5C–F ), implying poor model performance. Also, the regions falling in the south and southwest had a greater magnitude of F⁻ content than the original F⁻ levels, indicating an overestimation of excess contaminant levels. These spatial distribution maps of the proposed ML models, in comparison to the original F⁻ distribution, are evidently in accordance with the performance evaluation metrics, i.e., R² , RMSE, and MAE ( Table 3 , Figure 4 ). Consequently, it can be stated that ELM outperformed the remaining four models in groundwater fluoride levels in the study region.

The predictor (input) variables govern the robustness and stability of the prediction models (82, 83), and therefore, the relative importance ranking of these variables aid in determining the significant variables or factors influencing the contamination. The ranking of the variables was found to be consistent in ELM and SVM ( Figure 6 ), while other models displayed certain variations. In both ELM and SVM, the variables contributing the most to model prediction were TDS, EC, Cl⁻, Na⁺, and Ca²⁺, each with relative importance greater than 15% ( Figure 6 ). These variables also displayed a significant correlation with F⁻. Based on the relative importance scores, TDS, EC, and Na⁺ were the top three variables in these two models in the order TDS > EC > Na⁺, indicating their potential role in mobilizing and enhancing F⁻ in the groundwater of the study domain, particularly in the shallow aquifers. The variable importance of the abovementioned factors was also observed to be significant in the remaining three models; however, it was not in the same order as ELM and SVM. Based on the relative importance scores, the order of the variables was EC > TDS > Na⁺ in Xgboost, TDS > Na⁺ > EC in RF, whereas EC > Ca²⁺ > Na⁺ > TDS in MLP was observed. The relative importance of Cl⁻ and Ca²⁺ ranked 4^th and 5^th in ELM, SVM, and RF, while Xgboost (Ca²⁺ and Cl⁻ ranking 4^th and 5^th, respectively) and MLP (Ca²⁺ and Cl⁻ ranking 2^nd and 5^th, respectively) displayed slight variation. The variable SO₄ ²⁻ attained the least importance in ELM (0.5), SVM (0.2), Xgboost (1.1), RF (0.3), and MLP (0.3) ( Figure 6 ). The ranking discrepancies between these variables in all the five models could have resulted due to the differences in model algorithms (24).

The occurrence of elevated F⁻ concentrations in the groundwater of the Punjab basin is attributed mainly to geogenic origin (84), which results from the interplay of multiple complex interdependent hydrogeochemical processes (85). Naturally occurring F⁻ in minerals and rocks are nearly insoluble in water. However, favorable conditions facilitate the dissolution of these minerals, further releasing F⁻ into the groundwater (86). Dissolution of fluoride-bearing minerals, in particular fluorite, under suitable conditions, such as alkaline pH with excess EC and TDS, as observed in the aquifers of Punjab, favors F⁻ enrichment. The alkaline nature of the groundwater could be attributed to the presence of sediments containing abundant carbonate minerals (85). Elevated EC and TDS could be a consequence of the rapid and greater degree of rock and mineral weathering, waterlogging, and dissolution of salts (17, 56). Furthermore, anthropogenic inputs from agricultural and industrial practices and groundwater recharge through shallow aquifers also contribute to high EC and TDS (87). An increase in EC and TDS elevates the water’s ionic strength and major ion concentration. This results in a greater competitive effect between ions and F⁻ from soil exchange sites and mineral surfaces through the ion exchange process, thus reducing the adsorption potential of F⁻ and enhancing their mobilization (84, 88). Sodium and chloride ions, responsible for TDS/EC, also significantly correlated with F⁻. Calcite precipitation (decrease in Ca²⁺) enriches F⁻ in the groundwaters, hence the negative correlation between the two. Also, F⁻ and all other hydrochemical parameters in shallow aquifers surpassed their respective permissible limits compared to the deeper water samples. Shallow waters are easily accessible for human consumption and other activities, and therefore, raises concern over affecting human well-being.

Besides the groundwater hydrogeochemical conditions, the prevailing arid and semi-arid climate in the study region increases evaporation rates relative to humid areas. High rainfall inputs and subsequent dilution effect in humid climatic zones result in lower groundwater fluoride levels compared to drier environments. Also, the groundwater movement in arid/semi-arid regions is generally slow, thereby increasing the contact time between the water and rock, which further causes F⁻ enrichment in water (85, 89).

Excess groundwater fluoride incidence in arid and semi-arid regions is a common phenomenon (90, 91), as observed in the aquifers of Punjab. Furthermore, due to its geogenic origin, the concentration of F⁻ depends directly on the hydrogeochemical conditions. Also, relatively limited studies have forecast F⁻ levels in arid and semi-arid locations using hydrochemical characteristics. Therefore, developing a predictive model for determining F⁻ levels using water quality variables in locations lacking monitoring assessments is essential. This study proposed five different ML models (RF, MLP, SVM, ELM, and Xgboost) and determined the best-performing model based on the evaluation metrics (R² , RMSE, and MAE). Out of all the models, MLP trained extremely poorly for the dataset and is, therefore, unsuitable for making reliable predictions of groundwater fluoride concentration in our study area. This finding is contrary to Nafouanti et al. (28) and Gupta and Maiti (30), where MLP performed accurately in predicting F⁻ concentrations in the Datong basin (China) and Maharashtra (India), respectively. The poor performance of MLP in our study could be due to its incapability to extrapolate beyond the data used for training, which further leads to overfitting issues during the training phase (92, 93). The worst prediction performing model in both the training and testing stages was MLP, which contradicts Bui et al. (94). The MLP is based on neural network architecture that can generate more accurate results on a badly structured dataset than on tree-based models such as RF (94). On the contrary, the RF model overcomes overfitting issues by combining many trees, thereby free from bias, resulting in enhanced prediction performance (32). Regardless, RF performed poorly in the training phase as well. Both MLP and RF generated values that deviated greatly from the original values, implying unsatisfactory performance. MLP tends to overfit from the training data, interfering with its ability to infer the remaining data (test/cross-validation dataset) (28). Gupta and Maiti (30), in their work, also stated that MLP, RF, and SVM are less effective in uncovering the intricate non-linear association between the target and predictor variables. Unsatisfactory training values of MLP, RF, SVM, and Xgboost could have resulted from a very wide variability in the range of both the target (output) and predictor (input) variables within the compiled dataset. This adds a limitation to model fitting in our study, resulting in inaccurate prediction results. Data pre-treatment involving outlier suppression and logarithmic transformation can be a possible solution to further improve the prediction accuracies (30). However, these pre-processing steps on the raw dataset and their influence on the model performance need further evaluation. Gupta and Maiti (30) also emphasized on the limited prediction efficiency of the ELM model due to its design and direct inverse in estimating the bias and weights. Irrespective of this fact, the ELM model with relatively higher R² value and low RMSE and MAE values in the training phase in comparison to the remaining models implies good generalization capability for our dataset without undergoing much pre-processing, unlike other models tested in this study.

Classification-based ML models have been commonly used in groundwater fluoride level prediction (10, 26, 27, 95, 96), with very few studies on regression-based prediction modeling of the same (30, 31). The complexity and accuracy of the datasets, diverse algorithm architectures, and type and number of input parameters significantly influence the performance of the models, and therefore, there is no universal agreement on which ML model performs the best for all prediction-related studies (94). For instance, groundwater fluoride concentration in the Datong basin, China, was modelled using RF, Linear regression (LR), and MLP-based Artificial neural network (ANN), where RF proved to be the best prediction model (28). Similarly, the RF model displayed higher prediction accuracy for other contaminants, such as nitrate, than the enhanced regression tree, classified regression tree, and multiple linear regression (97). Khosravi et al. (98), instead, reported that the M5P model had the highest predictive power than Instance Based Learner (IBK), KStar, Locally Weighted Learning (LWL), and Regression by discretization (RBD) that were tested for predicting F⁻ in the aquifers of Maku plain in Iran. Similarly, F⁻ levels in the groundwater of Sindhudurg district in Maharashtra, India, were predicted using six different models, out of which ELM yielded the most unsatisfactory results (30). On the contrary, Barzegar et al. (31) compared the performance of three different models and determined ELM to be the best for forecasting F⁻ in the Maku Valley of Iran, which is in accordance with the findings of the current study. Therefore, it is advisable to test diverse algorithms with the same dataset and assess their performance in terms of prediction before selecting the best.

The ELM model has a simple architecture and an uncomplicated training process and is generally known for its efficient computational power, requiring fewer hyperparameters for model tuning and training. The parameters of the hidden layer in this model do not require manual adjustments and are also independent of the input data. It only determines the weights of the output analytically and thus has rapid learning speed and lower computation complexity (99) than the other models proposed in this work. Additionally, ELM also has good generalization capability for high dimensional datasets by initialization of weights and biases stochastically to avoid overfitting problems and thus making the model more robust (100). The superior predictive performance of ELM over other models was also confirmed in other works (69, 101–104).

It is worth nothing that most of the studies conducted in the Indian subcontinent have reported RF and MLP to be the best models for predicting F⁻ concentrations in groundwater. However, these models did not take into consideration the groundwater-physicochemical parameters and used only continuous variables such as climate, soil, geological, and topological parameters as predictors (10, 27, 96). Machine learning algorithms are designed to perform both classification and regression-related tasks. Classification-based ML models have been commonly applied for studies in India that facilitated in forecasting the contamination-risk prone areas (105–107). Nevertheless, it is equally essential to develop models for predicting the concentration of the contaminants based on the driving factors that directly influence its enhancement and mobility. In this context, regression-based modeling will prove to be much more beneficial than classification models. This study attempted to achieve this goal; therefore, such contrasting results could be due to these reasons.

The different variables influencing groundwater fluoride contamination in any region are complex and require an in-depth understanding to identify the potential parameters for proper groundwater resource management. The variable importance ranking in the present work highlighted that TDS, EC, Na⁺, Cl⁻, and Ca²⁺ were the most crucial factors and were highly correlated with F⁻ content in the study region. The increase in TDS and EC results in increased ionic strength and higher concentration of major ions dissolved in water. These factors enhance competition between ions and F⁻ from mineral surfaces and soil exchange sites through the ion-exchange process, which further minimizes the adsorption of F⁻ and makes them more mobile (84, 88). Sodium, one of the important parameters responsible for EC and TDS, forms compounds with F⁻, such as NaF, which further dissolves in water and becomes more mobile (108). Other factors, such as Cl⁻ and Ca²⁺, contributed significantly to the model performance at varying degrees. The primary source of F⁻ in this region is fluorite mineral (CaF₂), which undergoes dissolution further releasing F⁻ and Ca²⁺ and the latter precipitates in the presence of excess bicarbonate (HCO₃ ⁻), thereby resulting in free F⁻ ions (109). On the contrary, the chloride ion undergoes ionic exchange with F⁻ from the aquifer substrate, bringing about the discharge of F⁻ ions from these surfaces. Furthermore, the significant contributions of the top 5 variables, i.e., TDS, EC, Na⁺, Cl⁻, and Ca²⁺, in all the five models, irrespective of their prediction accuracies, indicate their potential in determining F⁻ levels in regions lacking groundwater quality monitoring practices. However, the variability in the accuracy post-tuning of the different models might have impacted the outcomes. In other words, data quality, input parameters, hyperparameter tuning process, and varying algorithm architecture play a significant role in the prediction of the target variable. Moreover, the top contributors and their influence on groundwater F⁻ concentration are clearly highlighted from their relative importance scores in all the models, particularly in the ELM model ( Figure 6 ). In addition to this, the maximum prediction accuracy of ELM relative to other proposed models ( Table 3 , Figures 4 , 5 ) makes it an acceptable method for groundwater quality assessment investigations.