The residual stresses and scalar fluxes are required to close the momentum and scalar transport equations in simulations of turbulence that are not fully resolved in space. In stratified turbulence, the stress and fluxes are statistically anisotropic unless the smallest resolved length scale is smaller than the Ozmidov scale and the buoyancy Reynolds number is sufficiently high for there to exist a range of scales that is statistically isotropic. In this work, a tensorial basis set is derived analytically that potentially contains sufficient information about the anisotropic interaction between resolved and residual scales. The residual stress tensor is evaluated by filtering data from direct numerical simulations of homogeneous stratified turbulence, with unity Prandtl number, resolved on up to 8192 x 8192 x 4096 grid points along with an isotropic homogeneous case resolved on 81923 grid points. Five approximations for the residual stress tensor are derived by projecting the tensor onto different tensorial basis sets. It is found that an approximation for the residual stress which is based on a combination of the rate-of-strain and a tensor related to energy redistribution is a very good starting point for models involving just two coefficients. An extension for eddy viscosity models based on a non-linear tensor is proposed and tested on flow cases with different levels of density stratification. Different eddy viscosity models are mixed with the proposed extension on LES simulations. Results confirm the inability of pure eddy viscosity models to produce the expected anisotropy. It is shown that the proposed extension improves results in all cases, especially the correct reproduction of anisotropy.