  In the past few years, approximate Bayesian Neural Networks (BNNs) have
demonstrated the ability to produce statistically consistent posteriors on a
wide range of inference problems at unprecedented speed and scale. However, any
disconnect between training sets and the distribution of real-world objects can
introduce bias when BNNs are applied to data. This is a common challenge in
astrophysics and cosmology, where the unknown distribution of objects in our
Universe is often the science goal. In this work, we incorporate BNNs with
flexible posterior parameterizations into a hierarchical inference framework
that allows for the reconstruction of population hyperparameters and removes
the bias introduced by the training distribution. We focus on the challenge of
producing posterior PDFs for strong gravitational lens mass model parameters
given Hubble Space Telescope (HST) quality single-filter, lens-subtracted,
synthetic imaging data. We show that the posterior PDFs are sufficiently
accurate (i.e., statistically consistent with the truth) across a wide variety
of power-law elliptical lens mass distributions. We then apply our approach to
test data sets whose lens parameters are drawn from distributions that are
drastically different from the training set. We show that our hierarchical
inference framework mitigates the bias introduced by an unrepresentative
training set's interim prior. Simultaneously, given a sufficiently broad
training set, we can precisely reconstruct the population hyperparameters
governing our test distributions. Our full pipeline, from training to
hierarchical inference on thousands of lenses, can be run in a day. The
framework presented here will allow us to efficiently exploit the full
constraining power of future ground- and space-based surveys.
