# Compactness peaks: An astrophysical interpretation of the mass distribution of merging binary black holes

Shanika Galaudage, Astrid Lamberts

► To cite this version:

Shanika Galaudage, Astrid Lamberts. Compactness peaks: An astrophysical interpretation of the mass distribution of merging binary black holes. Astronomy & Astrophysics - A&A, 2025, 694, pp.A186. 10.1051/0004-6361/202451654. [hal-04668246](https://hal.science/hal-04668246v1)

HAL Id: hal-04668246

https://hal.science/hal-04668246v1

Submitted on 20 Mar 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Compactness peaks: An astrophysical interpretation of the mass distribution of merging binary black holes

Shanika Galaudage<fnref n="1" />,<fnref n="2" />, * and Astrid Lamberts<fnref n="1" />,<fnref n="2" />, †

Received 25 July 2024 / Accepted 22 October 2024

## ABSTRACT

With the growing number of detections of binary black hole (BBH) mergers, we are beginning to probe structure in the distribution of mass. A recent study proposes that the isolated binary evolution of stripped stars naturally gives rise to the peaks at $M \sim 8 M_\odot$ and $14 M_\odot$ in the chirp-mass distribution and explains the dearth of black holes (BHs) in the mass range of $M \approx 10-12 M_\odot$. The gap in chirp mass results from an apparent gap in the component-mass distribution within $m_1, m_2 \approx 10-15 M_\odot$ and the specific pairing of these BHs. This component-mass gap results from variation in the core compactness of the progenitor, where a drop in compactness as a function of carbon-oxygen core mass means that BHs are no longer formed from core collapse. We develop a population model motivated by this scenario to probe the structure of the component-mass distribution of two populations of BBHs: one population consisting of two peak components, representing BHs formed in the compactness peaks, and another population with a power-law component to account for any polluting events, that is, binaries that may have formed from different channels (e.g. dynamical). We perform hierarchical Bayesian inference to analyse the events from the third gravitational-wave transient catalogue (GWTC-3) with our population model. We find that there is a preference for the lower-mass peak to drop off sharply at $\sim 11 M_\odot$ and the upper mass peak to turn on at $\sim 13 M_\odot$, in line with predictions in the literature. However, we find no clear evidence for a gap. We also find mild support for a scenario where the two populations have different spin distributions. In addition to these population results, we highlight observed events of interest that differ from the expected population distribution of compact objects formed from stripped stars.

**Key words.** gravitational waves – binaries: close – stars: black holes – stars: evolution

## 1. Introduction

Gravitational-wave astronomy is a rapidly growing field with approximately 90 detections (Abbott et al. 2023a) reported by the LIGO-Virgo-KAGRA collaboration (LVK). This catalogue of events was produced with data from the Advanced LIGO (Aasi et al. 2015) and Advanced Virgo (Acernese et al. 2015) gravitational-wave observatories. The majority of these detections are of merging binary black holes (BBHs). With population studies of these sources, the structure of the mass distribution of BBHs is slowly being revealed. This structure has implications for stellar physics processes, such as the explosion mechanisms of supernovae and mass-transfer processes in binaries (see Mapelli 2020; Mandel & Farmer 2022; Mandel & Broekgaarden 2022; Spera et al. 2022, for reviews).

The majority of the focus so far in terms of the astrophysical interpretation of these merging binaries has been in relation to the component masses of the compact binaries, with statistically significant peaks at $\approx 10 M_\odot$ and $\approx 35 M_\odot$ in the primary mass distribution (e.g. Talbot & Thrane 2018; Sadiq et al. 2022; Edelman et al. 2022; Wong & Cranmer 2022; Farah et al. 2023; Abbott et al. 2023b). The former peak may be associated with the binary mass-transfer processes (e.g. van Son et al. 2022), while the latter could be related to stellar physics, such as that regarding pulsational pair-instability supernovae (PPISNe; Heger & Woosley 2002; Croon & Sakstein 2023; Rahman et al. 2022), but there are suggestions that the fraction of events in this latter peak is too high for them to be PPISNe (e.g.

Stevenson et al. 2019), that the location of the peak due to PPISNe should be higher (e.g. Farmer et al. 2020; Farag et al. 2022; Hendriks et al. 2023; Golomb et al. 2024), and that the peak is populated by merging binaries in globular clusters (Antonini et al. 2023). We note that there is also a possible peak at $\sim 16-20 M_\odot$ and an under-density at $\sim 14 M_\odot$ (Tiwari & Fairhurst 2021; Tiwari 2022, 2023; Abbott et al. 2023b), although these features are not statistically significant (e.g. Farah et al. 2023). Recent studies have, however, started to look beyond the component-mass distributions to investigate and interpret the chirp-mass distribution in more detail (e.g. Abbott et al. 2023b; Tiwari & Fairhurst 2021; Tiwari 2022, 2023; Schneider et al. 2023).

The chirp mass (Finn & Chernoff 1993; Poisson & Will 1995; Blanchet et al. 1995) determines the leading-order orbital evolution of a binary system from gravitational-wave emission and is therefore more accurately measured than the component masses,

$$M = \frac{(m_1 m_2)^{3/5}}{(m_1 + m_2)^{1/5}} \quad (1)$$

It is a combination of the individual component masses $m_1$ and $m_2$ and is directly related to gravitational-wave frequency and frequency evolution. The chirp-mass distribution of BBH mergers shows statistically significant peaks at 8, 14, and 28 $M_\odot$, with an apparent lack of BBH mergers with $M \approx 10-12 M_\odot$ (e.g. Abbott et al. 2023b; Tiwari 2023), which is hereafter referred to as the ‘chirp-mass gap’.

<footnote n="1">¹ Laboratoire Lagrange, Université Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Bd de l'Observatoire, 06300 Nice, France</footnote>

<footnote n="2">² Laboratoire Artemis, Université Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Bd de l'Observatoire, 06300 Nice, France</footnote>

<footnote>* Corresponding authors; shanika.galaudage@oca.eu</footnote>
Schneider et al. (2023) recently proposed that merging BBHs formed via isolated binary evolution, where the progenitors lose their hydrogen-rich envelopes through mass-transfer processes, naturally give rise to the peaks at $M \sim 8 M_\odot$ and $14 M_\odot$ and a dearth of BHs within the range of $M \approx 10-12 M_\odot$. These authors find that the bimodality in the BBH chirp-mass distribution is due to the pairing of ‘compactness peaks’ in the component-mass distribution, which is linked to stellar physics and supernova mechanisms.

The explodability of stars is dependent on the compactness of the stellar core, $\xi = (M/M_\odot)/(R(M)/1000 \text{ km})$, where M is the mass of the core, and $R(M)$ is the radius at the mass coordinate M, where $M = 2.5 M_\odot$ is used in Schneider et al. (2023). The compactness is generally considered to not steadily increase with core mass, but rather to have sharp peaks and dips due to variations in burning stages in the core (e.g. O'Connor & Ott 2011; Sukhbold & Woosley 2014; Patton & Sukhbold 2020; Chieffi & Limongi 2020; Schneider et al. 2021, 2023; Burrows et al. 2023, 2024); however, there are studies in which a monotonic increase in compactness is considered (Mapelli et al. 2020). Considering the binary stripped stars in the Schneider et al. (2023) models, for 10% solar metallicity, the compactness level peaks due to neutrino-dominated burning of carbon-oxygen (CO) at core masses ($M_{CO}$) of $7.5 M_\odot$, producing BHs from core collapse supernova at $\approx 9-10 M_\odot$; these are the low-mass BHs, which we refer to as BH$_L$, which is the same notation as that used in Schneider et al. (2023). Above this core mass point, the compactness in the core drops, and this is a region where neutron stars are formed from core collapse rather than BHs (see Figure 3 in Schneider et al. 2023). At around $M_{CO} \approx 13 M_\odot$, burning becomes neutrino dominated in the oxygen-neon (ONe) stage as well, producing BHs of greater than $15 M_\odot$ from core collapse supernova; these are the high-mass BHs, which we refer to as BH$_H$. This variation in compactness level results in no BH masses being formed within $\approx 10-15 M_\odot$. These boundaries can shift due to uncertainties in the stellar physics, including that related to the metallicity-dependant mass-loss rate (Vink & de Koter 2005; Belczynski et al. 2010; Spera et al. 2015); convective core boundary mixing during core hydrogen and core helium burning (e.g. Schneider et al. 2023; Temaj et al. 2024), and $^{12}$C($\alpha$, $\gamma$)$^{16}$O nuclear reaction rates (e.g. Farmer et al. 2019, 2020; Costa et al. 2021)$^1$. However, the magnitude of the shifts at the lower-mass end are not expected to exceed $\approx 1 M_\odot$.

The peaks in chirp mass are produced by a mixture of pairings between BH$_L$ and BH$_H$. The peak below the chirp-mass gap is formed by BH$_L$ + BH$_L$ and the structures above the chirp-mass gap are formed by BH$_H$ + BH$_H$. Any BH$_L$ + BH$_H$ binaries are expected to fill this chirp-mass gap, but the pairing of BH$_L$ + BH$_H$ is less likely due to the fact that isolated binary evolution is more likely to produce systems with closer-to-equal mass ratios (Schneider et al. 2021). Assuming this model hypothesis is correct, the BH$_L$ + BH$_H$ population may be completely suppressed according to current observations (Abbott et al. 2023a).

In the present work, we develop a parametric population model motivated by the results of Schneider et al. in order to obtain an astrophysical interpretation of the mass distribution of BBH mergers. In Sect. 2 we provide a description of our population model and analysis method. In Sect. 3 we present our findings based on this model and compare them to results derived

from models commonly used by the LVK. Finally, in Sect. 4 we discuss the astrophysical implications of our results and avenues for future work.

## 2. Methodology

In this work, we explore whether the findings of Schneider et al. (2023) can be used to explain the population of observed BBH mergers. To this end, we first investigate the properties of the events individually in Sect. 2.1. We then develop an astrophysically informed population model in Sect. 2.2 to investigate the properties of the population and to investigate whether or not there is a gap-like feature in the component-mass distribution.

### 2.1. Identifying binary black hole systems of interest

We consider the BBH mergers analysed in previous population studies (e.g. Abbott et al. 2023b), which include 69 BBH mergers from the third LVK gravitational-wave transient catalogue (GWTC-3; Abbott et al. 2023a) detected with a FAR < 1yr$^{-1}$.

As an initial step, we find the positions of these events in the chirp-mass distribution and sort them according to whether they reside in the low- or high-mass region, as illustrated in the upper panel of Figure 1. We assume that those events in the low-mass region ($M < 11 M_\odot$) are BH$_L$ + BH$_L$ mergers and that those in the high-mass region ($M > 11 M_\odot$) are BH$_H$ + BH$_H$ mergers. We then show where these objects lie in the component-mass space. We find that from this first inspection, while there is a clear gap in the chirp-mass distribution, the component-mass distribution does not show a gap-like feature.

We cannot robustly determine whether or not there exists a gap-like structure from visual inspection of the posteriors alone, especially given that the uncertainties on the mass posteriors are wider than the width of the gap. A population study, as outlined in Sect. 2.2, is required to investigate this structure. However, it is interesting to consider whether there are binaries ‘polluting’ the $\approx 10-15 M_\odot$ region in component mass that were formed via other channels.

We identify the events where the median of the mass posterior distribution (for one or both components) lies above $10 M_\odot$ for events with $M < 11 M_\odot$ and below $15 M_\odot$ for events with $M > 11 M_\odot$. We find 20 events that pollute the component-mass gap and list some of the properties in Table 1, including the mass ratio ($q = m_2/m_1$) and the effective inspiral spin (Santamaria et al. 2010; Ajith et al. 2011),

$$ \chi_{\text{eff}} = \frac{\chi_1 \cos \theta_1 + q\chi_2 \cos \theta_2}{1 + q}, \qquad (2) $$

a mass-weighted spin parameter that is projected along the orbital angular momentum, where $\chi$ is the spin magnitude and $\theta$ is the spin orientation, and the effective precession spin (Hannam et al. 2014; Schmidt et al. 2015)

$$ \chi_p = \max \left[ \chi_1 \sin \theta_1, q \frac{(3q+4)}{(4q+3)} \chi_2 \sin \theta_2 \right], \qquad (3) $$

which is a measure of in-plane spin, and enables the parameterisation of the rate of relativistic precession of the orbital plane. These parameters may suggest the binary formed via a different channel. Mergers with an unequal mass ratio are unlikely to occur in isolated binary evolution (e.g. Kruckow et al. 2018; Giacobbo & Mapelli 2018; Broekgaarden et al. 2022; Iorio et al. 2023; Schneider et al. 2021). Negative values of $\chi_{\text{eff}}$ and large

<footnote n="1">¹ See Figure 6 in Schneider et al. (2023) for an illustration of how physical processes shift the edges of the $10-15 M_\odot$ gap in component mass.</footnote>
<figure><img src="image_2.png" /><figcaption>Fig. 1. Events in GWTC-3 with FAR < 1yr⁻¹ sorted according to the chirp-mass peak in which they lie, assuming the Schneider et al. (2023) predictions regarding their classification and the locations of the gaps. Top panel: Chirp-mass posteriors showing the BH<sub>L</sub> + BH<sub>L</sub> (blue) and the BH<sub>H</sub> + BH<sub>H</sub> (pink). The grey shaded region represents the chirp-mass gap region (M ≈ 10–12 M<sub>☉</sub>). Bottom panel: Component-mass posteriors showing the BH<sub>L</sub> (blue) and the BH<sub>H</sub> (pink) BHs. The grey shaded region represents the expected component-mass gap region (m ≈ 10–15 M<sub>☉</sub>). The vertical ticks illustrate the median of each posterior.</figcaption></figure>

values of $\chi_p$ are generally considered to be a signature of dynamical formation given the support for misaligned and in-plane spin orientations (Mandel & O'Shaughnessy 2010; Rodriguez et al. 2016; Zevin et al. 2017; Rodriguez et al. 2018; Doctor et al. 2020).

We also indicate whether or not the events have been considered in other studies to have formed via a channel that is different from isolated binary evolution. To obtain a visual sense of these potential ‘pollutants’ (i.e. events that may not have formed via isolated binary evolution), we also provide the $q$, $\chi_{eff}$, and $\chi_p$ distributions for these events in Figure 2 along with the median of the posterior for the other 49 events, which are illustrated by the black dots.

There are a few events in the component-mass gap (Table 1) that appear to differ from the typical properties of systems formed via isolated evolution with more extreme mass ratios and/or higher values of $\chi_p$ (e.g. GW190412, GW151226, GW200225) compared to the rest of the BBH mergers populating the gap (see the blue and pinks dots in Figure 2). However, there are also events not in this mass gap that have these differing properties (black dots in Figure 2). Given this information, we need to build a population model that can account for potential pollutants in the population in addition to finding the location and edges of the compactness peaks. Finding this structure will help determine whether or not there is a gap-like feature within ≈10–15 M<sub>☉</sub> in the component-mass distribution of BBH mergers.

## 2.2. Building an astrophysically motivated population model

We developed an astrophysically informed population model incorporating the features predicted in Schneider et al. (2023) for the component-mass distributions while accounting for ‘pol-

lutant’ events. We performed hierarchical Bayesian inference in order to measure the population hyper-parameters using gravitational-wave data from the LVK<fnref n="2" />. To this end, we used GWPopulation (Talbot et al. 2019), a hierarchical Bayesian inference package that employs Bilby (Ashton et al. 2019; Romero-Shaw et al. 2020). We use the nested sampler DYNESTY (Speagle 2020) for our analyses. Our dataset contains the 69 BBH observations from GWTC-3 (Abbott et al. 2023a) that were considered reliable for population analyses (events with a false-alarm rate of <1 yr⁻¹; Abbott et al. 2023b). We used the posterior samples from parameter estimation results used in Abbott et al. (2023b) along with the injection sets to account for selection effects (Essick & Farr 2022). There are three components to our model with various hyperparameters (Λ):

- A low-mass narrow Gaussian (BH<sub>L</sub> Peak): this models the compactness peak where BHs are formed from increased compactness of CO cores of stars from neutrino-dominated burning,

$$ \pi(m_1|\Lambda)_{\text{BH}_L} = \mathcal{N}(m_1|\mu_{\text{BH}_L}, \sigma_{\text{BH}_L}, m_{\text{min}}^{\text{BH}_L}, m_{\text{max}}^{\text{BH}_L}), \quad (4) $$

where $\mu_{\text{BH}_L}$ and $\sigma_{\text{BH}_L}$ are the mean and width of the Gaussian and $m_{\text{min}}^{\text{BH}_L}$ and $m_{\text{max}}^{\text{BH}_L}$ are its lower and upper limits.

- A high-mass broad Gaussian (BH<sub>H</sub> Peak): this models the BHs formed from the increased compactness from neutrino release in the ONe burning stage,

$$ \pi(m_1|\Lambda)_{\text{BH}_H} = \mathcal{N}(m_1|\mu_{\text{BH}_H}, \sigma_{\text{BH}_H}, m_{\text{min}}^{\text{BH}_H}, m_{\text{max}}^{\text{BH}_H}), \quad (5) $$

<footnote n="2">² See Thrane & Talbot (2019) for a detailed review of parameter estimation and hierarchical inference in gravitational wave astronomy; the methods used for our analysis are equivalent to those described in their Section 5.</footnote>
<table><caption>**Table 1.** Properties of events 'polluting' the component-mass gap ($m = 10-15 M_\odot$).</caption><thead><tr><th>Event; [REF]</th><th>q</th><th>χ<sub>eff</sub></th><th>χ<sub>p</sub></th><th>Suggestions for other formation scenarios</th></tr></thead><tbody><tr><td>GW190412_053044; [1]</td><td>0.28<sup>+0.12</sup><sub>-0.06</sub></td><td>0.25<sup>+0.08</sup><sub>-0.11</sub></td><td>0.31<sup>+0.2</sup><sub>-0.16</sub></td><td>e.g. Safarzadeh & Hotokezaka (2020); Gerosa et al. (2020); Liu & Lai (2021)</td></tr><tr><td>GW151226; [2]</td><td>0.57<sup>+0.37</sup><sub>-0.34</sub></td><td>0.18<sup>+0.2</sup><sub>-0.12</sub></td><td>0.5<sup>+0.39</sup><sub>-0.33</sub></td><td>e.g. Chatterjee et al. (2017)</td></tr><tr><td>GW190828_065509; [3]</td><td>0.43<sup>+0.38</sup><sub>-0.16</sub></td><td>0.08<sup>+0.15</sup><sub>-0.16</sub></td><td>0.29<sup>+0.41</sup><sub>-0.21</sub></td><td>–</td></tr><tr><td>GW200316_215756; [4]</td><td>0.6<sup>+0.34</sup><sub>-0.39</sub></td><td>0.12<sup>+0.28</sup><sub>-0.1</sub></td><td>0.29<sup>+0.4</sup><sub>-0.21</sub></td><td>–</td></tr><tr><td>GW190728_064510; [3]</td><td>0.68<sup>+0.27</sup><sub>-0.31</sub></td><td>0.12<sup>+0.13</sup><sub>-0.06</sub></td><td>0.28<sup>+0.32</sup><sub>-0.2</sub></td><td>–</td></tr><tr><td>GW190720_000836; [3]</td><td>0.63<sup>+0.31</sup><sub>-0.27</sub></td><td>0.17<sup>+0.13</sup><sub>-0.11</sub></td><td>0.3<sup>+0.35</sup><sub>-0.21</sub></td><td>–</td></tr><tr><td>GW191204_171526; [4]</td><td>0.7<sup>+0.26</sup><sub>-0.26</sub></td><td>0.16<sup>+0.08</sup><sub>-0.05</sub></td><td>0.39<sup>+0.36</sup><sub>-0.26</sub></td><td>–</td></tr><tr><td>GW190708_232457; [3]</td><td>0.75<sup>+0.22</sup><sub>-0.28</sub></td><td>0.02<sup>+0.1</sup><sub>-0.08</sub></td><td>0.3<sup>+0.43</sup><sub>-0.24</sub></td><td>–</td></tr><tr><td>GW190930_133541; [3]</td><td>0.66<sup>+0.28</sup><sub>-0.37</sub></td><td>0.14<sup>+0.2</sup><sub>-0.14</sub></td><td>0.33<sup>+0.38</sup><sub>-0.23</sub></td><td>–</td></tr><tr><td>GW191216_213338; [4]</td><td>0.63<sup>+0.31</sup><sub>-0.28</sub></td><td>0.11<sup>+0.13</sup><sub>-0.06</sub></td><td>0.23<sup>+0.33</sup><sub>-0.16</sub></td><td>–</td></tr><tr><td>GW190707_093326; [3]</td><td>0.73<sup>+0.22</sup><sub>-0.22</sub></td><td>-0.03<sup>+0.08</sup><sub>-0.07</sub></td><td>0.25<sup>+0.38</sup><sub>-0.2</sub></td><td>–</td></tr><tr><td>GW190512_180714; [3]</td><td>0.54<sup>+0.37</sup><sub>-0.18</sub></td><td>0.03<sup>+0.13</sup><sub>-0.13</sub></td><td>0.23<sup>+0.35</sup><sub>-0.18</sub></td><td>–</td></tr><tr><td>GW191103_012549; [4]</td><td>0.67<sup>+0.29</sup><sub>-0.36</sub></td><td>0.21<sup>+0.16</sup><sub>-0.1</sub></td><td>0.4<sup>+0.4</sup><sub>-0.28</sub></td><td>–</td></tr><tr><td>GW170608; [5]</td><td>0.7<sup>+0.27</sup><sub>-0.37</sub></td><td>0.03<sup>+0.19</sup><sub>-0.07</sub></td><td>0.35<sup>+0.44</sup><sub>-0.27</sub></td><td>–</td></tr><tr><td>GW190725_174728; [3]</td><td>0.57<sup>+0.38</sup><sub>-0.31</sub></td><td>-0.04<sup>+0.26</sup><sub>-0.14</sub></td><td>0.37<sup>+0.46</sup><sub>-0.29</sub></td><td>–</td></tr><tr><td>GW200225_060421; [4]</td><td>0.73<sup>+0.24</sup><sub>-0.28</sub></td><td>-0.11<sup>+0.16</sup><sub>-0.28</sub></td><td>0.53<sup>+0.34</sup><sub>-0.38</sub></td><td>e.g. Antonelli et al. (2023); Zhang et al. (2023)</td></tr><tr><td>GW191105_143521; [4]</td><td>0.71<sup>+0.25</sup><sub>-0.3</sub></td><td>-0.02<sup>+0.13</sup><sub>-0.09</sub></td><td>0.29<sup>+0.46</sup><sub>-0.24</sub></td><td>–</td></tr><tr><td>GW151012; [6]</td><td>0.59<sup>+0.36</sup><sub>-0.35</sub></td><td>0.05<sup>+0.3</sup><sub>-0.21</sub></td><td>0.34<sup>+0.45</sup><sub>-0.26</sub></td><td>–</td></tr><tr><td>GW191129_134029; [4]</td><td>0.63<sup>+0.31</sup><sub>-0.3</sub></td><td>0.06<sup>+0.16</sup><sub>-0.08</sub></td><td>0.26<sup>+0.37</sup><sub>-0.2</sub></td><td>–</td></tr><tr><td>GW200202_154313; [4]</td><td>0.72<sup>+0.24</sup><sub>-0.31</sub></td><td>0.04<sup>+0.12</sup><sub>-0.06</sub></td><td>0.28<sup>+0.41</sup><sub>-0.21</sub></td><td>–</td></tr></tbody></table>

**Notes.** We provide the 90% credible intervals for the mass ratio ($q$), effective inspiral spin ($\chi_{\text{eff}}$), and effective precession spin ($\chi_{\text{p}}$) of events for which the median of the mass posterior distribution (for one or both components) lies above $10 M_\odot$ for events with $M < 11 M_\odot$ and below $15 M_\odot$ for events with $M > 11 M_\odot$. The events are ordered from having the most to least posterior support in the gap. We include the publication where the event is first reported: [1] Abbott et al. (2020); [2] Abbott et al. (2016b); [3] Abbott et al. (2021a); [4] Abbott et al. (2023a); [5] Abbott et al. (2017); [6] Abbott et al. (2016a); and some studies that have suggested a formation scenario different from isolated binary evolution for a given system.

where $\mu_{BH_H}$ and $\sigma_{BH_H}$ are the mean and width of the Gaussian and $m_{\text{min}}^{BH_H}$ and $m_{\text{max}}^{BH_H}$ are its lower and upper limits.

- A power-law component (Power Law or PL): this models everything that does not fit within these peaks, possibly from BHs formed via other channels,

$$ \pi(m_1|\Lambda)_{\text{PL}} \propto \begin{cases} m_1^{-\alpha} & m_{\text{min}}^{\text{PL}} < m_1 < m_{\text{max}}^{\text{PL}} \\ 0 & \text{otherwise,} \end{cases} \quad (6) $$

where $\alpha$ is the slope index of the power law and $m_{\text{min}}^{\text{BL}}$ and $m_{\text{max}}^{\text{BL}}$ are its lower and upper truncation points.

Each of these components is paired with different mass-ratio distributions, where

$$ \pi(q|m_1, \Lambda)_X \propto \begin{cases} q^{\beta_X} & m_{\text{min}}^X < m_2 < m_1 \\ 0 & \text{otherwise,} \end{cases} \quad (7) $$

where $X = BH_L$ is the mass-ratio distribution for the low-mass Gaussian component, $X = BH_H$ is for the high-mass Gaussian component, and $X = PL$ is for the power-law component. For each component, the mass ratio is constrained such that $m_2$ also has to be larger than the minimum mass of that component $m_{\text{min}}$.

We consider the fraction of BBHs in each component, where $f_{\text{peaks}}$ is the fraction of BBHs in the $BH_L$ and $BH_H$ peaks, and $f_1$ is the fraction of BBHs in the $BH_L$ peak compared to the total number of BBHs in the two peaks. Therefore, the fraction of BBHs in the $BH_L$ peak, in the $BH_H$ peak, and in the power-law component are $f_{\text{BH_L}} = f_{\text{peaks}} \times f_1$, $f_{\text{BH_H}} = f_{\text{peaks}} \times (1 - f_1)$, and $f_{\text{PL}} = 1 - f_{\text{peaks}}$, respectively. Details on the prior ranges of these parameters are provided in Table A.1.

Considering our model, we may naively expect that the peaked distributions may have different spin properties compared to the power-law component. The peaked distributions

are from BHs that form via isolated evolution, and therefore this population would generally be expected to form BBH systems with spins preferentially aligned with the orbital angular momentum. Supernova kicks may lead to misalignment of the spins, but the angle of this misalignment is generally considered to be modest (O'Shaughnessy et al. 2017; Stevenson et al. 2017; Gerosa et al. 2018; Rodriguez et al. 2016; Bavera et al. 2020). The BHs in the power-law component may have formed through different channels, and some may have formed dynamically. These systems are expected to form BH binaries with isotropically distributed spin orientations (Mandel & O'Shaughnessy 2010; Rodriguez et al. 2016; Zevin et al. 2017; Rodriguez et al. 2018; Doctor et al. 2020). If our sample contains hierarchical mergers, we may also expect to see greater spin magnitudes in this population (e.g. Yang et al. 2019; Gerosa & Berti 2017).

We consider two spin populations, one for the BBHs in the peaks ($X = \text{peaks}$), and one for the BBHs in the power law $X = \text{PL}$. The spin magnitude distribution is modelled by a truncated Gaussian:

$$ \pi(\chi_{1,2}|\Lambda)_X = \prod_{i=1}^{2} \mathcal{N}(\chi_i^X, \sigma_i^X), \quad (8) $$

where $\mu_X$ and $\sigma_X$ are the mean and width of the Gaussian, and the distribution is truncated at $\chi = 0$ and $\chi = 1$. The spin orientation is modelled as a mixture model of a uniform distribution and a truncated Gaussian (Talbot & Thrane 2017),

$$ \pi(\cos \theta_{1,2} | \Lambda)_X = \frac{1 - \zeta^X}{4} + \zeta^X \prod_{i=1}^{2} \mathcal{N}(\cos \theta_i | \sigma_i^X), \quad (9) $$

where $\theta$ is the angle of misalignment of a component with respect to the orbital angular momentum, $\sigma_i$ is the width of the
<figure><img src="image_3.png" /><figcaption>Fig. 2. Two-dimensional distributions of the mass ratio (q), effective inspiral spin (χ<sub>eff</sub>), and effective precession spin (χ<sub>p</sub>) of events polluting the m = 10–15 M<sub>☉</sub> gap (pink and blue shaded region; 90% credibility) and the median (pink and blue dots). The medians of the posteriors of the events outside the gap are also shown (black dots). Events polluting the gap (GW190412, GW151226, and GW200225), which may have formed via a channel different from isolated binary evolution, are highlighted with stars.</figcaption></figure>

Gaussian, and the mean is fixed to 1. The distribution is truncated at cos θ = −1 and cos θ = 1.

Hereafter, we refer to the model with these mass and spin prescriptions as COMPACTNESS PEAKS. The full COMPACTNESS PEAKS model is given by

$$
\begin{align}
& \pi(m_1, q, \chi_1, \chi_2, \cos \theta_1, \cos \theta_2 | \Lambda) \nonumber \\
& = (1 - f_{\text{peaks}}) \times \pi(m_1 | \Lambda)_{\text{PL}} \times \pi(q | m_1, \Lambda)_{\text{PL}} \nonumber \\
& \quad \times \pi(\chi_{1,2} | \Lambda)_{\text{PL}} \times \pi(\cos \theta_{1,2} | \Lambda)_{\text{PL}} \nonumber \\
& \quad + \pi(\chi_{1,2} | \Lambda)_{\text{peaks}} \times \pi(\cos \theta_{1,2} | \Lambda)_{\text{peaks}} \nonumber \\
& \quad \times [f_{\text{peaks}} f_1 \times p(m_1 | \Lambda)_{\text{BH}_L} \times \pi(q | m_1, \Lambda)_{\text{BH}_L} \times \nonumber \\
& \quad + f_{\text{peaks}} (1 - f_1) \times p(m_1 | \Lambda)_{\text{BH}_H} \times \pi(q | m_1, \Lambda)_{\text{BH}_H}]. \tag{10}
\end{align}
$$

In addition to masses and spins, we also fit for the redshift distri-
bution using the POWERLAW redshift evolution (Fishbach et al.
2018) used in Abbott et al. (2023b) for ease of comparison with
results from LVK analyses. We assume the same redshift distri-
bution for each component of the model (π(z|Λ)), and hence do
not discuss these analysis results.

3. Results

We analyse the population of BBH mergers and provide popula-
tion predictive distributions (PPDs). The PPD, for a given model,
is our best guess for the distribution of some source parameter
θ, averaged over the posterior for population parameters Λ (see
Table A.1 for descriptions):

$$
p_{\Lambda}(\theta|d) = \int d\Lambda p(\Lambda|d) \pi(\theta|\Lambda), \quad (11)
$$

<figure><img src="image_4.png" /><figcaption>Fig. 3. Population distributions for primary mass (m<sub>1</sub>, top panel) and mass ratio (q, bottom panel) using the COMPACTNESS PEAKS model for BBH mergers scaled by the fraction of events in the BH<sub>L</sub> Peak (blue), BH<sub>H</sub> Peak (pink), and POWER LAW (black) components. The solid curve is the median of the distribution and the shaded regions represent the 90% credible interval. The dashed line represents the median of the distribution for the MULTI PEAK model, which has the same mass ratio distribution for all components.</figcaption></figure>

where $p(\Lambda|d)$ represents the hyperposterior and $\pi(\theta|\Lambda)$ is the population prior.

We plot and highlight the components of COMPACTNESS PEAKS and compare these results to the MULTI PEAK model, which is the parametric model with the strongest preference used in Abbott et al. (2023b)<fnref n="3" />. This model is similar to COMPACTNESS PEAKS in that it has a power law and two peaks in the primary mass component; however, it does not have different mass-ratio or spin distributions for the different primary mass components, nor does it allow for the peaks to truncate at locations that are different from the minimum and maximum values of the power-law components. We calculate the Bayes factor (B) to provide a comparison of the models and determine which model best fits (or is ‘preferred by’) the data (Table 2). We find that the COMPACTNESS PEAKS model is strongly preferred over the MULTI PEAK model by log₁₀ B = 4.63 (i.e. preferred a by factor of ~40 000), supporting the hypothesis in Schneider et al. (2023).

In Figure 3 we show the astrophysical distributions for $m_1$ in the top panel and $q$ in the bottom panel. We find that by allowing

<footnote n="3">³ We use a variation of the MULTI PEAK model without smoothing ($\delta_m$ = 0) and different priors on the parameters. For details of this model, see Section B.4 in Abbott et al. (2021b); the priors used are provided in Table A.2. We use the same spin and redshift model as those used for COMPACTNESS PEAKS, but the same spin distribution is used for all components.</footnote>
<figure><img src="image_5.png" /><figcaption>Fig. 4. Fraction of binary BHs in the BH<sub>L</sub> Peak (blue), BH<sub>H</sub> Peak (pink), and POWER LAW (black) components using the COMPACTNESS PEAKS (solid lines) and MULTI PEAK (dashed lines) model analyses.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Fig. 5. Slope index of the mass-ratio distribution of the BH<sub>L</sub> Peak (blue), BH<sub>H</sub> Peak (pink), and POWER LAW (black) components using the COMPACTNESS PEAKS (solid lines) and COMPACTNESS PEAKS + NO PEAK PAIRING (dashed lines), which has relaxed conditions on the pairing of m<sub>1</sub> and m<sub>2</sub>.</figcaption></figure>

different mass-ratio distributions in each component and pairing of the peaked regions in m<sub>1</sub> and m<sub>2</sub>, the power-law component of the primary mass distribution becomes less steep – and shows a slope index of α = 1.82<sup>+1.01</sup><sub>-1.5</sub> – in comparison to the MULTI PEAK model, which exhibits a slope with an index of α = 2.82<sup>+0.93</sup><sub>-0.98</sub> at 90% credibility. We find that the edges of the component-mass gap, defined by the maximum of the lower-mass peak and the minimum of the upper-mass peak are m<sub>1</sub> = 10.92<sup>+1.62</sup><sub>-0.56</sub> M<sub>☉</sub> and m<sub>1</sub> = 13.09<sup>+3.02</sup><sub>-2.67</sub> M<sub>☉</sub>, respectively, at 90% credibility (see the posterior distributions in Figure B.1). We note that the higher mass peak in the COMPACTNESS PEAKS analysis is broader than the peak from the MULTI PEAK analysis, which is consistent with the broad BH<sub>H</sub> Peak distribution predicted in Schneider et al. (2023).

We also find shifts in the fraction of BBHs in each component with events in the power-law component moving to the peaked components when comparing the COMPACTNESS PEAKS and MULTI PEAK models, as shown in Figure 4. The 90% credible intervals for the fraction of binaries in each component are f<sub>PL</sub> = 0.15<sup>+0.29</sup><sub>-0.11</sub>, f<sub>BH<sub>L</sub></sub> = 0.73<sup>+0.13</sup><sub>-0.26</sub>, and f<sub>BH<sub>H</sub></sub> = 0.1<sup>+0.08</sup><sub>-0.05</sub> for COMPACTNESS PEAKS and f<sub>PL</sub> = 0.35<sup>+0.44</sup><sub>-0.28</sub>, f<sub>BH<sub>L</sub></sub> = 0.56<sup>+0.25</sup><sub>-0.42</sub>, and f<sub>BH<sub>H</sub></sub> = 0.08<sup>+0.1</sup><sub>-0.04</sub> for MULTI PEAK.

In Figure 5, we consider the slope index, β, of the mass ratio distribution for each component. We find that the power-law component has no preference for q ~ 1, but rather a uniform distribution, whereas the BH<sub>H</sub> Peak component strongly prefers more positive slopes (i.e. q ~ 1), which is consistent with findings from Li et al. (2022, 2024a). The BH<sub>L</sub> Peak component appears to return the prior for β. This result is because the truncation in the mass-ratio distribution – due to the strict pairing of the peaks in m<sub>1</sub> and m<sub>2</sub> – is driving its shape. However, if we relax this restriction by letting

$$ \pi(q|m_1, \Lambda)_{\text{BH}_L} \propto \begin{cases} q^{\beta_{\text{BH}_L}} & m_{\text{min}}^{\text{PL}} < m_2 < m_1 \\ 0 & \text{otherwise,} \end{cases} \quad (12) $$

and

$$ \pi(q|m_1, \Lambda)_{\text{BH}_H} \propto \begin{cases} q^{\beta_{\text{BH}_H}} & m_{\text{min}}^{\text{PL}} < m_2 < m_1 \\ 0 & \text{otherwise,} \end{cases} \quad (13) $$

such that m<sub>2</sub> is not confined to be within the component defined in m<sub>1</sub> space, we find the BH<sub>L</sub> Peak components also prefer positive slopes for the mass-ratio distribution. This analysis also reveals a preference for values of q ~ 1 for the peak components. This model variation is hereafter referred to as COMPACTNESS PEAKS + NO PEAK PAIRING. We also find that this model, with a relaxed condition on the mass-ratio model, is slightly more preferred than the COMPACTNESS PEAKS model, by log<sub>10</sub> B = 0.3. While this preference suggests that the strict pairing of m<sub>1</sub> and m<sub>2</sub> in the peaked component regions (as illustrated in Equation 7) may not be necessary to describe the population at present, the preference of q ~ 1 in the peaked components still remains and can be accounted for by the positive values of β, as shown in Figure 5. The preference for an equal mass ratio in these peaked regions is also supported by findings in Sadiq et al. (2024).

Figure 6 shows that the populations in the power-law and compactness peaks have different spin properties, where the peaked components have lower spins than the power-law component and are generally more aligned with respect to the orbital angular momentum. However, the result for the power-law component is somewhat prior driven (see Figure B.2) given the small fraction of BBHs in the power-law component of namely f ~ 0.15. We also observe that the power-law component has accounted for events with higher spins and very low spins, with reduced support for these values in the peaked distributions in comparison to the MULTI PEAK model.

We also considered a few variations of the COMPACTNESS PEAKS and COMPACTNESS PEAKS + NO PEAK PAIRING without separate spin distributions (i.e. the same values for the spin parameters μ<sub>χ</sub>, σ<sub>χ</sub>, σ<sub>t</sub>, ζ) for the power law and peak populations (Table 2). We find that support for these models is significantly driven by allowing for separate spin populations (log<sub>10</sub> B ~ 1–3). Interestingly, the degree of increase in model support with and without separate spin populations varies depending on the mass-ratio distribution of the model. This result suggests that mass and spin are not strictly independent parameters and this may be a consequence of the χ<sub>eff</sub> – q correlation (Callister et al. 2021; Adamcewicz & Thrane 2022).

## 4. Discussion

Motivated by the hypothesis presented in Schneider et al. (2023) for the presence of a gap in the component-mass spectrum of merging BBHs, we analysed the events in Abbott et al. (2023b) using a population designed to capture the structure of the compactness peaks and the degree of pairing of the peak regions in
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Population distributions for spin magnitude ($\chi$, top panel) and spin orientation (cos $\theta$, bottom panel) using the COMPACTNESS PEAKS model for BBHs mergers scaled by the fraction of events in the BH$_L$ and BH$_H$ PEAKS (purple), and POWER LAW (black) components. The solid curve is the median of the distribution and the shaded regions represent the 90% credible interval. The dashed line represents the median of the distribution obtained from the MULTI PEAK model analysis, which has the same spin distribution for all components.</figcaption></figure>

$m_1$ and $m_2$. We also explored whether there is a subpopulation of events that have properties that are different from those of the BBHs we expect from binary stripped stars (e.g. formed via dynamical evolution).

We find that there is a preference for the lower-mass peak to drop off sharply at $m_1 = 10.92_{-0.56}^{+1.62} M_{\odot}$ and the upper-mass peak to turn on at $m_1 = 13.09_{-2.67}^{+3.02} M_{\odot}$ (90% credible intervals). The locations of the possible edges of the peaks are consistent with predictions from Schneider et al. (2023), with a possible gap, but there is no clear evidence for a gap-like feature given that the posteriors of the edges of the peaks overlap. If we consider the case where there is no power-law component to sweep up the ‘pollutants’ (i.e. only the compactness peaks), support for a gap-like structure decreases, with the edges at $m_1 = 10.66_{-0.44}^{+1.28} M_{\odot}$ and $m_1 = 12.1_{-1.95}^{+1.78} M_{\odot}$ (90% credible intervals).

We also find mild support for the peaked populations to have a different spin distribution to the population defined by the power-law component. This result may hint at the possibility that the power-law component events were formed via alternative pathways to that leading to binary stripped stars in isolation, with 4–44% of binaries in the power-law component compared to 56–96% in the peaked components, and the majority of binaries, 48–87%, being found in the low-mass peak (90% credible intervals). However, more events are needed to constrain these results. With events from the fourth LVK observing run (O4),

which includes the KAGRA (Akutsu et al. 2019) observatory in the gravitational-wave network, we expect our sample size to grow by a factor of about 3 (Abbott et al. 2018; Iacovelli et al. 2022; Kiendrebeogo et al. 2023). This increase may help us to constrain the parameters defining the power-law component and the edges of the compactness peaks.

Recent studies with a focus on looking for hierarchical mergers in the population find evidence that the spin properties of observed BBHs change at $m \approx 40–50 M_{\odot}$ (e.g. Li et al. 2024b; Antonini et al. 2025; Pierra et al. 2024). In our work, we find that the upper mass compactness peak turns off at $m \approx 37–42 M_{\odot}$, which appears to be consistent with these findings. This result is not driven by the spin information for each component. When fixing each mass component to have the same spin distribution, the turn off is at $m \approx 37–41 M_{\odot}$. However, there is a preference for the power-law component to have a different spin distribution to the compactness peak populations. Given that $\chi_{\text{eff}}$ is a more tightly constrained parameter, it would be interesting to investigate the variation of $\chi_{\text{eff}}$ across these different populations and how it compares to our analysis with component spins.

In the later stages of the preparation of this paper, the authors became aware of a similar study exploring the proposed gap in the component masses (Adamcewicz et al. 2024). In their work, the authors aim to answer the question of whether there is a lack of BHs within $m \approx 10–15 M_{\odot}$. To do this, they directly model a dip-like feature in the component-mass distribution using a notch filter – the same as that proposed in Farah et al. (2022) – to explore the gap between the masses of neutron stars and BHs. Similar to our model, they also use two peaks that surround this gap. They do not find evidence for a dip-like feature in the distribution. Considering there are events ‘polluting the gap’, and given the mass posterior uncertainties are generally larger than the width of the gap, this result is not unexpected. The authors also conclude that the dip-like feature will likely not reveal itself with an increased catalogue of events from O4. In our work, we aim to answer a slightly different question, which pertains to whether or not there is a lack of BHs within $m \approx 10–15 M_{\odot}$ given the expectation of pollution in this region from BBHs formed via other channels. These two approaches lead to slightly different conclusions, with support for the edges of the possible gap-like feature more pronounced in our work, and strong support for our Schneider et al. (2023) inspired model ($\log_{10} \mathcal{B} = 4.63$). However, both studies conclude that there is no clear evidence for a gap-like feature in component mass.

In addition to this, the requirement of $q = m_2/m_1 \le 1$ may also be impacting our results. This prior restriction results in the skewed posteriors we observe for some events, where the mass ratio is close to unity. This effect is particularly evident in the events with masses of $<15 M_{\odot}$, where slight shifts in mass lead to much more unequal mass ratios. To mitigate this, future analyses extending this prior range to the individual events may help resolve the structure and edges of the compactness peaks. Extending this, it may be interesting to study this same population with a spin-sorting approach rather than a mass-sorting approach, where the objects are ordered according to their spin magnitude rather than their mass. This approach can also constrain spin information and has the greatest impact for events of $q \sim 1$ (e.g. Biscoveanu et al. 2021).

We highlight that in recent years, there has been a move towards data-driven models (e.g. Edelman et al. 2022; Tiwari 2023; Callister & Farr 2024) for the analysis of compact binary populations. These data-driven methods are valuable for finding structure in distribution without strong assumptions on the shape of the underlying distribution. However, astrophysically
<table><caption>**Table 2.** Population model comparison.</caption><thead><tr><th>Model</th><th>log<sub>10</sub> B</th><th>Δ log<sub>10</sub> L<sub>max</sub></th></tr></thead><tbody><tr><td>COMPACTNESS PEAKS</td><td>4.63</td><td>6.84</td></tr><tr><td>COMPACTNESS PEAKS (w/ same μ<sub>χ</sub>, σ<sub>χ</sub>, σ<sub>t</sub>, ζ)</td><td>3.78</td><td>6.02</td></tr><tr><td>COMPACTNESS PEAKS + NO PEAK PAIRING</td><td>4.93</td><td>7.49</td></tr><tr><td>COMPACTNESS PEAKS + NO PEAK PAIRING (w/ same β)</td><td>4.78</td><td>6.22</td></tr><tr><td>COMPACTNESS PEAKS + NO PEAK PAIRING (w/ same μ<sub>χ</sub>, σ<sub>χ</sub>, σ<sub>t</sub>, ζ)</td><td>2.88</td><td>5.88</td></tr><tr><td>COMPACTNESS PEAKS + NO PEAK PAIRING (w/ same β, μ<sub>χ</sub>, σ<sub>χ</sub>, σ<sub>t</sub>, ζ)</td><td>1.91</td><td>1.51</td></tr><tr><td>MULTI PEAK (δ<sub>m</sub> = 0)</td><td>0.0</td><td>0.0</td></tr><tr><td>POWERLAW PEAK (δ<sub>m</sub> = 0)</td><td>-1.97</td><td>-3.44</td></tr></tbody></table>

**Notes.** Log Bayes factors and maximum log likelihood differences for models compared to the MULTI PEAK model (Abbott et al. 2021b). We also include results for the POWERLAW PEAK (Talbot & Thrane 2018), a commonly used model for the population analysis of BBHs (e.g. Abbott et al. 2023b).

informed parametric models provide a clear method for astro-physical interpretations of the features observed and are ideal for model comparison. In this data-rich era, both data-driven and phenomenological models are needed: one to identify the possible new structure and one to interpret the astrophysical implications for the features we observe.

## Data availability

Supplementary material including analysis inputs, posterior samples and additional plots are available here github.com/shanikagalaudage/bbh_compactness_peaks

**Acknowledgements.** We thank Giuliano Iorio and the anonymous reviewer for their helpful comments on this work. Galaudage and Lamberts are supported by the ANR COSMERGE project, grant ANR-20-CE31-001 of the French Agence Nationale de la Recherche. The authors are grateful for computational resources provided by the LIGO Laboratory and supported by National Science Foundation Grants PHY-0757058 and PHY-0823459. This material is based upon work supported by NSF's LIGO Laboratory which is a major facility fully funded by the National Science Foundation

## References

Aasi, J., Abbott, B. P., Abbott, R., et al. 2015, CQG, 32, 074001

Abbott, B.P., Abbott, R., Abbott, T.D., et al. 2016a, PRX, 6, 041015; Erratum: PRX 8, 039903 (2018)

Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2016b, PRL, 116, 241103

Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2017, ApJ, 851, L35

Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2018, Liv. Rev. Rel., 21, 3

Abbott, R., Abbott, T. D., Abraham, S., et al. 2020, PRD, 102, 043015

Abbott, R., Abbott, T. D., Abraham, S., et al. 2021a, PRX, 11, 021053

Abbott, R., Abbott, T. D., Abraham, S., et al. 2021b, ApJ, 913, L7

Abbott, R., Abbott, T. D., Acernese, F., et al. 2023a, PRX, 13, 041039

Abbott, R., Abbott, T. D., Acernese, F., et al. 2023b, PRX, 13, 011048

Acernese, F., Agathos, M., Agatsuma, K., et al. 2015, CQG, 32, 024001

Adamcewicz, C., & Thrane, E. 2022, MNRAS, 517, 3928

Adamcewicz, C., Lasky, P. D., Thrane, E., & Mandel, I. 2024, ApJ, 975, 253

Ajith, P., Hannam, M., Husa, S., et al. 2011, PRL, 106, 241101

Akutsu, T., Ando, M., Arai, K., et al. 2019, Nat. Astron., 3, 35

Antonelli, A., Kritos, K., Ng, K. K. Y., Cotesta, R., & Berti, E. 2023, PRD, 108, 084044

Antonini, F., Gieles, M., Dosopoulou, F., & Chattopadhyay, D. 2023, MNRAS, 522, 466

Antonini, F., Romero-Shaw, I. M., & Callister, T. 2025, Phys. Rev. Lett., 134, 011401

Ashton, G., Hübner, M., Lasky, P. D., et al. 2019, ApJS, 241, 27

Bavera, S. S., Fragos, T., Qin, Y., et al. 2020, A&A, 635, A97

Belczynski, K., Bulik, T., Fryer, C. L., et al. 2010, ApJ, 714, 1217

Biscoveanu, S., Isi, M., Vitale, S., & Varma, V. 2021, PRL, 126, 171103

Blanchet, L., Damour, T., Iyer, B. R., Will, C. M., & Wiseman, A. G. 1995, PRL, 74, 3515

Broekgaarden, F. S., Berger, E., Stevenson, S., et al. 2022, MNRAS, 516, 5737

Burrows, A., Vartanyan, D., & Wang, T. 2023, ApJ, 957, 68

Burrows, A., Wang, T., & Vartanyan, D. 2024, ApJ, 964, L16

Callister, T. A., & Farr, W. M. 2024, PRX, 14, 021005

Callister, T. A., Haster, C.-J., Ng, K. K. Y., Vitale, S., & Farr, W. M. 2021, ApJ, 922, L5

Chatterjee, S., Rodriguez, C. L., Kalogera, V., & Rasio, F. A. 2017, ApJ, 836, L26

Chieffi, A., & Limongi, M. 2020, ApJ, 890, 1

Costa, G., Bressan, A., Mapelli, M., et al. 2021, MNRAS, 501, 4514

Croon, D., & Sakstein, J. 2023, arXiv e-prints [arXiv:2312.13459]

Doctor, Z., Wysocki, D., O'Shaughnessy, R., Holz, D. E., & Farr, B. 2020, ApJ, 893, 35

Edelman, B., Doctor, Z., Godfrey, J., & Farr, B. 2022, ApJ, 924, 101

Essick, R., & Farr, W. 2022, arXiv e-prints [arXiv:2204.00461]

Farag, E., Renzo, M., Farmer, R., Chidester, M. T., & Timmes, F. X. 2022, ApJ, 937, 112

Farah, A. M., Fishbach, M., Essick, R., Holz, D. E., & Galaudage, S. 2022, ApJ, 931, 108

Farah, A. M., Edelman, B., Zevin, M., et al. 2023, ApJ, 955, 107

Farmer, R., Renzo, M., de Mink, S. E., Marchant, P., & Justham, S. 2019, ApJ, 887, 53

Farmer, R., Renzo, M., de Mink, S., Fishbach, M., & Justham, S. 2020, ApJ, 902, L36

Finn, L. S., & Chernoff, D. F. 1993, PRD, 47, 2198

Fishbach, M., Holz, D. E., & Farr, W. M. 2018, ApJ, 863, L41

Gerosa, D., & Berti, E. 2017, PRD, 95, 124046

Gerosa, D., Berti, E., O'Shaughnessy, R., et al. 2018, PRD, 98, 084036

Gerosa, D., Vitale, S., & Berti, E. 2020, PRL, 125, 101103

Giacobbo, N., & Mapelli, M. 2018, MNRAS, 480, 2011

Golomb, J., Isi, M., & Farr, W. 2024, ApJ, 976, 121

Hannam, M., Schmidt, P., Bohé, A., et al. 2014, PRL, 113, 151101

Heger, A., & Woosley, S. E. 2002, ApJ, 567, 532

Hendriks, D. D., van Son, L. A. C., Renzo, M., Izzard, R. G., & Farmer, R. 2023, MNRAS, 526, 4130

Iacovelli, F., Mancarella, M., Foffa, S., & Maggiore, M. 2022, ApJ, 941, 208

Iorio, G., Mapelli, M., Costa, G., et al. 2023, MNRAS, 524, 426

Kiendrebeogo, R. W., Farah, A. M., Foley, E. M., et al. 2023, ApJ, 958, 158

Kruckow, M. U., Tauris, T. M., Langer, N., Kramer, M., & Izzard, R. G. 2018, MNRAS, 481, 1908

Li, Y.-J., Wang, Y.-Z., Tang, S.-P., & et al. 2022, ApJ, 933, L14

Li, Y. J., Tang, S. P., Gao, S. J., Wu, D. C., & Wang, Y. Z. 2024a, ApJ, 977, 67

Li, Y.-J., Wang, Y.-Z., Tang, S.-P., & Fan, Y.-Z. 2024b, PRL, 133, 051401

Liu, B., & Lai, D. 2021, MNRAS, 502, 2049

Mandel, I., & Broekgaarden, F. S. 2022, Liv. Rev. Rel., 25, 1

Mandel, I., & Farmer, A. 2022, Phys. Rept., 955, 1

Mandel, I., & O'Shaughnessy, R. 2010, CQG, 27, 114007

Mapelli, M. 2020, Proc. Int. Sch. Phys. Fermi, 200, 87

Mapelli, M., Spera, M., Montanari, E., et al. 2020, ApJ, 888, 76

O'Connor, E., & Ott, C. D. 2011, ApJ, 730, 70

O'Shaughnessy, R., Gerosa, D., & Wysocki, D. 2017, PRL, 119, 011101

Patton, R. A., & Sukhbold, T. 2020, MNRAS, 499, 2803

Pierra, G., Mastrogiovanni, S., & Perriès, S. 2024, arXiv e-prints [arXiv:2406.01679]

Poisson, E., & Will, C. M. 1995, PRD, 52, 848
Rahman, N., Janka, H.-T., Stockinger, G., & Woosley, S. 2022, MNRAS, 512, 4503

Rodriguez, C. L., Zevin, M., Pankow, C., Kalogera, V., & Rasio, F. A. 2016, ApJ, 832, L2

Rodriguez, C. L., Amaro-Seoane, P., Chatterjee, S., & Rasio, F. A. 2018, PRL, 120, 151101

Romero-Shaw, I. M., Talbot, C., Biscoveanu, S., et al. 2020, MNRAS, 499, 3295

Sadiq, J., Dent, T., & Wysocki, D. 2022, PRD, 105, 123014

Sadiq, J., Dent, T., & Gieles, M. 2024, ApJ, 960, 65

Safarzadeh, M., & Hotokezaka, K. 2020, ApJ, 897, L7

Santamaria, L., Ohme, F., Ajith, P., et al. 2010, PRD, 82, 064016

Schmidt, P., Ohme, F., & Hannam, M. 2015, PRD, 91, 024043

Schneider, F. R. N., Podsiadlowski, P., & Müller, B. 2021, A&A, 645, A5

Schneider, F. R. N., Podsiadlowski, P., & Laplace, E. 2023, ApJ, 950, L9

Speagle, J. S. 2020, MNRAS, 493, 3132

Spera, M., Mapelli, M., & Bressan, A. 2015, MNRAS, 451, 4086

Spera, M., Trani, A. A., & Mencagli, M. 2022, Galaxies, 10, 76

Stevenson, S., Berry, C. P. L., & Mandel, I. 2017, MNRAS, 471, 2801

Stevenson, S., Sampson, M., Powell, J., et al. 2019, ApJ, 882, 121

Sukhbold, T., & Woosley, S. 2014, ApJ, 783, 10

Talbot, C., & Thrane, E. 2017, PRD, 96, 023012

Talbot, C., & Thrane, E. 2018, ApJ, 856, 173

Talbot, C., Smith, R., Thrane, E., & Poole, G. B. 2019, PRD, 100, 043030

Temaj, D., Schneider, F. R. N., Laplace, E., Wei, D., & Podsiadlowski, P. 2024, A&A, 682, A123

Thrane, E., & Talbot, C. 2019, PASA, 36, e010; Erratum: PASA 37, e036 (2020)

Tiwari, V. 2022, ApJ, 928, 155

Tiwari, V. 2023, MNRAS, 527, 298

Tiwari, V., & Fairhurst, S. 2021, ApJ, 913, L19

van Son, L. A. C., de Mink, S. E., Renzo, M., et al. 2022, ApJ, 940, 184

Vink, J. S., & de Koter, A. 2005, A&A, 442, 587

Wong, K. W. K., & Cranmer, M. 2022, in 39th International Conference on Machine Learning Conference

Yang, Y., Bartos, I., Gayathri, V., et al. 2019, PRL, 123, 181101

Zevin, M., Pankow, C., Rodriguez, C. L., et al. 2017, ApJ, 846, 82

Zhang, R. C., Fragione, G., Kimball, C., & Kalogera, V. 2023, ApJ, 954, 23
## Appendix A: Population model details

In this section we provide details about the population models described above in Section 2.2. We include a summary of the parameters for the COMPACTNESS PEAKS and MULTI PEAK models and the prior distribution used for each parameter in Tables A.1 and A.2 respectively. The prior distributions are indicated using abbreviations: for example, U(0, 1) translates to uniform on the interval (0, 1).

<table><caption>**Table A.1.** Summary of parameters in COMPACTNESS PEAKS population model.</caption><thead><tr><th>Parameter</th><th>Description</th><th>Prior</th></tr></thead><tbody><tr><td>f<sub>peaks</sub></td><td>Fraction of BHs in the Gaussian components of the m<sub>1</sub> distribution.</td><td>U(0, 1)</td></tr><tr><td>f<sub>1</sub></td><td>Fraction of BHs in the lower-mass Gaussian component of the Gaussian components of the m<sub>1</sub> distribution.</td><td>U(0, 1)</td></tr><tr><td>α</td><td>Spectral index for the power-law of the m<sub>1</sub> distribution.</td><td>U(-4, 8)</td></tr><tr><td>m<sup>PL</sup><sub>min</sub> [M<sub>☉</sub>]</td><td>Minimum mass of the power-law component of the m<sub>1</sub> distribution.</td><td>U(2, 10)</td></tr><tr><td>m<sup>PL</sup><sub>max</sub> [M<sub>☉</sub>]</td><td>Maximum mass of the power-law component of the m<sub>1</sub> distribution.</td><td>U(60, 100)</td></tr><tr><td>β<sub>PL</sub></td><td>Spectral index for the power-law of the q distribution paired with the power law component of the m<sub>1</sub> distribution.</td><td>U(-4, 12)</td></tr><tr><td>μ<sub>BH<sub>L</sub></sub> [M<sub>☉</sub>]</td><td>Mean of the lower-mass Gaussian component in the m<sub>1</sub> distribution.</td><td>U(8, 12)</td></tr><tr><td>σ<sub>BH<sub>L</sub></sub> [M<sub>☉</sub>]</td><td>Width of the lower-mass Gaussian component in the m<sub>1</sub> distribution.</td><td>U(1, 4)</td></tr><tr><td>m<sup>BH<sub>L</sub></sup><sub>min</sub> [M<sub>☉</sub>]</td><td>Minimum mass of the lower-mass Gaussian component of the m<sub>1</sub> distribution.</td><td>U(5, 8)</td></tr><tr><td>m<sup>BH<sub>L</sub></sup><sub>max</sub> [M<sub>☉</sub>]</td><td>Maximum mass of the lower-mass Gaussian component of the m<sub>1</sub> distribution.</td><td>U(10, 15)</td></tr><tr><td>β<sub>BH<sub>L</sub></sub></td><td>Spectral index for the power-law of the q distribution paired with the power law component of the m<sub>1</sub> distribution.</td><td>U(-4, 12)</td></tr><tr><td>μ<sub>BH<sub>H</sub></sub> [M<sub>☉</sub>]</td><td>Mean of the upper-mass Gaussian component in the m<sub>1</sub> distribution.</td><td>U(12, 50)</td></tr><tr><td>σ<sub>BH<sub>H</sub></sub> [M<sub>☉</sub>]</td><td>Width of the upper-mass Gaussian component in the m<sub>1</sub> distribution.</td><td>U(4, 20)</td></tr><tr><td>m<sup>BH<sub>H</sub></sup><sub>min</sub> [M<sub>☉</sub>]</td><td>Minimum mass of the upper-mass Gaussian component of the m<sub>1</sub> distribution.</td><td>U(10, 25)</td></tr><tr><td>m<sup>BH<sub>H</sub></sup><sub>max</sub> [M<sub>☉</sub>]</td><td>Maximum mass of the upper-mass Gaussian component of the m<sub>1</sub> distribution.</td><td>U(30, 50)</td></tr><tr><td>β<sub>BH<sub>H</sub></sub></td><td>Spectral index for the power-law of the q distribution paired with the upper-mass Gaussian component of the m<sub>1</sub> distribution.</td><td>U(-4, 12)</td></tr><tr><td>μ<sup>X</sup><sub>χ</sub></td><td>Mean of the spin magnitude distribution for a given component.</td><td>U(0, 1)</td></tr><tr><td>σ<sup>X</sup><sub>χ</sub></td><td>Width of the spin magnitude distribution for a given component</td><td>U(0.1, 2)</td></tr><tr><td>σ<sup>X</sup><sub>t</sub></td><td>Width of the spin orientation distribution for the subpopulation aligned with respect to the orbital angular momentum for a given component</td><td>U(0.1, 4)</td></tr><tr><td>ζ<sup>X</sup></td><td>Fraction of BBHs in the subpopulation aligned with respect to the orbital angular momentum for a given component</td><td>U(0, 1)</td></tr></tbody></table>

**Notes.** The prior ranges are motivated by the results in Schneider et al. (2023)
<table><caption>**Table A.2.** Summary of parameters in MULTI PEAK model.</caption><thead><tr><th>Parameter</th><th>Description</th><th>Prior</th></tr></thead><tbody><tr><td>α</td><td>Spectral index for the power-law of the m<sub>1</sub> distribution.</td><td>U(-4, 12)</td></tr><tr><td>m<sub>min</sub> [M<sub>☉</sub>]</td><td>Minimum mass of the power-law component of the m<sub>1</sub> distribution.</td><td>U(2, 10)</td></tr><tr><td>m<sub>max</sub> [M<sub>☉</sub>]</td><td>Maximum mass of the power-law component of the m<sub>1</sub> distribution.</td><td>U(30, 100)</td></tr><tr><td>β</td><td>Spectral index for the power-law of the q distribution paired with the power law component of the m<sub>1</sub> distribution.</td><td>U(-4, 12)</td></tr><tr><td>f<sub>peaks</sub> (or λ)</td><td>Fraction of BHs in the Gaussian components of the m<sub>1</sub> distribution.</td><td>U(0, 1)</td></tr><tr><td>f<sub>1</sub> (or λ<sub>1</sub>)</td><td>Fraction of BHs in the lower-mass Gaussian component of the Gaussian components of the m<sub>1</sub> distribution.</td><td>U(0, 1)</td></tr><tr><td>μ<sub>m,1</sub> [M<sub>☉</sub>]</td><td>Mean of the lower-mass Gaussian component in the m<sub>1</sub> distribution.</td><td>U(5, 20)</td></tr><tr><td>σ<sub>m,1</sub> [M<sub>☉</sub>]</td><td>Width of the lower-mass Gaussian component in the m<sub>1</sub> distribution.</td><td>U(1,5)</td></tr><tr><td>μ<sub>m,2</sub> [M<sub>☉</sub>]</td><td>Mean of the upper-mass Gaussian component in the m<sub>1</sub> distribution.</td><td>U(20, 50)</td></tr><tr><td>σ<sub>m,2</sub> [M<sub>☉</sub>]</td><td>Width of the upper-mass Gaussian component in the m<sub>1</sub> distribution.</td><td>U(1, 10)</td></tr></tbody></table>

**Notes.** This model, with similar prior ranges, was the most preferred parametric model used in Abbott et al. (2023b); see Appendix D. The minimum and maximum of the Gaussian components are equivalent to the minimum and maximum of the power law component.
## Appendix B: Hyperparameter posteriors

In this section we provide the hyperparameter posteriors for the parameters governing the peak components of the mass distribution (Figure B.1) and the spin parameters for the peaked and powerlaw components (Figure B.2).

<figure><img src="image_8.png" /><figcaption>**Fig. B.1.** Mass hyperparameter posteriors of the peak components in the COMPACTNESS PEAKS model. Left panel: the lower mass (BH<sub>L</sub>) peak where μ<sub>BH<sub>L</sub></sub> is the mean of the peak, σ<sub>BH<sub>L</sub></sub> is the width of the peak and m<sup>BH<sub>L</sub></sup><sub>min</sub> and m<sup>BH<sub>L</sub></sup><sub>max</sub> are the minimum and maximum mass of the peak. Right panel: the upper mass (BH<sub>H</sub>) peak where μ<sub>BH<sub>H</sub></sub> is the mean of the peak, σ<sub>BH<sub>H</sub></sub> is the width of the peak and m<sup>BH<sub>H</sub></sup><sub>min</sub> and m<sup>BH<sub>H</sub></sup><sub>max</sub> are the minimum and maximum mass of the peak. The credible intervals on the 2D posterior distributions are at 1σ, 2σ and 3σ, using increasingly light shading. The intervals on the 1D posterior distributions are at 1σ.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Fig. B.2.** Spin hyperparameter posteriors of the COMPACTNESS PEAKS model, where μ<sub>χ</sub> is the mean of the spin magnitude distribution, σ<sub>χ</sub> is the width of the spin magnitude distribution, ζ is the fraction of BBHs in the subpopulation aligned with respect to the orbital angular momentum and σ<sub>t</sub> is the width of the spin orientation distribution for the subpopulation aligned with respect to the orbital angular momentum for a given component. Left panel: The spin parameters paired with the BH<sub>L</sub> and BH<sub>H</sub> peaks Right panel: The spin parameters paired with the powerlaw (PL) component. The credible intervals on the 2D posterior distributions are at 1σ, 2σ and 3σ, using increasingly light shading. The intervals on the 1D posterior distributions are at 1σ.</figcaption></figure>