In this paper, we couple a numerical method aimed at simulation of flow and heat transfer of nanofluids with stochastic modelling of input and material parameters. In order to simulate nanofluids, an in-house numerical method was developed, based on the solution of 3D velocity–vorticity formulation of Navier–Stokes equations. A fast Boundary-Domain Integral Method has been employed to solve the governing equations and set up the deterministic flow and heat transfer solver. The developed algorithm is used to simulate natural convection of a nanofluid in a closed cavity. The uncertainty present in the input parameters is propagated to the output of interest via the Stochastic Collocation Method. The stochastic mean, variance, and higher-order moments of the output values are presented. The non-intrusive nature of the Stochastic Collocation Method facilitates the previously validated deterministic code to remain unchanged. The stochastic analysis reveals that the uncertainty of input parameters influences the output results most in the areas where high flow field gradients appear.