In the present study, we perform a set of numerical simulations in a moderately stable boundary layer with four types of subgrid-scale parameterization schemes and attempt to evaluate the error of the vertical flux for these schemes in terms of self-consistency on the basis of the Germano identity. If the effects of grid-scale components in higher wavenumbers are excluded from the analysis, the error estimated by the Germano identity is insensitive to the reference data utilized. The subgrid-scale flux, evaluated by the Smagorinsky model, tends to excessively weaken the positive temperature gradient at the top of the boundary layer with decreasing model resolution. The Deardorff and two-part models overestimate the subgrid-scale temperature flux at a coarser resolution, and the dynamic Smagorinsky model tends to underestimate both the subgrid-scale momentum and the temperature fluxes throughout the entire boundary layer. The underestimation of the subgrid-scale flux found in the dynamic Smagorinsky model could be attributed to a low correlation between the resolved and the parameterized components.