  Accurate redshift calibration is required to obtain unbiased cosmological
information from large-scale galaxy surveys. In a forward modelling approach,
the redshift distribution n(z) of a galaxy sample is measured using a
parametric galaxy population model constrained by observations. We use a model
that captures the redshift evolution of the galaxy luminosity functions,
colours, and morphology, for red and blue samples. We constrain this model via
simulation-based inference, using factorized Approximate Bayesian Computation
(ABC) at the image level. We apply this framework to HSC deep field images,
complemented with photometric redshifts from COSMOS2020. The simulated
telescope images include realistic observational and instrumental effects. By
applying the same processing and selection to real data and simulations, we
obtain a sample of n(z) distributions from the ABC posterior. The photometric
properties of the simulated galaxies are in good agreement with those from the
real data, including magnitude, colour and redshift joint distributions. We
compare the posterior n(z) from our simulations to the COSMOS2020 redshift
distributions obtained via template fitting photometric data spanning the
wavelength range from UV to IR. We mitigate sample variance in COSMOS by
applying a reweighting technique. We thus obtain a good agreement between the
simulated and observed redshift distributions, with a difference in the mean at
the 1$\sigma$ level up to a magnitude of 24 in the i band. We discuss how our
forward model can be applied to current and future surveys and be further
extended. The ABC posterior and further material will be made publicly
available at https://cosmology.ethz.ch/research/software-lab/ufig.html.
