Managing soil phosphorus is essential for agricultural production and environmental protection. This requires information on the phosphorus sorption capacity (PSC) of the soil. In this study, we map the PSC for Danish soils in four depths. Measuring PSC directly is expensive and time-consuming, and we therefore used a pedotransfer function based on oxalate-extractable aluminum (Al_o) and iron (Fe_o). We mapped Al_o, Fe_o, and their uncertainties using Quantile Regression Forests (QRF). We then calculated the uncertainties for PSC with a quasi-Monte Carlo complete combinatorial convolution (CCC) of the prediction quantiles for Al_o and Fe_o. The main factors for predicting Al_o were the parent material, topography and precipitation. In many areas, podzolization also affected the vertical distribution of Al_o. The main factors for Fe_o were the soil texture, organic matter and wetland areas. The average predicted PSC was 36 ± 9 mmol kg⁻¹ (±1 SD) in topsoil, but it was generally highest at 25–50 cm depth with a mean value of 37 ± 9 mmol kg⁻¹. The weighted root mean squared error of the mapped properties was 14 mmol kg⁻¹ for Al_o, 32 mmol kg⁻¹ for Fe_o and 19 mmol kg⁻¹ for PSC. The prediction accuracies were moderate at best, but the prediction quantiles were generally reliable. The mapped uncertainties were largest in wetland areas, while they were smallest in young loamy moraine deposits.