An accurate analysis of casing stress distribution and its variation regularities present several challenges during hydraulic fracturing of shale gas wells. In this paper, a new analytical mechanical-thermal coupling method was provided to evaluate casing stress. For this new method, the casing, cement sheath, and formation (CCF) system was divided into three parts such as initial stress field, wellbore disturbance field, and thermal stress field to simulate the processes of drilling, casing, cementing, and fracturing. The analytical results reached a good agreement with a numerical approach and were in-line with the actual boundary condition of shale gas wells. Based on this new model, the parametric sensitivity analyses of casing stress such as mechanical and geometry properties, operation parameters, and geostress were conducted during multifracturing. Conclusions were drawn from the comparison between new and existing models. The results indicated that the existing model underestimated casing stress under the conditions of the geostress heterogeneity index at the range of 0.5–2.25, the fracturing pressure larger than 25 MPa, and a formation with large elastic modulus or small Poisson’s ratio. The casing stress increased dramatically with the increase of in situ stress nonuniformity degree. The stress decreased first and then increased with the increase of fracturing pressure. Thicker casing, higher fluid temperature, and cement sheath with small modulus, large Poisson’s ratio, and thinner wall were effective to decrease the casing stress. This new method was able to accurately predict casing stress, which can become an alternative approach of casing integrity evaluation for shale gas wells.