We present a calculation of the hindered M$1$ ${\rm Y}(2S)\to \eta b(1S)\gamma$ decay rate using lattice non-relativistic QCD. The calculation includes spin-dependent relativistic corrections to the NRQCD action through $O(v6)$ in the quark's relative velocity, relativistic corrections to the leading order current which mediates the transition through the quark's magnetic moment, radiative corrections to the leading spin-magnetic coupling and for the first time a full error budget. We also use gluon field ensembles at multiple lattice spacing values, all of which include $u$, $d$, $s$ and $c$ quark vacuum polarisation. Our result for the branching fraction is $\mathcal{B}(\Upsilon(2S)\to\eta_b(1S)\gamma) = 5.4(1.8)\times 10^{-4} $, which agrees with the current experimental value.