Gully erosion is one of the most extreme land degradation processes that exhibit spatial and temporal variation depending on topography, soil, climate, and land use and management characteristics. This study investigated the impact of changes in rainfall, land use/land cover (LULC), and land use management (LUM) practices on gully erosion in two midland watersheds (treated Kecha and untreated Laguna) in the Upper Blue Nile basin of Ethiopia by using the LANDPLANER model in combination with intensive field measurements and remote sensing products. We simulated gully erosion under past (in 2005), present (in 2021), and three potential future curve number conditions, each time under four rainfall scenarios (10, 30, 60, and 100 mm) using the dynamic erosion index (e), static topographic (esp), and erosion channel (esp_channel) thresholds. Density plot analyses showed that gullies frequently occur in low-lying gentle slope areas with relatively higher curve number values. The best predictions of gullies identified through true positive rates (TPR) and true negative rates (TNR) were achieved considering the static esp_channel > 1 for Kecha (TPR = 0.667 and TNR = 0.544) and the dynamic e > 0.1 for 60 mm of rainfall in Laguna (TPR = 0.769 and TNR = 0.516). Despite the 10 mm rainfall having negligible erosion-triggering potential in both watersheds, the 60 and 100 mm rainfall scenarios were 4–5 and 10–17 times, respectively, higher than the 30 mm rainfall scenario. While the LULC change in the untreated Laguna watershed increased the impact of rainfall on gully initiation by only 0–2% between 2005 and 2021, the combination of LULC and LUM significantly reduced the impact of rainfall in the treated Kecha watershed by 64–79%. Similarly, the gully initiation area in Kecha was reduced by 28% (from 33% in 2005 to 5% in 2021) due to changes in LULC and LUM practices, whereas Laguna showed little increment by only 1% (from 42% in 2005 to 43% in 2021) due to LULC change. In addition, the future predicted alternative land use planning options showed that gully initiation areas in Laguna could be reduced by 1% with only LULC conversion; 39% when only LUM practices were implemented; and 37% when both were combined. These results indicate that LUM practices outweigh the impact of LULC on gully erosion in the studied paired watersheds. Overall, LANDPLANER successfully simulated the spatio-temporal variation of gully erosion with scenario-based analyses and hence can be used to predict gullies in the study area and other data-scarce regions with similar agro-ecological settings.