Open Access
Issue
A&A
Volume 710, June 2026
Article Number L39
Number of page(s) 8
Section Letters to the Editor
DOI https://doi.org/10.1051/0004-6361/202659498
Published online 30 June 2026

© The Authors 2026

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. This email address is being protected from spambots. You need JavaScript enabled to view it. to support open access publication.

1. Introduction

Scaling relations between the visible mass of galaxies and their observed kinematics provide fundamental insights into the interplay between baryons, dark matter (DM), and standard gravity. A classic example is the Faber-Jackson relation (FJR), which describes a correlation between the luminosity and stellar velocity dispersion of elliptical galaxies (Faber & Jackson 1976). An extension of the classic FJR is the baryonic Faber-Jackson relation (BFJR), which links the total baryonic mass (Mbar; stars plus gas) to the line-of-sight velocity dispersion (σlos; e.g., Sanders 2010; Famaey & McGaugh 2012). This is analogous to superseding the classic Tully-Fisher relation with the baryonic Tully-Fisher relation (McGaugh et al. 2000; Lelli et al. 2016, 2019). For elliptical galaxies, the FJR is generally thought to be a projection of the so-called fundamental plane (FP; Djorgovski & Davis 1987; Dressler et al. 1987), which adds the stellar effective radius (Re) as a third variable. The FP is expected from the Newtonian virial theorem (M ∝ R ⋅ σlos2), but its observed parameters show a “tilt” with respect to the virial parameters (e.g., Ciotti et al. 1996), which may be driven by variations in the stellar mass-to-light ratio, inner DM fractions, orbital anisotropy, or structural non-homology (e.g., Cappellari et al. 2006; Bolton et al. 2007; Cappellari et al. 2013a).

In the Λ cold dark matter (ΛCDM) cosmological paradigm, scaling relations must emerge from a combination of complex stochastic physical processes, including the hierarchical growth of DM halos, accretion and cooling of gas, star formation, feedback from supernovae and active galactic nuclei, and the feedback effects on the structural properties of DM halos (Pillepich et al. 2018). Within the DM framework, this leads to a “fine-tuning” problem: a precise coupling between baryons and DM is required to place DM-dominated dwarf galaxies and baryon-dominated giant galaxies on the same relations (e.g., Famaey & McGaugh 2012; Desmond & Wechsler 2017; Lelli 2022).

The main alternative to particle dark matter is modified Newtonian dynamics (MOND; or MilgrOmiaN Dynamics; Milgrom 1983). MOND postulates that the nonrelativistic laws of dynamics (gravity or inertia) are modified at accelerations below a characteristic scale, a0 = 1.2 × 10−10m s−2 (see Famaey & McGaugh 2012 and Banik & Zhao 2022 for reviews). In particular, the MOND virial theorem (Milgrom 1984, 2014a) implies a universal scaling for isolated, self-gravitating, virialized systems in the deep-MOND regime (internal accelerations g ≪ a0):

M bar σ 3 D 4 / ( G a 0 ) , Mathematical equation: $$ \begin{aligned} M_{\mathrm{bar} }\propto \sigma _{\rm {3D}}^4/(Ga_0), \end{aligned} $$(1)

where σ3D is the 3D mass-weighted velocity dispersion of the system. Equation (1) is markedly different from the Newtonian virial theorem because it does not contain any dependence on the characteristic size of the system. It is expected to hold only for systems in which g ≪ a0, while systems in which g ≫ a0 should follow the Newtonian relation with no DM. The general MOND paradigm predicts the slope of the BFJR to be exactly four. The value of the intercept depends on the specific MOND theory, but differences are of order O(1) (Milgrom 2014b, 2025).

To study the properties of the BFJR and FP, we compiled an unprecedented dataset spanning eight orders of magnitude in baryonic mass, from dwarf galaxies to massive ellipticals and galaxy groups. We examined how the parameters of the BFJR and FP vary with internal acceleration, providing a stringent test for both ΛCDM theories of galaxy formation and MOND.

2. Data analysis

We built a comprehensive sample covering different types of pressure-supported systems, including (1) 63 galaxy groups in the Local Supercluster (Makarov & Karachentsev 2011; Milgrom 2019; Sadhu & Tian 2024); (2) 1218 elliptical galaxies from the Mapping Nearby Galaxies at APO survey (MaNGA; e.g., Bundy et al. 2015; Duann et al. 2023); (3) 26 elliptical galaxies from the ATLAS3D survey (Cappellari et al. 2011); (4) 34 dwarf ellipticals in the Virgo cluster (Toloba et al. 2014); (5) 31 dwarf ellipticals in the Fornax cluster (Eftekhari et al. 2022); and (6) 28 dwarf spheroidals in the Local Group (Lelli et al. 2017). This combined sample spans the ranges Mbar ≃ 105 − 1013M and σlos ≃ 10 − 300 km s−1. A detailed description of each subsample and processing methods is provided in Appendix A.

We homogenized the definitions of key quantities across all datasets. The total baryonic mass is Mbar = Mstar + Mgas. For dwarf galaxies with negligible gas content, Mbar ≃ Mstar (Kroupa initial mass function). For ellipticals and galaxy groups, we included the hot gas mass with a median gas-to-baryon mass fraction of about 8% (Sadhu & Tian 2024). The line-of-sight velocity dispersion (σlos) was taken as the stellar velocity dispersion measured within one effective radius (σe) for individual galaxies, and as the velocity dispersion of member galaxies with available redshifts for galaxy groups, using the robust biweight estimator. The same member galaxies were also used to estimate the mean projected radius of the group.

For each system, we computed the Newtonian gravitational acceleration due to the observed distribution of baryons, gbar(r), assuming spherical symmetry. For elliptical and dwarf galaxies, the radial variation of gbar was computed by deprojecting a Sérsic profile, using the observed Sérsic index, effective radius, and stellar mass (Domínguez Sánchez et al. 2022; Krajnović et al. 2013; Toloba et al. 2014; Eftekhari et al. 2022). Then, we considered the median baryonic acceleration within one effective radius, ⟨gbar⟩, as the characteristic internal acceleration of each system. Comparisons with nonspherical estimates for available datasets suggest that deviations remain within ∼30%, which we consider acceptable for this homogeneous analysis. For galaxy groups, the median acceleration was estimated within the mean projected radius using the spatially resolved distribution of member galaxies. This characteristic acceleration, ⟨gbar⟩, served as the principal parameter for selecting low-acceleration subsamples throughout our study.

We modeled the BFJR in logarithm space as a linear relation: y = mx + b with y = log10(Mbar/M) and x = log10(σe/km s−1). We fitted the data using the BayesLineFit software (Lelli et al. 2019), which implements a Markov chain Monte Carlo (MCMC) method that takes errors on both variables as well as intrinsic scatter (σint) into account. The intrinsic scatter is assumed to be Gaussian and can be defined either in the vertical direction (along the y variable) or in the orthogonal one (perpendicular to the best-fit line); we explored both options.

3. Results: Acceleration-dependent relations

3.1. Baryonic Faber-Jackson relation

The BFJR is depicted for the full sample in Fig. 1 (left panel); the points are color-coded according to the internal acceleration, ⟨gbar⟩. The data points are all broadly correlated across eight orders of magnitude in mass, but high-acceleration systems are clearly shifted toward higher σe for a given baryonic mass.

Thumbnail: Fig. 1. Refer to the following caption and surrounding text. Fig. 1.

BFJR in galaxy groups, elliptical galaxies, and dwarf galaxies. Left: Total baryonic mass (Mbar) versus velocity dispersion within the effective radius (σe) for the full sample. The data points are color-coded by the internal median baryonic acceleration, ⟨gbar⟩. Middle: BFJR for the low-acceleration subsample only (⟨gbar⟩ < 0.6a0). In both panels, the dashed green line shows the MOND prediction in the low-acceleration regime, the solid black line is the best fit from the orthogonal MCMC, and the orange region is its 1σ credible interval. Right: Variation in the fitted parameters with the acceleration cutoff value, ⟨gbar⟩/a0: slope (m; top) and intercept (b; bottom). Orange diamonds are the result from orthogonal MCMC fitting, blue circles from vertical MCMC fitting. The number of objects at each cutoff is listed in the upper-right panel. The horizontal dashed lines mark the theoretical expectations from MOND modified gravity theories (m = 4, b = 3.1).

This visual trend suggests that the properties of the BFJR are not universal but instead depend on ⟨gbar⟩. To test this hypothesis, we fitted various subsamples for which ⟨gbar⟩ < Xa0, where X ranges from 0.1 to 20. Figure 1 (right panel) illustrates how the fitted parameters (slope m and intercept b) vary as a function of the ⟨gbar⟩/a0 threshold. We do not show the intrinsic scatter because the various subsamples have different sizes and heterogeneous error estimates, so comparing the intrinsic scatters can be misleading. As we restricted the sample to progressively lower accelerations, the slope (m) steadily converges toward ∼4 and the intercept (b) becomes stable. This convergence suggests that low-acceleration and high-acceleration systems follow different BFJRs. This is confirmed in Appendix B, in which we perform the same exercise for subsamples for which ⟨gbar⟩ > Xa0.

We selected ⟨gbar⟩ < 0.6 a0 as a practical threshold to compromise between sample statistics (153 systems) and convergence of the best-fit results. The resulting BFJR is shown in Fig. 1 (middle panel). The best-fit relation is

log 10 ( M bar M ) = ( 4.19 ± 0.10 ) log 10 ( σ e km s 1 ) + ( 2 . 55 0.16 + 0.16 ) , Mathematical equation: $$ \begin{aligned} \log _{10}\left(\frac{M_{\mathrm{bar} }}{M_\odot }\right) = (4.19 \pm 0.10) \log _{10}\left(\frac{\sigma _{\rm e}}{\mathrm{{km\,s}^{-1}}}\right) + (2.55^{+0.16}_{-0.16}), \end{aligned} $$(2)

The corner plots of both vertical and orthogonal MCMC analyses are shown in Appendix C. We explore possible residual correlations in Appendix D. Importantly, the residuals show no correlation with effective radius, indicating that it is not possible to decrease the observed scatter with a third structural variable, contrary to the case of the BFJR of high-acceleration systems (see Fig. D.2); see Sect. 4 for further discussion.

3.2. Fundamental plane

Figure 2 shows the baryonic FP considering Mbar versus the Newtonian virial estimator 5Reσe2/G. The factor of 5 is adopted following Cappellari et al. (2013b). If no DM is present, the data should lie on the line of unity according to the Newtonian virial theorem. Massive ellipticals with high ⟨gbar⟩ follow the Newtonian expectation, while systems with low ⟨gbar⟩ systematically deviate at both low masses (dwarf galaxies) and high masses (galaxy groups). Performing a Bayesian orthogonal regression on the subsample with ⟨gbar⟩ > 6 a0, we find

log 10 ( M bar M ) = ( 0.99 ± 0.01 ) log 10 ( 5 R e σ e 2 G M ) + ( 0.04 ± 0.15 ) . Mathematical equation: $$ \begin{aligned} \log _{10}\left(\frac{M_{\mathrm{bar} }}{M_\odot }\right) = (0.99 \pm 0.01) \log _{10}\left(\frac{5R_{\mathrm{e} }\sigma ^2_{\mathrm{e} }}{GM_\odot }\right) + (0.04 \pm 0.15)\,. \end{aligned} $$(3)

Thumbnail: Fig. 2. Refer to the following caption and surrounding text. Fig. 2.

FP for pressure-supported systems, including galaxy groups, ellipticals, and dwarf galaxies. The x-axis shows the expected Newtonian dynamical mass, log10(5Reσe2/G), while the y-axis gives the observed baryonic mass, log10(Mbar). The symbols are color-coded by the median baryonic acceleration within the effective radius. The Newtonian expectations (dashed line) are followed only by high-acceleration systems with ⟨gbar⟩ > a0, while low-acceleration systems systematically depart from it. The inset presents the MCMC analysis for the subsample restricted to systems with ⟨gbar⟩ > 6a0.

4. Discussion

4.1. Consistency with the MOND paradigm

Our studies empirically corroborate three predictions of MOND:

  1. Acceleration dependence: High-acceleration systems (elliptical galaxies) follow the Newtonian FP with no need for DM, while low-acceleration systems (dwarf galaxies and galaxy groups) deviate from it and require large amounts of DM. Yet, low-acceleration systems define a linear BFJR with small scatter. This confirms the MOND prediction that a0 marks a transition scale below which the dynamical behavior of galaxies and galaxy systems fundamentally changes.

  2. Slope: Our measured slopes of 4.19 ± 0.10 (orthogonal fit) and 3.92 ± 0.10 (vertical fit) are consistent with the MOND prediction of m = 4 within the 95%(∼2σ) and 68%(∼1σ) confidence intervals, respectively.

  3. Normalization: In the specific cases of the nonrelativistic modified gravity theories AQUAL (Bekenstein & Milgrom 1984) and QUMOND (Milgrom 2010), the proportionality factor in Eq. (1) is exactly 9/4. Assuming orbital isotropy, we have σ 3 D = 3 σ los Mathematical equation: $ \sigma_\mathrm{{3D}} = \sqrt{3}\sigma_\mathrm{{los}} $ and the proportionality factor becomes 81/4. The best-fit intercept of the low-acceleration BFJR is statistically consistent with this predicted value.

The dwarf galaxies in our sample are either satellites of massive spirals (the Milky Way and Andromeda) or in galaxy clusters (Virgo and Fornax), so they may be affected by the MOND external field effect (EFE; Bekenstein & Milgrom 1984), the Newtonian external field experienced by these dwarfs is around ∼0.01 − 0.1a0, less than our acceleration cut (∼0.6a0). Similarly, galaxy groups may be affected by the EFE due to the large-scale structure of the Universe, which becomes relevant for ⟨gbar⟩ ≲ 0.01a0 (Chae et al. 2021; Kelleher & Lelli 2024). The fact that the low-acceleration BFJR agrees with the MOND prediction for isolated systems suggests that the EFE must play a secondary, subtle role. Interestingly, the residuals around the BFJR show a weak correlation with ⟨gbar⟩ (see Fig. D.1). This is qualitatively consistent with the EFE because the systems with the lowest internal accelerations could be more Newtonian than MONDian and so display a lower σe at fixed Mbar (or higher Mbar at fixed σe). Ideally, one would like to plot the BFJR residuals against ⟨gbar⟩/gext, where gext is the baryonic external field felt by each system. This requires a more accurate study of the environment of each dwarf galaxy and each galaxy group.

4.2. Implications for galaxy formation models

In a ΛCDM context, it is surprising that dwarf galaxies and galaxy groups lie on the same BFJR despite being totally different systems that are shaped by very different physical processes. The small scatter (∼0.11 dex) observed across ∼8 orders of magnitude in mass leaves little room for stochastic variation in these processes (Desmond & Wechsler 2017). Even if we focus only on galaxies, the tightness of the BFJR demands finely tuned feedback processes (e.g., from supernovae and active galactic nuclei) to regulate star formation and set a precise DM fraction as a function of mass. Future work could compare our results to state-of-the-art cosmological simulations, such as EAGLE (Crain et al. 2015), BAHAMAS (McCarthy et al. 2017), and IllustrisTNG (Nelson et al. 2019). In general, it is unclear why the formation and evolution of galaxy groups and dwarf galaxies in ΛCDM should conspire to resemble the a priori prediction of MOND.

5. Conclusions

By compiling a comprehensive sample of pressure-supported systems spanning about eight orders of magnitude in mass, we find that the properties of the BFJR and of the FP systematically change with the mean internal acceleration, ⟨gbar⟩. In the low-acceleration regime (⟨gbar⟩ < 0.6a0), the BFJR converges to a tight power law, M bar σ e 4.19 ± 0.10 Mathematical equation: $ {M_{\text{bar}}}\propto \sigma_{\mathrm{e}}^{4.19\pm0.10} $. In the high-acceleration regime (⟨gbar⟩ > 6a0), the FP offers a superior description of the data compared to the BFJR. These findings are in excellent agreement with the predictions of MOND, with the characteristic acceleration scale (a0) naturally accounting for the observed behavior.

Overall, the BFJR and FP emerge as fundamental scaling relations of galaxies and galaxy groups. Their tightness and acceleration dependence pose a significant fine-tuning challenge for models within the standard cosmological framework but provide a powerful empirical testbed for distinguishing between competing theories of gravity and galaxy formation.

Data availability

Datasets are available at the CDS via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/710/L39

Acknowledgments

We thank the referee for the feedback and suggestions. We also thank Moti Milgrom and Pradyumna Sadhu. YT acknowledges the Taiwan National Science and Technology Council (NSTC) grants 110-2112-M-008-015-MY3 and 114-2112-M-008-024-MY3. YT and KHC acknowledge the National Research Foundation of Korea (grant no. NRF-2022R1A2C1092306). MSP acknowledges funding via a Leibniz-Junior Research Group (project number J94/2020). SSM is supported in part by NASA ADAP grant 80NSSC19k0570 and also acknowledges support from NSF PHY-1911909. YD is supported by the Postdoctoral Research Abroad Program (PRAP) grant NSTC 114-2917-I-564-044 and NSTC 114-2124-M-008-003. EDT is supported by the European Research Council (ERC) under grant agreement No. 101040751. MHK was supported by NSTC grant 110-2112-M-008-015-MY3. CMK is supported by the Taiwan NSTC 114-2112-M-008-018.

References

  1. Banik, I., & Zhao, H. 2022, Symmetry, 14, 1331 [NASA ADS] [CrossRef] [Google Scholar]
  2. Bekenstein, J., & Milgrom, M. 1984, ApJ, 286, 7 [NASA ADS] [CrossRef] [Google Scholar]
  3. Bolton, A. S., Burles, S., Treu, T., Koopmans, L. V. E., & Moustakas, L. A. 2007, ApJ, 665, L105 [NASA ADS] [CrossRef] [Google Scholar]
  4. Bundy, K., Bershady, M. A., Law, D. R., et al. 2015, ApJ, 798, 7 [Google Scholar]
  5. Cappellari, M. 2013, ApJ, 778, L2 [NASA ADS] [CrossRef] [Google Scholar]
  6. Cappellari, M., Bacon, R., Bureau, M., et al. 2006, MNRAS, 366, 1126 [Google Scholar]
  7. Cappellari, M., Emsellem, E., Krajnović, D., et al. 2011, MNRAS, 413, 813 [Google Scholar]
  8. Cappellari, M., McDermid, R. M., Alatalo, K., et al. 2013a, MNRAS, 432, 1862 [NASA ADS] [CrossRef] [Google Scholar]
  9. Cappellari, M., Scott, N., Alatalo, K., et al. 2013b, MNRAS, 432, 1709 [Google Scholar]
  10. Chae, K.-H., Desmond, H., Lelli, F., McGaugh, S. S., & Schombert, J. M. 2021, ApJ, 921, 104 [NASA ADS] [CrossRef] [Google Scholar]
  11. Ciotti, L., Lanzoni, B., & Renzini, A. 1996, MNRAS, 282, 1 [Google Scholar]
  12. Crain, R. A., Schaye, J., Bower, R. G., et al. 2015, MNRAS, 450, 1937 [NASA ADS] [CrossRef] [Google Scholar]
  13. Desmond, H., & Wechsler, R. H. 2017, MNRAS, 465, 820 [Google Scholar]
  14. Djorgovski, S., & Davis, M. 1987, ApJ, 313, 59 [Google Scholar]
  15. Domínguez Sánchez, H., Margalef, B., Bernardi, M., & Huertas-Company, M. 2022, MNRAS, 509, 4024 [Google Scholar]
  16. Dressler, A., Lynden-Bell, D., Burstein, D., et al. 1987, ApJ, 313, 42 [Google Scholar]
  17. Duann, Y., Tian, Y., & Ko, C.-M. 2023, RAS Tech. Instrum., 2, 649 [Google Scholar]
  18. Eftekhari, F. S., Peletier, R. F., Scott, N., et al. 2022, MNRAS, 517, 4714 [NASA ADS] [CrossRef] [Google Scholar]
  19. Faber, S. M., & Jackson, R. E. 1976, ApJ, 204, 668 [Google Scholar]
  20. Famaey, B., & McGaugh, S. S. 2012, Liv. Rev. Rel., 15, 10 [Google Scholar]
  21. Kauffmann, G., Heckman, T. M., White, S. D. M., et al. 2003, MNRAS, 341, 33 [Google Scholar]
  22. Kelleher, R., & Lelli, F. 2024, A&A, 688, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Krajnović, D., Emsellem, E., Cappellari, M., et al. 2011, MNRAS, 414, 2923 [Google Scholar]
  24. Krajnović, D., Alatalo, K., Blitz, L., et al. 2013, MNRAS, 432, 1768 [Google Scholar]
  25. Kroupa, P. 2001, MNRAS, 322, 231 [NASA ADS] [CrossRef] [Google Scholar]
  26. Lelli, F. 2022, Nat. Astron., 6, 35 [NASA ADS] [CrossRef] [Google Scholar]
  27. Lelli, F., McGaugh, S. S., & Schombert, J. M. 2016, ApJ, 816, L14 [Google Scholar]
  28. Lelli, F., McGaugh, S. S., Schombert, J. M., & Pawlowski, M. S. 2017, ApJ, 836, 152 [Google Scholar]
  29. Lelli, F., McGaugh, S. S., Schombert, J. M., Desmond, H., & Katz, H. 2019, MNRAS, 484, 3267 [NASA ADS] [CrossRef] [Google Scholar]
  30. Makarov, D., & Karachentsev, I. 2011, MNRAS, 412, 2498 [NASA ADS] [CrossRef] [Google Scholar]
  31. McCarthy, I. G., Schaye, J., Bird, S., & Le Brun, A. M. C. 2017, MNRAS, 465, 2936 [Google Scholar]
  32. McGaugh, S. S., Schombert, J. M., Bothun, G. D., & de Blok, W. J. G. 2000, ApJ, 533, L99 [Google Scholar]
  33. McGaugh, S. S., Mistele, T., Duey, F., et al. 2026, ApJ, 1001, 65 [Google Scholar]
  34. Milgrom, M. 1983, ApJ, 270, 365 [Google Scholar]
  35. Milgrom, M. 1984, ApJ, 287, 571 [CrossRef] [Google Scholar]
  36. Milgrom, M. 2010, MNRAS, 403, 886 [NASA ADS] [CrossRef] [Google Scholar]
  37. Milgrom, M. 2014a, Phys. Rev. D, 89, 024016 [NASA ADS] [CrossRef] [Google Scholar]
  38. Milgrom, M. 2014b, MNRAS, 437, 2531 [NASA ADS] [CrossRef] [Google Scholar]
  39. Milgrom, M. 2019, Phys. Rev. D, 99, 044041 [NASA ADS] [CrossRef] [Google Scholar]
  40. Milgrom, M. 2025, arXiv e-prints [arXiv:2510.16520] [Google Scholar]
  41. Nelson, D., Springel, V., Pillepich, A., et al. 2019, Comput. Astrophys. Cosmol., 6, 2 [Google Scholar]
  42. Pillepich, A., Nelson, D., Hernquist, L., et al. 2018, MNRAS, 475, 648 [Google Scholar]
  43. Sadhu, P., & Tian, Y. 2024, MNRAS, 528, 5612 [Google Scholar]
  44. Sanders, R. H. 2010, MNRAS, 407, 1128 [CrossRef] [Google Scholar]
  45. Scott, N., Eftekhari, F. S., Peletier, R. F., et al. 2020, MNRAS, 497, 1571 [Google Scholar]
  46. Serra, P., Oosterloo, T., Morganti, R., et al. 2012, MNRAS, 422, 1835 [Google Scholar]
  47. Taylor, E. N., Hopkins, A. M., Baldry, I. K., et al. 2011, MNRAS, 418, 1587 [Google Scholar]
  48. Tian, Y., Cheng, H., McGaugh, S. S., Ko, C.-M., & Hsu, Y.-H. 2021, ApJ, 917, L24 [CrossRef] [Google Scholar]
  49. Toloba, E., Guhathakurta, P., Peletier, R. F., et al. 2014, ApJS, 215, 17 [Google Scholar]
  50. Vazdekis, A., Sánchez-Blázquez, P., Falcón-Barroso, J., et al. 2010, MNRAS, 404, 1639 [NASA ADS] [Google Scholar]

Appendix A: Samples

The data analyzed in this paper are drawn from six public catalogs. A description of each catalog is provided below. Figure A.1 shows the baryonic mass against the effective radius (or the equivalent mean radius for galaxy groups). Our sample spans ∼8 dex in baryonic mass (Mbar ≃ 105 − 1013 M) and ∼4 dex in characteristic size (from ∼100 pc to ∼1 Mpc). The complete data table is available at the CDS.

Thumbnail: Fig. A.1. Refer to the following caption and surrounding text. Fig. A.1.

Relation between the baryonic mass and characteristic size (effective radius for galaxies, mean radius for galaxy groups) of our sample. Different symbols represent different datasets: galaxy groups (large diamonds), MaNGA ellipticals (circles), ATLAS3D ellipticals (squares), Virgo dwarfs (upward triangles), Fornax dwarfs (downward triangles), and Local Group dwarfs (small diamonds). The color coding represents the median baryonic gravitational acceleration, log10(⟨gbar⟩), within the effective radius (or the mean radius for galaxy groups).

Galaxy groups: We used data for 13 galaxy groups from Sadhu & Tian (2024) and 50 groups from Milgrom (2019), with original measurements primarily based on Makarov & Karachentsev (2011). The baryonic mass is given by the stellar mass of the member galaxies plus the hot gas mass from X-ray observations. Stellar masses were estimated from the K-band absolute magnitude (MK) of each member galaxy using the relation from Cappellari (2013), as adopted in Sadhu & Tian (2024):

log 10 M star = 10.58 0.44 ( M K + 23 ) . Mathematical equation: $$ \begin{aligned} \log _{10} M_{\mathrm{star} }= 10.58 - 0.44(M_{K} + 23)\,. \end{aligned} $$(A.1)

This relation is consistent with a Kroupa initial mass function (IMF). According to these measurements, the stellar mass dominates over the hot gas mass. In principle, there may be a warm-hot intergalactic medium (with temperatures of 105 − 106 K) that is not observed in the X-rays. Given the considerable uncertainties on the amount of such warm-hot gas component, we do not include it in our baryonic mass estimate. Velocity dispersions are calculated using a robust biweight estimator applied to the line-of-sight velocities of group members within the mean radius, defined as the average projected radius of all member galaxies.

In the ΛCDM context, there should be missing baryons in galaxy groups, possibly in the form of warm-hot gas (e.g., McGaugh et al. 2026). This does not need to hold in other paradigms, such as MOND. In this work, we consider only the baryonic mass that is directly observed: stars and X-ray gas. A substantial amount of missing baryons in galaxy groups will lead to a steeper slope and smaller intercept of the low-acceleration BFJR. It will also move galaxy groups close to the Newtonian expectation in the baryonic FP.

MaNGA elliptical galaxies: We selected galaxies from the Sloan Digital Sky Survey IV (SDSS-IV) MaNGA survey (Bundy et al. 2015). MaNGA provides two-dimensional spectroscopic maps for nearly 10,000 galaxies, enabling precise measurements of stellar kinematics. From this dataset, 2,632 elliptical galaxies were identified based on morphological classifications using a convolutional neural network, as presented in Domínguez Sánchez et al. (2022). To ensure that the velocity dispersion is a faithful tracer of the equilibrium gravitational potential, we excluded galaxies exhibiting ascending (5.6%) or irregular (2.1%) velocity dispersion profiles, as classified by Duann et al. (2023). Only galaxies with declining or flat velocity dispersion profiles were included in this study. Finally, we restrict to galaxies with velocity dispersion measurements out to at least one effective radius (Re), so we can calculate σe (the average line-of-sight velocity dispersion within Re). Our final sample comprises 1218 elliptical galaxies. The effective radius and Sérsic index are provided by the MaNGA PyMorph photometric Value Added Catalogue (MPP-VAC-DR17) in Domínguez Sánchez et al. (2022).

Stellar masses are derived from fitting the spectral energy distribution assuming a Kroupa IMF (Kauffmann et al. 2003; Kroupa 2001). The hot gas mass in elliptical galaxies was estimated using the scaling relation from Chae et al. (2021):

log ( M g , hot / M ) = 1.47 log ( M star / M ) 5.414 . Mathematical equation: $$ \begin{aligned} \log (M_{\rm g,hot}/M_\odot ) = 1.47\log (M_{\mathrm{star} }/M_\odot )-5.414\,. \end{aligned} $$(A.2)

Given that the cold gas mass of ellipticals may be (at most) a few percent of the stellar mass (Serra et al. 2012), we assume Mbar ≈ Mstar + Mg, hot.

ATLAS3D elliptical galaxies: This sample comprises 260 local early-type galaxies from the ATLAS3D survey (Cappellari et al. 2011). We restrict the sample to 68 elliptical galaxies using the morphological classifications from Krajnović et al. (2011). For consistency with other samples, we limit our analysis to 26 galaxies with data extending out to at least one effective radius, enabling measurements of σe. Stellar masses and kinematic data are taken from Cappellari et al. (2013b), while the Sérsic indices and effective radii are adopted from Krajnović et al. (2013).

Stellar masses are estimated using the mass-to-light ratio (M/L) at the effective radius, as determined with the Salpeter IMF in Cappellari et al. (2013a). To ensure consistency with other samples, we convert these values to a Kroupa IMF using Mstar = MSalp/1.6. The hot gas mass is estimated using Eq. (A.2), yielding a median value of 25% of Mstar in our sample. Cold gas masses, as measured by Serra et al. (2012), contribute only a few percent relative to Mstar. The total baryonic mass is thus defined as Mbar = Mstar + Mgas + Mg, hot.

Virgo dwarfs: We selected 34 out of 39 Virgo dwarf galaxies from Toloba et al. (2014), requiring that the velocity dispersion within the effective radius exceed the rotational velocity. Stellar masses were estimated assuming a constant mass-to-light ratio of M/L = 0.73 in the H band at one effective radius by adopting the Kroupa IMF (Vazdekis et al. 2010). The Sérsic index ranges from 1.0 to 2.2, as reported in Toloba et al. (2014); for a few galaxies without measured indices, we adopted n = 1 for the calculation of baryonic acceleration. Given their cluster environment, these dwarfs contain minimal gas, so we assumed Mbar ≃ Mstar.

Fornax dwarfs: Data for dwarf galaxies in Fornax are drawn from the SAMI-Fornax Dwarf Galaxy Survey (Scott et al. 2020; Eftekhari et al. 2022). We selected 31 galaxies morphologically classified as dwarf ellipticals. The SAMI integral-field spectrograph provides kinematic data extending beyond the half-light radius. The velocity dispersion is consistently measured within one effective radius. The Sérsic parameters used for the baryonic acceleration calculations are reported in Table 1 of Eftekhari et al. (2022). Stellar masses were estimated by Taylor et al. (2011) and Eftekhari et al. (2022) using the Chabrier IMF, based on r-band measurements and g − i and r − i colors:

log ( M star / M ) e = 1.15 + 0.75 ( g i ) 0.4 M r , e + 0.4 ( r i ) . Mathematical equation: $$ \begin{aligned} \log (M_{\mathrm{star} }/M_\odot )_{e} = 1.15 + 0.75(g-i) - 0.4M_{r,e} + 0.4(r-i). \end{aligned} $$(A.3)

The Chabrier IMF is virtually equivalent to the Kroupa IMF assumed for the other samples. Given their cluster environment, these dwarfs contain very little gas, so we assume Mbar ≈ Mstar.

Local Group dwarf spheroidals: To extend the BFJR to the lowest mass regime, we include dwarf spheroidal (dSph) satellite galaxies of the Milky Way and Andromeda. Data are taken from the compilation by Lelli et al. (2017). We select 28 dSph galaxies with luminosities greater than 105M (so we do not consider the so-called "ultra-faint dwarfs" that have substantially more uncertain data) and with velocity dispersions measured from more than 20 member stars. Stellar masses are estimated assuming a mass-to-light ratio of M/L = 2 in the V band (for a Kroupa IMF), and a Sérsic index of n = 1 is adopted for all systems. As these satellites are almost entirely devoid of gas, their baryonic mass is effectively equal to their stellar mass.

Finally, we note that compiling data across eight orders of magnitude in mass inherently requires combining systems with different observational constraints. Consequently, the adopted velocity dispersions are not strictly identical observables across all subsamples. For individual elliptical and dwarf galaxies, we use the stellar line-of-sight velocity dispersion integrated within one effective radius (σe). In contrast, for galaxy groups, we use the velocity dispersion derived from the discrete line-of-sight velocities of member galaxies. While our homogenization is intended to provide a characteristic measure of pressure support across vastly different classes of systems, these quantities represent physically distinct tracers. These systematic differences between subsamples should be kept in mind when interpreting the absolute parameters of the scaling relations, as they naturally contribute additional intrinsic scatter to the overarching empirical trend.

Appendix B: The baryonic Faber-Jackson relation for high-acceleration subsamples

To characterize the BFJR in the high-acceleration regime, we fitted a linear relation to subsamples with ⟨gbar⟩ > Xa0, where X is a number from 1 to 50. Similarly to Sect. 3, we performed both orthogonal and vertical MCMC fits, explicitly accounting for measurement uncertainties in log10(Mbar) and log10(σe), as well as intrinsic scatter. This methodology provides robust estimates of the slope and intercept even when measurement errors affect both axes.

The behavior of the BFJR under increasingly stringent acceleration thresholds is summarized in Fig. B.1. The left panel reproduces the full sample, while the middle panel shows the high-acceleration subsample (⟨gbar⟩ > 6a0) together with the best-fitting orthogonal MCMC relation, y = 4.32x + 1.36. The right panel compares orthogonal and vertical MCMC fits for subsamples selected above increasing thresholds of ⟨gbar⟩/a0, complementary to the presentation in the right panel of Fig. 1. The two fitting approaches yield systematically different results in this regime, indicating that the inferred BFJR parameters depend on the adopted fitting scheme. This is common for linear relations with steep slopes (e.g., Lelli et al. 2019; Tian et al. 2021).

Thumbnail: Fig. B.1. Refer to the following caption and surrounding text. Fig. B.1.

Left: BFJR for the full sample. Middle: BFJR for the high-acceleration subsample, (⟨gbar⟩ > 6a0). In both panels, the black line shows the best-fitting relation from the orthogonal MCMC fit, and the orange region denotes the 1σ credible interval. Right: Variation in the fitted slope and intercept as a function of the acceleration threshold ⟨gbar⟩/a0.

To revisit the existence of the FP, Fig. D.2 presents both the orthogonal and vertical residuals of the high-acceleration subsample as a function of internal acceleration and effective radius. It is evident by eye that the residuals correlate with Re (as expected due to the FP) and more weakly with ⟨gbar⟩ (also expected because the internal baryonic acceleration depends on the baryonic surface density). Indeed, the Pearson’s test gives a correlation coefficient r ≃ 0.1 with p ≃ 2 × 10−4 for the orthogonal residuals against Re, and r ≃ 0.58 with p ≪ 10−5 for the vertical ones. This analysis confirms that an additional structural parameter can reduce the scatter around the BFJR of the high-acceleration regime, contrary to the case of the BFJR for the low-acceleration sample (see Fig. D.1).

Appendix C: Posterior probability distributions

Figure C.1 shows the “corner plots” for the orthogonal and vertical MCMC fits to the low-acceleration subsample (gbar < 0.6a0). The 1D posterior probability distributions are single-peaked and close to a Gaussian function, so the best-fit parameters and their uncertainties are well defined. As always happens in linear fits, the slope and intercept are somewhat degenerate.

Thumbnail: Fig. C.1. Refer to the following caption and surrounding text. Fig. C.1.

Posterior probability distributions of BFJR fit parameters when selecting galaxies with gbar < 0.6a0. Top: Results from the orthogonal MCMC fit, showing the marginalized and joint posterior distributions for the slope, intercept, and intrinsic scatter of the BFJR. Bottom: Results from the vertical MCMC fit using the same data. The best-fit values and 1σ uncertainties are indicated by red lines and annotations. Contours correspond to the 1σ and 2σ credible regions. The comparison demonstrates the impact of the fitting method on the derived BFJR parameters and their uncertainties.

Appendix D: Residuals of the low-acceleration baryonic Faber-Jackson relation

To check for potential secondary correlations, Fig. D.1 shows the residuals of the low-acceleration BFJR against different properties of the systems. The orthogonal residuals show no systematic correlation with the characteristic size of the system (the effective radius of galaxies or the mean radius of galaxy groups). This indicates that the scatter around the BFJR of low-acceleration systems cannot be decreased by adding a third structural variable, contrary to the case of the FP of high-acceleration systems (see Appendix B). Indeed, the Pearson’s test gives a negligible correlation coefficient r ≃ 0.1 with p ≃ 0.14. The same test gives a potentially significant correlation between vertical residuals and Re (r ≃ 0.3 and p ≃ 5 × 10−5) but this may be a shortcoming of the vertical fit, which is known to underestimate the best-fit slope of steep linear relations (Lelli et al. 2019).

Thumbnail: Fig. D.1. Refer to the following caption and surrounding text. Fig. D.1.

BFJR residuals versus two different physical quantities for the low-acceleration subsample (⟨gbar⟩ < 0.6a0). Top panels: Orthogonal residuals versus the logarithm of the mean internal acceleration, log10gbar⟩ (left) and the logarithm of the effective radius, log10(Re) (right). Bottom panels: Same as the top panels but for the vertical residuals.

Intriguingly, we find a statistically significant, albeit weak, anti-correlation with the internal acceleration ⟨gbar⟩ for both the orthogonal fit (r ≃ −0.2 with p ≃ 5 × 10−3) and for the vertical fit (r ≈ −0.3 with p ≃ 6 × 10−4), which may possibly be driven by the EFE in MOND.

Thumbnail: Fig. D.2. Refer to the following caption and surrounding text. Fig. D.2.

BFJR residuals versus two different physical quantities for the high-acceleration subsample (⟨gbar⟩ > 6a0). Top panels: Orthogonal residuals versus the logarithm of the internal acceleration, log10(⟨gbar⟩) (left) and the logarithm of the effective radius, log10(Re) (right). Bottom panels: Same as the top panels but for the vertical residuals.

All Figures

Thumbnail: Fig. 1. Refer to the following caption and surrounding text. Fig. 1.

BFJR in galaxy groups, elliptical galaxies, and dwarf galaxies. Left: Total baryonic mass (Mbar) versus velocity dispersion within the effective radius (σe) for the full sample. The data points are color-coded by the internal median baryonic acceleration, ⟨gbar⟩. Middle: BFJR for the low-acceleration subsample only (⟨gbar⟩ < 0.6a0). In both panels, the dashed green line shows the MOND prediction in the low-acceleration regime, the solid black line is the best fit from the orthogonal MCMC, and the orange region is its 1σ credible interval. Right: Variation in the fitted parameters with the acceleration cutoff value, ⟨gbar⟩/a0: slope (m; top) and intercept (b; bottom). Orange diamonds are the result from orthogonal MCMC fitting, blue circles from vertical MCMC fitting. The number of objects at each cutoff is listed in the upper-right panel. The horizontal dashed lines mark the theoretical expectations from MOND modified gravity theories (m = 4, b = 3.1).

In the text
Thumbnail: Fig. 2. Refer to the following caption and surrounding text. Fig. 2.

FP for pressure-supported systems, including galaxy groups, ellipticals, and dwarf galaxies. The x-axis shows the expected Newtonian dynamical mass, log10(5Reσe2/G), while the y-axis gives the observed baryonic mass, log10(Mbar). The symbols are color-coded by the median baryonic acceleration within the effective radius. The Newtonian expectations (dashed line) are followed only by high-acceleration systems with ⟨gbar⟩ > a0, while low-acceleration systems systematically depart from it. The inset presents the MCMC analysis for the subsample restricted to systems with ⟨gbar⟩ > 6a0.

In the text
Thumbnail: Fig. A.1. Refer to the following caption and surrounding text. Fig. A.1.

Relation between the baryonic mass and characteristic size (effective radius for galaxies, mean radius for galaxy groups) of our sample. Different symbols represent different datasets: galaxy groups (large diamonds), MaNGA ellipticals (circles), ATLAS3D ellipticals (squares), Virgo dwarfs (upward triangles), Fornax dwarfs (downward triangles), and Local Group dwarfs (small diamonds). The color coding represents the median baryonic gravitational acceleration, log10(⟨gbar⟩), within the effective radius (or the mean radius for galaxy groups).

In the text
Thumbnail: Fig. B.1. Refer to the following caption and surrounding text. Fig. B.1.

Left: BFJR for the full sample. Middle: BFJR for the high-acceleration subsample, (⟨gbar⟩ > 6a0). In both panels, the black line shows the best-fitting relation from the orthogonal MCMC fit, and the orange region denotes the 1σ credible interval. Right: Variation in the fitted slope and intercept as a function of the acceleration threshold ⟨gbar⟩/a0.

In the text
Thumbnail: Fig. C.1. Refer to the following caption and surrounding text. Fig. C.1.

Posterior probability distributions of BFJR fit parameters when selecting galaxies with gbar < 0.6a0. Top: Results from the orthogonal MCMC fit, showing the marginalized and joint posterior distributions for the slope, intercept, and intrinsic scatter of the BFJR. Bottom: Results from the vertical MCMC fit using the same data. The best-fit values and 1σ uncertainties are indicated by red lines and annotations. Contours correspond to the 1σ and 2σ credible regions. The comparison demonstrates the impact of the fitting method on the derived BFJR parameters and their uncertainties.

In the text
Thumbnail: Fig. D.1. Refer to the following caption and surrounding text. Fig. D.1.

BFJR residuals versus two different physical quantities for the low-acceleration subsample (⟨gbar⟩ < 0.6a0). Top panels: Orthogonal residuals versus the logarithm of the mean internal acceleration, log10gbar⟩ (left) and the logarithm of the effective radius, log10(Re) (right). Bottom panels: Same as the top panels but for the vertical residuals.

In the text
Thumbnail: Fig. D.2. Refer to the following caption and surrounding text. Fig. D.2.

BFJR residuals versus two different physical quantities for the high-acceleration subsample (⟨gbar⟩ > 6a0). Top panels: Orthogonal residuals versus the logarithm of the internal acceleration, log10(⟨gbar⟩) (left) and the logarithm of the effective radius, log10(Re) (right). Bottom panels: Same as the top panels but for the vertical residuals.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.