The present paper aims to assess the effectiveness of four sub-grid stress models in large eddy simulation of natural convection in a buoyant square cavity. The four models used in the present study are, namely, dynamic Smagorinsky model, Vreman model, one equation eddy viscosity model and buoyancy modified one equation eddy viscosity model. We have considered a square cavity that is heated at the left end and kept cold at the right end with Rayleigh number (Ra) based on the height of the cavity equal to 1.58 x 10(9). A weakly compressible version of the Navier-Stokes equations, a non Boussinesq formulation, has been solved to simulate flow and heat transfer inside the buoyant cavity. The computational results have been compared with the experimental data reported in the literature. The present study shows that the Vreman model overpredicts the peak of the time-averaged velocity and underpredicts the peak of the root mean square fluctuations. The one equation model and buoyancy modified one equation model can predict the time-averaged field, root mean square field, turbulent shear stress and turbulent heat flux more accurately than those by the dynamic Smagorinsky and Vreman models. Some improvements in the results has been observed by applying the buoyancy modification to the one equation SGS model. An overall good agreement is observed between the computational and experimental results. Small deviations have been reported between the predicted second order turbulent statistics and experimental results. The physics of the mean flow field and turbulent structures in the buoyant cavity are discussed in detail. (C) 2016 Elsevier Ltd. All rights reserved.