Stress changes associated with reservoir depletion are often observed in the field. The four-dimensional stress evolution within and surrounding drainage areas can greatly affect completion of infill wells and refracturing. To accurately predict the four- dimensional stress distribution of shale gas reservoir, a coupled fluid- flow/geomechanics model considering the microscopic seepage mechanism of shale gas and the distribution of complex natural fractures (NFs) is derived based on the Biot's theory, the embedded discrete fracture model (DEFM) and finite volume method (FVM). Based on this model, the four-dimensional stress prediction can be realized considering the mechanism of adsorption, desorption, diffusion and slippage of shale gas and the random distribution of NFs. The results show that in the process of four- dimensional stress evolution, there will be extremes of σxx, σyy, σxy, Δσ, α and stress reversal area at some time, and the time of occurrence of extremes is different at different positions. The key to determine this law is the pore pressure gradient with spatio-temporal evolution effect. Different microscopic seepage mechanisms have great influence on the storage and transmission of shale gas, which leads to great differences in the distribution of reservoir pressure and four-dimensional stress. The influence of microscopic seepage mechanism should be considered in the process of four- dimensional stress prediction. The larger the initial stress difference is, the more difficult the stress reversal is. When the initial stress difference exceeds a certain limit value, the stress reversal phenomenon will not occur in the reservoir. This research is of great significance for understanding the four-dimensional stress evolution law of shale gas reservoir, guiding completion of infill wells and refracturing design.