The structure of soils is often heterogeneous with layered strata having distinct permeabilities. An advanced mathematical and numerical coupled model of elastic wave propagation in poroelastic multi-layered soils subjected to subsoil water infiltration is proposed in this study. The coupled model was based on the introduction of an inhomogeneous functionally graded fluid-saturation of the considered soil depending on the infiltration time, which was evaluated employing Richards’ equation. The time-harmonic solution was formulated in terms of the Fourier transform of Green’s matrix and the surface load that excites the vibration. The convergence and efficiency of the proposed approach are demonstrated. An example of dispersion curves for partially saturated porous strata made of loam, sand, and rock at different infiltration times is provided, and it is shown that the characteristics of the surface acoustic waves change with time, which can be further used for inverse problems’ solution.