  The hadronic contribution to the muon anomalous magnetic moment
$a_\mu=(g_\mu-2)/2$ has to be determined at the per-mille level for the
Standard Model prediction to match the expected final uncertainty of the
ongoing E989 experiment. That is 3 times better than the current precision from
the dispersive approach, and 5-15 times smaller than the uncertainty based on
the purely theoretical determinations from lattice QCD. So far the
stumbling-block is the large statistical error in the Monte Carlo evaluation of
the required correlation functions which can hardly be tamed by brute force. In
this talk we present our proposal to solve this problem by multi-level Monte
Carlo integration, a technique which reduces the variance of correlators
exponentially in the distance of the fields. We report the results of our
feasibility tests for the computation of the Hadronic Vacuum Polarization on a
lattice with a linear extension of 3~fm, a spacing of 0.065 fm, and a pion mass
of 270 MeV. Indeed the two-level integration makes the contribution to the
statistical error from long-distances de-facto negligible by accelerating its
inverse scaling with the cost of the simulation. These findings establish
multi-level Monte Carlo as a solid and efficient method for a precise lattice
determination of the hadronic contribution to $a_\mu$.
