Gully head erosion significantly contributes to land degradation, and affects gully dynamics on the Loess Plateau of China. Modeling with a gully head erosion susceptibility map (GHEM) is an essential step toward appropriate mitigation measures. This study evaluates the effectiveness of two varieties of advanced data mining techniques -a bivariate statistical model (certainty factor (CF)) and a machine learning model (random forest (RF)) for the accurate mapping of gully head erosion susceptibility taking the Dongzhi Loess Tableland of China as an example at a regional scale. A database comprising 11 geographic and environmental parameters was extracted with 415 spatially distributed gully heads, of which 70% (291) were selected for model training and 30% (124) were used for validation. An accuracy evaluation using the area under the curve (AUC) value revealed that the CF model was 84.1% accurate, while the AUC value of the RF model map was 88.8% accurate. According to the RF method used to assess the relative significance of predictor variables, the most significant factors influencing the spatial distribution of the GHEM were the slope angle, slope aspect, topographic wetness index, and slope length. The GHEM can ultimately aid in decision making with respect to soil planning and management and sustainable development of the study area, which can be applied to other similar loess regions. (c) 2021 Elsevier B.V. All rights reserved.