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 $(Zalpha)^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.