Based on the upper bound limit analysis theorem and the pseudo-dynamic approach, a modified horizontal slice method is proposed to calculate the work rates of soil gravity and seismic inertia forces, both characterized by evident nonlinear distribution features. Traditional plasticity theory is extended to unsaturated soils to account for the influence of soil suction on the resistance forces acting on piles. The proposed method is validated by comparing it with analytical solutions. A series of parametric analysis are conducted regarding the impacts of soil suction and pile design parameters on soil slope stability. The results show that the slope stability is underestimated by 20% ~ 50% when soil suction is not considered in practice. The suction effect depends on seismic excitations and is significantly enhanced when the seismic wave approaches its peak. The horizontal, vertical acceleration coefficients, and the soil shear modulus all have a great influence on seismic slope stability. When the pile spacing D1/dp ≤ 4, pile reinforcement becomes obvious and increases rapidly as the internal friction angle increases. The optimum pile location commonly lies in the middle and upper part of the slope. The optimum pile location gradually decreases, and the pile reinforcement is slightly reduced as the internal friction angle increases.