Geostatistical simulation techniques are used to quantify uncertainty of spatial attributes of interest describing mineral deposits, petroleum reservoirs, hydrogeological horizons, environmental contaminants and so on. The majority of existing methods consider second-order spatial statistics and Gaussian processes, while the more advanced multiple point-based simulation approaches are algorithmic and do not consistently account for the high-order spatial relations in data. Recently, simulation techniques for complex and non-Gaussian spatially distributed variables have been developed, based on high-order spatial cumulants, and make no assumptions on data distribution or require data transformations. In this paper, the previous developments are extended and a new approach for the joint simulation of multiple correlated variables using high-order spatial statistics is proposed. The technique is based on a new algorithm for the decorrelation of correlated variables into factors, using the so-termed diagonal domination condition of high-order cumulants. The decorrelated factors are then simulated using high-order simulation and back-transformed into the initial correlated variables. The decorrelation using diagonal domination of high-order statistics is tested with a data set from a multielement iron ore deposit, and then a fully known two-dimensional data set with two correlated variables is used to demonstrate the practical intricacies of the proposed method.CITATION:Minniakhmetov, I and Dimitrakopoulos, R, 2014. Joint simulation of correlated variables using high-order spatial statistics for orebody modelling, in Proceedings Orebody Modelling and Strategic Mine Planning Symposium 2014 , pp 53-60 (The Australasian Institute of Mining and Metallurgy: Melbourne).