Gully erosion susceptibility mapping is an essential land management tool to reduce soil
erosion damages. This study investigates gully susceptibility based on multiple diagnostic analysis,
support vector machine and random forest algorithms, and also a combination of these models,
namely the ensemble model. Thus, a gully susceptibility map in the Kondoran watershed of Iran
was generated by applying these models on the occurrence and non-occurrence points (as the target
variable) and several predictors (slope, aspect, elevation, topographic wetness index, drainage density,
plan curvature, distance to streams, lithology, soil texture and land use). The Boruta algorithm was
used to select the most effective variables in modeling gully erosion susceptibility. The area under
the receiver operating characteristic curve (AUC), the receiver operating characteristics, and true skill
statistics (TSS) were used to assess the model performance. The results indicated that the ensemble
model had the best performance (AUC = 0.982, TSS = 0.93) compared to the others. The most effective
factors in gully erosion susceptibility mapping of the study region were topological, anthropogenic,
and geological. The methodology and variables of this study can be used in other regions to control
and mitigate the gully erosion phenomenon by applying biophilic and regenerative techniques at the
locations of the most influential factors.