Three dimensional soil textural structure in a township was conditionally simulated using a transition probability-based indicator geostatistical method based on 270 soil texture samples from 27 profiles. Additionally the distribution of soil profiles lacking clay interlayers (indicating high irrigation water and nutrient leaching risk) was analyzed using 500 realizations from the simulation. The results indicated that the simulation could predict the soil texture distribution with low uncertainties using the existing data, and the predicted soil map (0-10 cm) formed by the maximum probable soil textures also exhibited a good agreement with the legacy soil survey map. For water and nutrient leaching risk analysis, the areas lacking clay interlayer could be located; however, their distribution was still highly uncertain if based only on the existing sampling data. That means supplementary sampling in future is required for the risk assessment, and the existing study can help to optimise the sampling points and their distribution. Generally, the transition probability-based geostatistical simulation, as a stochastic conditional simulation method, exhibited its potential in soil texture spatial reproduction and related risk assessment.