  The largest hyperfine interaction coefficients in the hydrogen molecular ion
HD$^+$, i.e. the electron-proton and electron-deuteron spin-spin scalar
interactions, are calculated with estimated uncertainties slightly below 1~ppm.
The $(Z\alpha)^2 E_F$ relativistic correction, for which a detailed derivation
is presented, QED corrections up to the order $\alpha^3 \ln^2 (\alpha)$ along
with an estimate of higher-order terms, and nuclear structure corrections are
taken into account. Improved results are also given for the electron-proton
interaction coefficient in H$_2^+$, in excellent agreement with RF spectroscopy
experiments. In HD$^+$, a 4$\sigma$ difference is found in the hyperfine
splitting of the $(v,L)=(0,3) \to (9,3)$ two-photon transition that was
recently measured with high precision. The origin of this discrepancy is
unknown.
