Action-based distribution function modelling for constraining the shape of the Galactic dark matter halo

We estimate the 3D density profile of the Galactic dark matter (DM) halo within r ≲ 30 kpc from the Galactic centre by using the astrometric data for halo RR Lyrae stars from Gaia DR2. We model both the stellar halo distribution function and the Galactic potential, fully taking into account the survey selection function, the observational errors, and the missing line-of-sight velocity data for RR Lyrae stars. With a Bayesian method, we infer the model parameters, including the density flattening of the DM halo q, which is assumed to be constant as a function of radius. We find that 99 per cent of the posterior distribution of q is located at q > 0.963, which strongly disfavours a flattened DM halo. We cannot draw any conclusions as to whether the Galactic DM halo at $r \lesssim 30 \, \mathrm{kpc}$ is prolate, because we restrict ourselves to axisymmetric oblate halo models with q ≤ 1. Our DM density profile might be biased especially in the inner few kpc, due to the uncertainty in the baryonic distribution. Our result is in tension with predictions from cosmological hydrodynamical simulations that advocate more oblate (<q> ~ 0.8 ± 0.15) DM haloes within ${\sim}15{{\ \rm per\ cent}}$ of the virial radius for Milky-Way-sized galaxies. An alternative possibility, based on our validation tests with a cosmological simulation, is that the true value q of the Galactic halo could be consistent with cosmological simulations but that disequilibrium in the Milky Way potential is inflating our measurement of q by 0.1-0.2. As a by-product, our model constrains the DM density in the Solar neighbourhood to be $\rho _{\mathrm{DM},\odot } = (9.01^{+0.18}_{-0.20})\times 10^{-3}{\,\rm M_\odot} \mathrm{pc}^{-3} = 0.342^{+0.007}_{-0.007}$ GeVcm^-3, consistent with other recent measurements.