Statistical analysis of the gravitational anomaly in {\it Gaia} wide binaries

The exploration of the low acceleration $a<a_{0}$ regime, where $a_{0}=1.2 \times 10^{-10}$m s$^{-2}$ is the acceleration scale of MOND around which gravitational anomalies at galactic scale appear, has recently been extended to the much smaller mass and length scales of local wide binaries thanks to the availability of the {\it Gaia} catalogue. Statistical methods to test the underlying structure of gravity using large samples of such binary stars and dealing with the necessary presence of kinematic contaminants in such samples have also been presented. However, an alternative approach using binary samples carefully selected to avoid any such contaminants, and consequently much smaller samples, has been lacking a formal statistical development. In the interest of having independent high quality checks on the results of wide binary gravity tests, we here develop a formal statistical framework for treating small, clean, wide binary samples in the context of testing modifications to gravity of the form $G \to \gamma G$. The method is validated through extensive tests with synthetic data samples, and applied to recent {\it Gaia} DR3 binary star observational samples of relative velocities and internal separations on the plane of the sky, $v_{2D}$ and $r_{2D}$, respectively. Our final results for a high acceleration $r_{2D}<0.01$pc region are of $\gamma=1.000 \pm 0.096$, in full accordance with Newtonian expectations. For a low acceleration $r_{2D}>0.01$pc region however, we obtain $\gamma=1.5 \pm 0.2$, inconsistent with the Newtonian value of $\gamma=1$ at a $2.6 \sigma$ level, and much more indicative of MOND AQUAL predictions of close to $\gamma=1.4$.