*4.2. Result of FR*

The results of the FR analysis in identifying the relationship of land subsidence occurrence with the conditioning factors are summarized in Table 3. Two out of five altitude classes had the highest probability (FR > 1.0), with 1119 to 1137 m being the most correlated class with land subsidence, followed by the altitude class lower than 1119 m. The results of the slope angle analysis showed that slopes ranging between 4.5 and 6.8 degrees had the highest FR (1.11). Further, a TWI class lower than 4.84, profile curvature higher than 0.0029, convex plan curvature class, and flat (F) slope aspect had the most influence on LSS for each corresponding factor. Land cover analysis results indicated that forest and urban classes had the highest probability of land subsidence occurrence, with FR values of 1.17 and 1.04, respectively. Distance to a stream of between 50 to 100 m had the highest FR, and the class of lower than 50 m had a considerable correlation; a stream density higher than 2.68 had the highest correlation with land subsidence and the 1.23 to 1.92 class, which also had a considerable FR value. For distance to road, the 0 to 100 m class had the highest FR followed by the 100 to 200 m classes. Finally, groundwater drawdown ranging from 55 to 83 m and from 28 to 55 m had higher impacts on land subsidence occurrence.


**Table 3.** Relationship between land subsidence occurrence and conditioning factors using the FR model.

### *4.3. Result of Hybrid Models*

In the course of implementation of the hybrid models, 70% of the land subsidence points (210 locations) were used for training with values 1, and the same number of randomly selected non-subsidence points were taken into account with 0 values in the training phase. For the test dataset, 30% of the subsidence inventory (90 locations) with a value of 1 was used, with 90 randomly assigned points with a value of 0. The training datasets were used to calibrate the weights of the membership functions. The testing dataset was used to evaluate the performance of the trained ANFIS ensemble models. Hybrid models were implemented in MATLAB 2017b software. The parameters used in meta-heuristic algorithms are presented in Table 4. The prediction power of ANFIS and the two hybrid models with the training dataset (target) along with the comparison of the output and target testing dataset is shown in Figure 6.


**Figure 6.** The target and output values for training and validation datasets of (**a**) ANFIS, (**b**) ANFIS-GWO, and (**c**) ANFIS-ICA.

The RMSEs of the training and testing phases were calculated and are shown in Table 5. The two ensemble models enhanced the ANFIS model, and the ANFIS-ICA outperformed the ANFIS-GWO with an RMSE of 0.276 in the training phase and 0.3199 in the validation and testing phase. The ANFIS-GWO yielded an RMSE of 0.313 and 0.3217 in training and validation phases, respectively. Finally, the ANFIS model resulted in 0.323 in training and 0.34 in the validation phase.

**Table 5.** The comparison of model performance.


The convergence results of the two ANFIS-ICA and ANFIS-GWO ensemble models up to 1000 iterations are shown in Figure 7. ANFIS-ICA had a better convergence value (0.276) than ANFIS-GWO (0.313). The lowest amount of the cost function (RMSE) indicates the best cost and thus the best performance in predicting the results.

**Figure 7.** The convergence graph of the objective functions.

#### *4.4. LSSM Using ANFIS and Its Optimized Models*

The original ANFIS and its optimized ensembles in this research were trained and used to estimate land subsidence susceptibility in the study area. Susceptibility modelling and estimation were all carried out in MATLAB 2017b and were then exported to ArcGIS 10.3 software to classify and generate the susceptibility maps. Land susceptibility index was classified into five classes, very high, high, moderate, low, and very low, based on a natural break classification scheme [22,41]. Figure 8 presents the generated classified subsidence susceptibility maps obtained from ANFIS, ANFIS-GWO, and ANFIS-ICA. As can be seen, all the output subsidence susceptibility maps are similar and consistent with each other, particularly the ones for ANFIS-ICA and ANFIS-GWO. Moreover, the map based on ANFIS-ICA is much smoother than the others.

#### *4.5. Validation*

The ROC curves were calculated for all LSS maps using the test data. Figure 9 demonstrates the comparison of AUC for all the models used. The results showed that the ANFIS-ICA had the highest prediction accuracy (0.932), followed by the ANFIS-GWO (0.926) and ANFIS (0.908). This proves that the combination of the ANFIS model with meta-heuristic algorithms such as GWO and ICA can significantly improve the output land subsidence susceptibility maps in comparison to ANFIS alone.

**Figure 8.** The LSS maps of the study area using (**a**) ANFIS, (**b**) ANFIS-ICA, and (**c**) ANFIS-GWO.

**Figure 9.** The ROC curves for the LSSMs of the models and their AUC.

#### **5. Discussion**

Land subsidence is the slow vertical lowering of the Earth's surface, posing a serious threat to both the environment and human life. Recently, there has been an increasing interest in land subsidence analysis and monitoring in Iran as it is one of the highest subsidence-prone countries [17,18,20,24]. Natural hazard vulnerability analysis using machine learning algorithms (i.e., ANFIS) has shown promising results. Therefore, in this research, the focus was to employ the ANFIS model in combination with meta-heuristics in land subsidence susceptibility mapping.

Land subsidence inventories are necessary for accurate subsidence susceptibility analysis. The use of remote sensing SAR data is suitable for providing subsidence inventories due to their wide availability, independence from fieldwork, time and cost efficiency, frequent repeatability over time, and, especially, high precision [6]. In this work, the PS-InSAR technique with its millimetric precision was employed to determine the subsided areas in the region of interest to form the inventory data used for training and testing the LSS models.

Important conditioning factors for determining land subsidence prone areas were identified and collected based on either the literature or availability of data. The FR model was used to evaluate the correlation and influence of the factors. The results showed that all the factors employed in this paper had a considerable effect on LSS in Shahryar County. Among all the factors, the flat (slope aspect) area had the highest FR value (3.7), indicating high subsidence susceptibility in flat areas. The slope angle is related to the hydro physiographic characteristics that can influence the water infiltration rate and the volume and velocity of the Earth's surface flow [13]. Altitude and groundwater drawdown were the best predictors of land subsidence in this study, followed by stream density and distance to stream. Rahmati et al. and Arabameri et al. [2,37] also found that groundwater drawdown had a greater impact on land subsidence. Other topo-hydrographic factors, such as stream density and distance to stream are indirectly related to LS as they impact groundwater recharge and infiltration [2,40] and, as can be seen in the results, the land

areas closer to streams and with a certain stream density were more susceptible to LS. In terms of altitude, the lower lands were more prone to subsidence as the class 1137 to 1119 m and those lower than 1119 m had the highest FR. TWI, plan curvature, and profile curvature are among secondary topographic derivatives indirectly influencing LS [2,40,71]. These factors were not among best predictors of subsidence in the study area, which may be due to smooth and low altitude changes in the study area. The FR analysis showed that a TWI lower than 4.84 was strongly correlated with LS. In a similar study [40], lower TWI values (i.e., 2.54 to 8) have been reported to be more prone to subsidence. Positive and convex plan and profile curvatures had the highest FR value, as reported in [40]. Cropland and urban land cover types exist in lower altitude and flat areas. The main water source of the area is groundwater; therefore, more extraction of water in recent years as a result of population increase has caused the subsidence rate to increase. Previous studies have stressed the impact of groundwater extraction on subsidence occurrence [72,73]. Regarding the distance to road factor, the closer to a road, the greater the land subsidence risk, which can be due to closeness to urban land cover and thus indirectly related to the subsidence phenomenon.

Two novel meta-heuristic algorithms, GWO and ICA, were used to optimize the rules and parameters of the ANFIS model. Both of these evolutionary algorithms belong to swarm intelligence. The results showed that ICA had a slower convergence rate than GWO; however, it had better performance. In order to evaluate the prediction power and accuracy of the models, RMSE and ROC criteria were used. The RMSE is simply based on error assessment, whereas ROC is based on true positive (TP), false positive (FP), true negative (TN), and false negative (FN), which is more appropriate for comparison [42]. ANFIS-ICA had the lowest RMSE in both the training (0.276) and testing (0.3199) phases, followed by ANFIS-GWO and ANFIS alone. According to the AUC-ROC results, the ANFIS-ICA model was more accurate (0.932), followed by the ANFIS-GWO model (0.926) and the ANFIS model (0.908). It can be seen that the use of machine learning algorithms resulted in higher prediction accuracy since ANFIS alone yielded a suitable performance compared to other statistical methods in other studies. It can also be concluded that optimization of the ANFIS algorithm by meta-heuristics improves its results considerably. This was also reported in cases of other applications [64,74]. The results showed that the ICA algorithm was more accurate than the GWO algorithm in optimization of the ANFIS model. The advantages of the ICA algorithm are high convergence speed and the ability to optimize functions with a large number of variables [75]. The GWO algorithm has a small number of disadvantages, including a low solving accuracy, poor local searching ability, and slow convergence rate [60].

The output land subsidence susceptibility maps of the three models were similar and in line with each other. However, the map produced by the ANFIS-ICA was smoother than that of the other two. As could be observed, the high-risk areas were predicted where the groundwater extraction was higher, elevation was lower, and agricultural land use was higher. This is because the main source of income in the study area is agriculture. Further, the population has increased; therefore, food production has stressed the groundwater, the main water source of the area, and thus the land subsidence risk has become higher in those areas. Further, it is evident that the subsidence trend is gradually reaching towards the urban part of the Shahryar County, posing a serious threat to settlements and human life. The generated LSSMs in this paper can benefit authorities and decision-makers to identify subsidence-prone areas regarding environmental and urban management.

#### **6. Conclusions**

Land subsidence is an important issue in Iran due to the semi-arid and arid climate and excessive groundwater extraction. Therefore, modeling, simulation, and risk mapping offer valuable knowledge of environmental geohazards. GIS-based predictions have proved to be essential for authorities in terms of planning and decision-making. In this work, we used remote sensing SAR data and the PSInSAR technique to create a land subsidence inventory of the study area as a high-precision tool with a low cost and frequent

reproducibility. Since machine learning tools have shown appropriate performance in modeling and mapping hazard susceptibility, the ANFIS model was used in this research to map the land subsidence risk in Shahryar County, Tehran province, Iran. Another objective of this paper was to investigate the effect of optimization of the ANFIS model through meta-heuristics. Two novel evolutionary algorithms, namely, GWO and ICA, were used to create ensemble models. The results of the three models in both training and testing phases were assessed by RMSE. In both phases, ANFIS-ICA had the lowest RMSE, followed by ANFIS-GWO and ANFIS alone. AUC-ROC analysis was also used for model evaluation, and its results indicated that ANFIS-ICA had the best prediction performance (0.932), followed by ANFIS-GWO (0.926) and ANFIS (0.908). To conclude, the results overall showed the applicability of the ANFIS machine learning algorithm in land subsidence susceptibility mapping and the effectiveness of its ensembles with metaheuristic algorithms. The methodology used is reproducible and can be applied to other regions with different environmental parameters to test the modelling performance. Further studies should be applied using other machine learning and deep learning algorithms to compare their prediction accuracy. In addition, future research can focus on developing risk monitoring and early-warning frameworks.

**Author Contributions:** Conceptualization, B.R. and S.V.R.-T.; Data curation, B.R. and S.V.R.-T.; methodology, S.V.R.-T., F.F. and B.R.; software, D.P., F.F. and S.V.R.-T.; validation, S.V.R.-T., F.F. and B.R.; formal analysis, B.R., S.V.R.-T., F.F., A.S.-N. and D.P.; investigation, B.R.; resources, F.F.; writing original draft preparation, B.R.; writing—review and editing, B.R., S.V.R.-T. and F.F.; visualization, S.V.R.-T. and B.R.; supervision, A.S.-N. and D.P.; project administration, F.F., A.S.-N. and D.P. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Acknowledgments:** The authors would like to thank the European Space Agency (ESA) for freely providing Sentinel-1 satellite imagery.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

