This study focusses on probabilistic modelling of the bladed disc system and numerical estimation of the distributions of the response quantities of the system. Stochastic finite element model of the system consisting of all the assemblies and the hub is developed and reported. The spatial inhomogeneity of mistuned structures is modelled as non-Gaussian random field. Experimentally, the system parameters can be measured at the specified locations of the bladed disk structure. In this analysis, a synthetic data is generated which represent this measured data set. Further, Nataf transformation is implemented to each component of the data set to get the polynomial chaos expansion framework of the system parameters. Since, the random field of the system parameter is approximated as correlated random variables, Spearman’s rank correlation coefficient is used in this manuscript to obtain that correlation among the random parameters across the domain. The approximated probability density function obtained through the aforementioned methodology is compared with the target probability density function of the parameter using Kullback–Liebler (KL) entropy as a metric. Also, the same KL entropy is used as a metric to check the convergence of polynomial chaos terms in the expansion. Next, the proposed polynomial chaos method is integrated with commercial finite element software to quantify the propagation of randomness associated with system parameters into the response quantities. Subsequently, the statistical processing helps in estimating the probabilistic measure of the required response quantities. The results obtained through the conventional Monte Carlo (MC) simulations have been used as the benchmark to compare the response characteristics obtained through the proposed algorithm.