  The growth of supermassive black holes is strongly linked to their galaxies.
It has been shown that the population mean black-hole accretion rate
($\overline{\mathrm{BHAR}}$) primarily correlates with the galaxy stellar mass
($M_\star$) and redshift for the general galaxy population. This work aims to
provide the best measurements of $\overline{\mathrm{BHAR}}$ as a function of
$M_\star$ and redshift over ranges of $10^{9.5}<M_\star<10^{12}~M_\odot$ and
$z<4$. We compile an unprecedentedly large sample with eight thousand active
galactic nuclei (AGNs) and 1.3 million normal galaxies from nine high-quality
survey fields following a wedding-cake design. We further develop a
semiparametric Bayesian method that can reasonably estimate
$\overline{\mathrm{BHAR}}$ and the corresponding uncertainties, even for
sparsely populated regions in the parameter space. $\overline{\mathrm{BHAR}}$
is constrained by X-ray surveys sampling the AGN accretion power and
UV-to-infrared multi-wavelength surveys sampling the galaxy population. Our
results can independently predict the X-ray luminosity function (XLF) from the
galaxy stellar mass function (SMF), and the prediction is consistent with the
observed XLF. We also try adding external constraints from the observed SMF and
XLF. We further measure $\overline{\mathrm{BHAR}}$ for star-forming and
quiescent galaxies and show that star-forming $\overline{\mathrm{BHAR}}$ is
generally larger than or at least comparable to the quiescent
$\overline{\mathrm{BHAR}}$.
