In this paper, an efficient yet accurate method for the thermal stress analysis using a first order shear deformation theory(FSDT) is presented. The main objective herein is to systematically modify transverse shear strain energy through the mixed variational theorem(MVT). In the mixed formulation, independent transverse shear stresses are taken from the efficient higher-order zigzag plate theory, and the in-plane displacements are assumed to be those of the FSDT. Moreover, a smooth parabolic distribution through the thickness is assumed in the transverse normal displacement field in order to consider a transverse normal deformation. The resulting strain energy expression is referred to as an enhanced first order shear deformation theory, which is obtained via the mixed variational theorem with transverse normal deformation effect(EFSDTM_TN). The EFSDTM_TN has the same computational advantage as the FSDT_TN(FSDT with transverse normal deformation effect) does, which allows us to improve the through-the-thickness distributions of displacements and stresses via the recovery procedure. The thermal stresses obtained by the present theory are compared with those of the FSDT_TN and three-dimensional elasticity.