The assessment and mapping of the risk of soil salinization can contribute to sustainable land planning aimed at mitigating soil degradation and increasing crop pro- duction. A probabilistic approach, based on multivariate geostatistics was used to model the spatial variation of soil salinization risk at the landscape scale and to delineate the areas at high risk. The study site is a citrus growing area in south-eastern Sardinia (Italy). Electrical conductivity (ECe), exchangeable sodium percentage (ESP), pH and ‘total clay + fine silt content’ (FIN), were measured in the topsoil (0–40 cm). The method requires indicator coding, which transforms measured data values into a binary variable according to critical thresholds. These latter were set to: 4 dS m-1 for ECe, 10% for ESP, 8 for pH, and 40% for ‘total clay + fine silt content’. To determine the probability of exceeding these critical values, multi-collocated indicator cokriging was used. Factorial kriging was also applied to identify one regionalized factor that summarizes the effects of the selected variables on soil salinization. Maps of each soil indicator and regionalized factor were produced to show the areas at risk of salinization. The results are valuable for planning the management of salinity.