  We use \texttt{GRUMPY}, a simple regulator-type model for dwarf galaxy
formation and evolution, to forward model the dwarf galaxy satellite population
of the Milky Way (MW) using the Caterpillar zoom-in simulation suite. We show
that luminosity and distance distributions of the model satellites are
consistent with the distributions measured in the DES, PS1 and SDSS surveys,
even without including a model for the orphan galaxies. We also show that our
model for dwarf galaxy sizes can simultaneously reproduce the observed {\it
distribution} of stellar half-mass radii, $r_{1/2}$, of the MW satellites and
the overall $r_{1/2}-M_\star$ relation exhibited by observed dwarf galaxies.
The model predicts that some of the observed faint stellar systems with
$r_{1/2}<10$ pc are ultra-faint dwarf galaxies. Scaling of the stellar mass
$M_\star$ and peak halo mass $M_{\rm peak}$ for the model satellites is not
described by a power law, but has a clear flattening of $M_\star-M_{\rm peak}$
scaling at $M_{\rm peak}<10^8\,M_\odot$ imprinted by reionization. As a result,
the fraction of low mass haloes ($M_{\rm peak} < 10^8\, M_\odot$) hosting
galaxies with $M_V<0$ is predicted to be 50% at $M_{\rm peak} \sim 3.6 \times
10^7\,M_\odot$. We find that such high fraction at that halo mass helps to
reproduce the number of dwarf galaxies discovered recently in the HSC-SSP
survey. Using the model we forecast that there should be the total of
$440^{+201}_{-147}$ (68\% confidence interval) MW satellites with $M_V < 0$ and
$r_{1/2} > 10$ pc within 300 kpc and make specific predictions for the HSC-SSP,
DELVE-WIDE and LSST surveys.
