The research on the dynamic response of slopes under seismic action currently focuses on a single mainshock without considering the effect of aftershocks, and the spatial variability of material parameters is generally ignored. In this paper, the spatial variability of parameters is fully considered, and a reliability analysis framework based on the Newmark method and probaility density evolution method (PDEM) is proposed to quantify the effects of aftershocks and spatial variability on the dynamic reliability. First, the physical random function model, Copula function and the narrowband harmonic group superposition method are combined to generate the mainshock-aftershock sequence (MAS). In addition, the random field is generated based on the spectral representation and parameters are assigned to the finite element model based on the corresponding coordinates. Then, the permanent displacement of the slope considering the spatial variability of parameters subjected to the MAS was batch calculated based on the Newmark method, and the effects of the coefficient of variation of cohesion and friction angle (COV_C and COV_F), aftershock, and peak ground acceleration (PGA) on the permanent displacement of the slope were analyzed by the mean of the displacement. Finally, based on PDEM, the effects of COV and aftershocks on the dynamic reliability of the slope are explained from a probabilistic point of view. The results of the study show that the mean of displacements shows a gradual increase with an increase of coefficient of variation (COV). In contrast, the COV_F has a more pronounced effect on the slope displacement . In addition, the mean of displacement of the slope subjected to the MAS is greater compared to that of the single mainshock. If the spatial variability of parameters and the influence of aftershocks are ignored, the seismic performance of the slope would be overestimated.