Oven-dry soil bulk density (BD) is a key soil parameter in biophysical models. Yet, its direct determination for large-scale modeling applications is limited by excessive efforts required for labor-demanding, time-consuming, expensive field campaigns and laboratory-based measurements. To circumvent these shortcomings, BD can be estimated using pedotransfer functions (PTFs) that, however, offer their optimal prediction capability if calibrated and validated within the area of interest. In this study, we exploited the availability of a dataset comprising 3,316 soil samples collected in the farmlands of Campania (a region of southern Italy) to develop regional PTFs for predicting BD using the Random Forest (RF) algorithm. RF was executed considering different combinations of seven soil and three terrain attributes with a 10-fold cross-validation approach to avoid performance overestimation. In light of the RF-based results, we further developed two new PTFs based on multiple linear regression equations. The first regression-based PTF was multiparametric and employed eight features (i.e., six soil properties and two terrain features as environmental covariates), whereas the second PTF was parsimonious and based on three easily available soil predictors. Both regression-based PTFs consistently outperformed 62 existing published PTFs. We also enhanced PTF prediction capabilities by employing regionalization through a clustering approach by grouping soil samples in ten land system classes. Finally, transferability of our models was tested using an external large independent dataset of 12,019 soil samples extracted from the European EU-HYDI database. The parsimonious PTF proved satisfactory prediction performance by corroborating results found in the Campania dataset.