We present STEEL a STatistical sEmi-Empirical modeL designed to probe the distribution of satellite galaxies in groups and clusters. Our fast statistical methodology relies on tracing the abundances of central and satellite haloes via their mass functions at all cosmic epochs with virtually no limitation on cosmic volume and mass resolution. From mean halo accretion histories and subhalo mass functions the satellite mass function is progressively built in time via abundance matching techniques constrained by number densities of centrals in the local Universe. By enforcing dynamical merging timescales as predicted by high-resolution N-body simulations, we obtain satellite distributions as a function of stellar mass and halo mass consistent with current data. We show that stellar stripping, star formation, and quenching play all a secondary role in setting the number densities of massive satellites above $M_*gtrsim 3times 10^{10}, M_{odot}$. We further show that observed star formation rates used in our empirical model over predict low-mass satellites below $M_*lesssim 3times 10^{10}, M_{odot}$, whereas, star formation rates derived from a continuity equation approach yield the correct abundances similar to previous results for centrals.