An integral model of the deep bed filtration process has been developed. It incorporates pore and particle size distributions, as well as the particle residence time distribution in the framework of the continuous time random walk theory. Numerical modeling is carried out to study the factors influencing breakthrough curves and deposition profiles for the deep bed filtration systems. Results are compared with a large set of experimental observations. Our findings show that highly dispersed breakthrough curves, e.g. those with early arrivals and large ending tails, correspond to large dispersion coefficients. For such cases the elliptic equation excels the advection dispersion equation in both fitting breakthrough curves and predicting deposition profiles related to natural or highly heterogeneous porous media. The deposition hyperexponentiality can be caused by the following three mechanisms: particle population in connection with the distribution of the filtration coefficients, heterogeneity in connection with non-Fickian transport, and heterogeneity in connection with the spatial distribution of the filtration coefficients. The influence and interaction of all three mechanisms have been analyzed in numerical computations and by comparison to several sets of experimental data.