(see Table \[Table:HFS\]) with theoretical values by @Dzuba2016a, we observe that our $A$ values match, while there is a roughly $10\%$ discrepancy in values of $B$. Itano (cited in Ref. [@Petrasiunas2012]) has previously computed the $A$ and $B$ constants for the $4f^{13}6s^2$ ($J=7/2$) state in $^{171}\mathrm{Yb}^+$ and $^{173}\mathrm{Yb}^+$. Itano has also used the MCDHF method, but his results are markedly different from ours. Since there are no details of calculations given in Ref. [@Petrasiunas2012], it is difficult to assess the reasons for this difference. We, however, point out that our MCDHF results are in a better agreement with experimental values. For example, the deviation is about 20% between his result and the experimental value $A(^{171}\textrm{Yb}) = 905$ MHz [@Taylor1999]. Multiplying our $A$ constant for $^{173}$Yb by the ratio $\mu(^{171}\textrm{Yb})I(^{173}\textrm{Yb})/\mu(^{173}\textrm{Yb})I(^{171}\textrm{Yb})$, we obtained $A=882$ MHz for $^{171}$Yb, which differs from the measurement [@Taylor1999] by only 3%.

Based on these comparisons we conservatively estimate the uncertainty of our MCDHF calculations to be $\sim 10\%$ for the magnetic dipole and electric quadrupole hyperfine interaction constants. This estimate is also consistent with that of Ref. [@Dzuba2016a], where they claimed a similar 10% theoretical uncertainty for these two constants using a different computational method. We assign 