# First constraints on compact binary environments from LIGO-Virgo data

Giada Caneva Santoro<EMAIL_ADDRESS>Institut de Física d’Altes Energies
(IFAE), The Barcelona Institute of Science and Technology, UAB Campus, E-08193
Barcelona, Spain  Soumen Roy<EMAIL_ADDRESS>Nikhef, Science Park 105,
1098 XG Amsterdam, The Netherlands  Institute for Gravitational and Subatomic
Physics (GRASP), Utrecht University, Princetonplein 1, 3584 CC Utrecht, The
Netherlands  Rodrigo Vicente<EMAIL_ADDRESS>Institut de Física d’Altes
Energies (IFAE), The Barcelona Institute of Science and Technology, UAB
Campus, E-08193 Barcelona, Spain  Maria Haney  Nikhef, Science Park 105, 1098
XG Amsterdam, The Netherlands  Ornella Juliana Piccinni  Institut de Física
d’Altes Energies (IFAE), The Barcelona Institute of Science and Technology,
UAB Campus, E-08193 Barcelona, Spain  Walter Del Pozzo  Dipartimento di Fisica
“Enrico Fermi”, Università di Pisa and INFN Sezione di Pisa, Pisa I-56127,
Italy  Mario Martinez  Institut de Física d’Altes Energies (IFAE), The
Barcelona Institute of Science and Technology, UAB Campus, E-08193 Barcelona,
Spain  Catalan Institution for Research and Advanced Studies (ICREA), E-08010
Barcelona, Spain

###### Abstract

The LIGO-Virgo analyses of signals from compact binary mergers observed so far
have assumed isolated binary systems in vacuum, neglecting the potential
presence of astrophysical environments. We present here the first
investigation of environmental effects on the events of GWTC-1. We find no
evidence for the presence of environmental effects. Most of the events
decisively exclude the scenario of dynamical fragmentation of massive stars as
their formation channel. GW170817 results in the most stringent upper bound on
the medium density ($\lesssim 21\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$). We
find that environmental effects can substantially bias the recovered
parameters in the vacuum model, even when these effects are not detectable. We
forecast that the Einstein Telescope and B-DECIGO will be able to probe the
environmental effects of accretion disks and superradiant boson clouds on
compact binaries.

Introduction. In-depth analyses of gravitational-wave (GW) data have been
routinely conducted by the LIGO-Virgo-KAGRA (LVK) collaboration, through
parameter estimation, population studies, cosmology, and tests of general
relativity (GR) [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14]. These
analyses relied on the assumption that the sources of GWs are in a _vacuum_
environment. However, there is a growing interest in exploring the potential
effects of astrophysical environments on these observations, particularly with
regard to binary black hole (BBH) formation in dense regions, such as star
clusters [15, 16, 17, 18, 19, 20], active galactic nuclei (AGN) accretion
disks [21, 22, 23], or through the dynamical fragmentation of very massive
stars [24]. The physical processes that take place in environments that
harbour black holes (BHs) significantly impact their formation, dynamics and
evolution, and it is therefore crucial to characterise them. These effects are
also relevant to multi-messenger astronomy, as electromagnetic counterparts
are expected to accompany GWs for BBH mergers in non-vacuum environments [24,
25, 26, 27].

Most studies of environmental effects on GW signals (e.g., [28, 29, 30, 31,
32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47]) have focused
on sources relevant to the frequency bandwidth of the LISA experiment [48],
like intermediate and extreme mass ratio inspirals. These works forecast that
LISA will be sensitive to the imprints of astrophysical environments,
potentially being able to identify the type of environment and reconstruct its
model parameters [41]. The reason for the generalized focus on LISA sources is
mainly twofold: (i) the long observation window of the inspiral stage, which
allows for the accumulation of usually small environmental effects, and (ii)
LISA’s sensitivity to low frequencies, which are the most affected by
environmental dephasing. LISA is expected to become operational only by 2037
[49] and so it is important to quantify the _current_ ability of LVK to “see”
environments.

Here, for the first time, a systematic Bayesian analysis is performed in the
framework of the existing LIGO-Virgo GWTC-1 data, testing for the presence of
environmental effects, trying to constrain the environmental properties, and
determining whether the inclusion of environmental effects in waveforms can
bias the estimations of the GW source parameters. A previous study using LIGO-
Virgo data focused on the first observed event GW150914 by comparing numerical
relativity waveforms in vacuum versus environment [50], which is
computationally challenging and out of scope for low-mass binaries. We
perform, for the first time, a parameterized post-Newtonian (PN) test of
environmental effects with LIGO-Virgo data. The test allows for modifications
to the in-vacuum GW phase of the binary that could be induced by environments
during the inspiral stage. Similar types of parameterized tests, considering
waveform deformations at PN orders higher than the ones induced by
environmental effects, are routinely performed by LVK to test GR [11, 12, 13].
The same framework can be used to detect non-GR signals [51]. However, these
methods do not specifically test scenarios of a binary immersed in an
astrophysical medium.

We use geometrized units with $G=c=1$. For reference, a supermassive BH
(SMBHs) with mass a $M_{\mathrm{SMBH}}$ can exist in environments with
densities as large as: $\sim
10^{-8}(10^{5}M_{\odot}/M_{\mathrm{SMBH}})\,\mathrm{g}/\mathrm{cm}^{3}$ for
thick accretion disks and $\sim
0.1(10^{5}M_{\odot}/M_{\mathrm{SMBH}})^{\frac{7}{10}}\,\mathrm{g}/\mathrm{cm}^{3}$
for thin ones [30], $\sim 10^{-6}\,\mathrm{g}/\mathrm{cm}^{3}$ for cold dark
matter spikes [52, 53], and $\sim
0.1(10^{5}M_{\odot}/M_{\mathrm{SMBH}})^{2}\,\mathrm{g}/\mathrm{cm}^{3}$ for
superradiant clouds of light bosons [54]. Compact binaries formed through the
dynamical fragmentation of a massive star [24] are expected to subsequently
merge in a gaseous environment with density $\gtrsim
10^{7}\,\mathrm{g}/\mathrm{cm}^{3}$ [50].

Environmental effects. When in astrophysical environments, the phase evolution
of compact binaries is expected to be slightly modified with respect to the
vacuum “chirp” by effects like accretion [55, 56, 57] or dynamical friction
(DF) [58, 59, 60, 61, 62] (for details, see Ref. [63]). While the specifics of
these effects and their relative impact on waveforms depend on the particular
environment and binary source and can only be completely captured by numerical
simulations, there are generic (agnostic) features that can be well captured
by semi-analytic expressions.

While environmental effects may lead to an increase in the BBH eccentricity
[64, 23], for small enough densities the GW radiation reaction dominates
leading to the circularization of these binaries [65]. For quasi-circular
orbits DF from the gravitational wake in the medium results in an energy loss

$\dot{E}_{\mathrm{DF}}\approx\frac{4\pi\rho
M^{2}\mathcal{I}(v,\eta)}{v}\Big{(}\frac{1-3\eta}{\eta}\Big{)},$ (1)

with $v\coloneqq(\pi M_{z}f)^{\frac{1}{3}}$, where $f/2$ denotes the orbital
frequency and we define $M_{z}\coloneqq(1+z)M$, with $M$ the binary’s total
mass (in the source frame) and $z$ the cosmological redshift to the source.
The symmetric mass ratio is defined as $\eta\coloneqq m_{1}m_{2}/M^{2}$,
$\rho$ is the (local) average mass density of the environment, and
$\mathcal{I}$ is a function of $v$ and $\eta$ which depends on the type of
environment; sources in the LVK mass range in a gaseous media with asymptotic
speed of sound $c_{\mathrm{s}}$ have $\mathcal{I}\sim
O\big{(}c_{\mathrm{s}}/10^{-3}\big{)}$, and so we treat $\mathcal{I}$ as a
constant in our analysis. Analytic expressions for $\mathcal{I}$ can be found,
e.g., in [66] for gaseous media, or [67, 68, 69, 70, 71] for (ultra-light)
dark matter. In this work we assume a static homogeneous density profile,
neglecting the feedback of the BBH on the environment (see, e.g., [72]).

Accretion causes the masses of the binary components to change through the
inspiral, which also affects the evolution of the orbital phase 111The
transfer of momentum via accretion is also responsible for an additional
(hydrodynamic) drag force which may be included within the factor
$\mathcal{I}$ in Eq. (1) [57, Barausse:2007dy, 70].. In collisional media,
numerical simulations (e.g., [74, 75, 76, 77, 78]) show that binary accretion
can be well approximated by Bondi-Hoyle-Lyttleton accretion (BHLA), with the
usual fudge parameter $\lambda\sim O(1)$ [79]. In media constituted by
particles with larger mean free path (like particle dark matter overdensities
[30, 80], or plasmas around BHs [81]), accretion is less efficient and it is
better described by collisionless accretion (CA) [82, 83].

The environment can alter the evolution of the GW phase also through other
effects (like the medium’s gravitational potential [84, 80, 85], acceleration
of the binary’s center of mass [33, 86], planet-like migration [28, 29, 39],
and others), but for LVK binaries the effects imparted by DF and accretion are
expected to be the most important [30, 34]. We consider diluted environments
causing only small corrections to the dominant “chirp” due to vacuum
radiation-reaction. As shown in _Supplemental Material_ (see also Ref. [30]),
at _linear order_ in $\rho M^{2}$, the GW phase in environments differs from
the vacuum one by additive terms $(3/128\eta)\widehat{\delta\Phi}_{k}v^{k-5}$,
where $k=-9$ ($-4.5$PN) for CA and $k=-11$ ($-5.5$PN) for DF and BHLA, with
the coefficients given in terms of physical parameters in Tab. 1. The same
$-5.5$PN effect on the GW phase has also been observed in numerical relativity
simulations of BBHs in gaseous media (see Fig. 2 in Ref. [50]). We have also
checked that, up to a $O(1)$ factor, the coefficients $\beta_{k}$ of Tab. 1
can reproduce the environment effect on the chirp observed in those
simulations.

Table 1: Dependence of environmental dephasing coefficients on the physical
parameters (see _Supplemental Material_). We write
$\widehat{\delta\Phi}_{k}=-\beta_{k}\tilde{\rho}M^{2}$, with $\tilde{\rho}$
denoting $\rho$ for CA, $\lambda\rho$ for BHLA, and $\mathcal{I\rho}$ for DF.
The effects enter at $k/2$-th PN order in the GW phase.

Effect |  $\beta_{k}$ |  $k$
---|---|---
CA |  $\dfrac{125\pi(1-3\eta)}{357\eta^{2}}$ |  $-9$
BHLA |  $\dfrac{125\pi[1-5\eta(1-\eta)]}{1824\eta^{4}}$ |  $-11$
DF |  $\dfrac{25\pi(1-3\eta)}{304\eta^{3}}$ |  $-11$

Analysis Setup. We consider an agnostic dephasing parameter $\delta\Phi_{k}$,
such that, in frequency domain, the GW phase with a specific environmental
effect is

$\phi^{\rm env}=\phi^{\rm{vac}}+\frac{3}{128\eta}\delta\Phi_{k}v^{k-5},$ (2)

where $k=-9$ for CA, and $k=-11$ for BHLA or DF. To generate the waveforms
with the environmental correction, we followed the model-agnostic framework of
parameterized tests of GR [87, 88, 89]. Throughout this work, we use the
phenomenological inspiral-merger-ringdown models IMRPhenomPv2 for binary black
hole (BBH) [90, 91, 92, 93] and IMRPhenomPv2_NRTidalv2 for binary neutron star
(BNS) systems [94], which incorporate spin-induced precession effects.

|  GW150914 |  GW151012 |  GW151226 |  GW170104 |  GW170608 |  GW170729 |  GW170809 |  GW170814 |  GW170817 |  GW170818 |  GW170823
---|---|---|---|---|---|---|---|---|---|---|---
$\delta\Phi_{-9}$ | $-2.09$ | $0.86$ | $-4.20$ | $-1.29$ | $-4.95$ | – | $-1.91$ | $-3.17$ | $-5.45$ | $-2.55$ | –
$\delta\Phi_{-11}$ | $-3.30$ | $0.77$ | $-5.52$ | $-3.33$ | $-6.17$ | – | $-2.59$ | $-2.81$ | $-6.47$ | $-2.57$ | –
Table 2: Logarithmic Bayes factor ($\log_{10}\mathcal{B}_{\rm vac}^{\rm env}$)
for GWTC-1 events. For the events GW170729 and GW170823, we could not find
informative $\delta\Phi_{k}$ posteriors even with a broad prior range due to
the low SNR of their inspiral.

We carry out a Bayesian parameter estimation analysis to measure the dephasing
parameters $\delta\Phi_{k}$ and quantify the evidence for the presence of an
environment. We compute the Bayes factor $\mathcal{B}_{\rm vac}^{\rm env}$ to
compare the two hypotheses: (i) the data $d$ is described by the environmental
model ${\mathcal{H}}_{\rm env}$ that allows nonzero values of the dephasing
parameter $\delta\Phi_{k}$, versus (ii) $d$ is described by the vacuum model
${\mathcal{H}}_{\rm vac}$ where no additional parameters are required. We
estimate the magnitude of the environmental effects by computing the
marginalized posterior probability distribution of the $\delta\Phi_{k}$
parameter within ${\mathcal{H}}_{\rm env}$. The posterior distribution is
obtained using Bayes’ rule for a given prior distribution of the model
parameters. The hypotheses share the same prior distribution for their common
parameters, and we consider a (uninformative) zero-centred uniform prior
distribution for $\delta\Phi_{k}$ in ${\mathcal{H}}_{\rm env}$. As we
discussed in the Supplemental Material, the mismatch ($\mathpzc{MM}$) between
vacuum waveform ($h_{\rm vac}$) and environment waveform ($h_{\rm env}$) shows
that the sensitivity of the dephasing parameters for low-mass systems is
significantly higher than for a high-mass system—a suitable prior range of
$\delta\Phi_{k}$ is crucial to ensure the sampler convergence on the global
maximum. Throughout this work, we choose to vary a single phase parameter at a
time to address specific environmental effects leading to phase deformation.

Density constraints from GWTC-1. We analyze individually each of the eleven
events of GWTC-1 to measure the evidence for environmental effects in the data
[95, 96, 97] from the LIGO [98] and Virgo [99] detectors. Table 2 shows the
values of $\log_{10}\mathcal{B}_{\rm vac}^{\rm env}$ for environmental effects
with $k=\\{-9,-11\\}$ relative to the vacuum model. We find negative values of
the log Bayes factor for all events, except for GW151012. Since this event was
observed with a low statistical significance [100, 1] and the Bayes factor
that we found is relatively small ($\log_{10}\mathcal{B}_{\rm vac}^{\rm
env}<1$), this result is indecisive. Therefore, we do not find any support for
environmental effects in the data of GWTC-1. On the other hand, we can not
derive any evidence for the vacuum model either. Environmental corrections are
more effective in the early inspiral phase of the coalescence, where LIGO-like
detectors might be less sensitive. Even if a binary evolves in a medium, the
detection of the environmental effect might be challenging because of noise
domination. Nevertheless, we can use the results of this analysis to estimate
the upper limits on the environment density using the posterior samples of
$\delta\Phi_{k}$, together with the ones of the binary masses.

The environmental effects are responsible for negative values of
$\delta\Phi_{k}$, but the sampler is not limited to those as we also intend to
conduct a model-independent (agnostic) vacuum GR test (see Supplemental
Material). To put upper bounds on the density, we select only the samples with
negative $\delta\Phi_{k}$ values, and calculate the density for the individual
environmental effects using $\widehat{\delta\Phi}_{k}=\delta\Phi_{k}$ and the
expressions in Tab. 1 222This procedure can have a non-trivial prior effect,
we re-weight the posterior to enforce a uniform density prior. We do not find
any noticeable effect in this process.. Figure 1 shows the 90% upper bound on
the density posterior obtained considering specific environmental effects. Our
results for the events GW150914, GW151226, GW170104, GW170608, GW170814, and
GW170817 show the density constraint
$\tilde{\rho}\lesssim(21\textup{--}2\times
10^{6})\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$, which decisively rules out
the binary formation scenario of dynamical fragmentation [24, 50]. The
remaining events’ low inspiral SNR with fewer inspiral cycles in the band led
to poor density constraints—inconclusive for that scenario. We find the
notable bound $\tilde{\rho}\lesssim 21\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$
from GW170817, which is the foremost constraint from GWTC-1 (roughly the
density of gold at room temperature on Earth).

Figure 1: 90% upper bounds on the environmental density obtained considering
the effect of CA (blue triangle), BHLA (red cross), and DF (green diamond).
Figure 2: Logarithmic Bayes factor as function of the medium density. This
shows the detectability of different environmental effects using the aLIGO
sensitivity curve, for GW150914, GW170608, and GW170817-like injections in
zero noise realization. We omit the BHLA curve since it follows very closely
the DF one.

Forecasts for future observations. To understand the detectability of
environmental effects with future detectors, we carry out a set of injection
analyses mimicking events like GW150914, GW170608, and GW170817 using the
design sensitivity of advanced LIGO (aLIGO) [102], and perform a mismatch-
based prediction for Einstein Telescope (ET) [103] and B-DECIGO [104, 105],
with the corresponding lower cutoff frequencies set to $15$ Hz, $5$ Hz, and
$0.1$ Hz, respectively. We select the maximum likelihood sample from the
standard vacuum model analysis on real data from the Bilby GWTC-1 rerun [106,
107]. Then, we add a set of $\delta\Phi_{k}$ values to manipulate the
individual environmental effects with a monotonically increasing density
parameter.

For aLIGO, we conduct Bayesian analysis with zero-noise realization
considering a fixed injection SNR ($\rho_{\rm inj}$) of 25 and compare the
environment versus vacuum models. Figure 2 shows the values of
$\log_{10}\mathcal{B}_{\rm vac}^{\rm env}$ for a range of medium densities
$0.1$–$10^{8}\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$. The effect of DF in
GW170817, GW170608, GW150914-like events is detectable when
$\tilde{\rho}\gtrsim(10,\,4\times
10^{4},\,10^{6})\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$, respectively. The
effect of CA is detectable only for media $\sim 10$–$100$ times denser. When
the environment density is set to be half of the (threshold) detectable value,
it can still be measurable even though it is not detectable. This is because
in comparing environment versus vacuum, the environment hypothesis is hit by
Occam’s razor factor due to the extra dephasing parameter in
${\mathcal{H}}_{\rm env}$. That also explains why almost all
$\log_{10}\mathcal{B}_{\rm vac}^{\rm env}$ values in Tab. 2 are negative.
Figure 3 shows that, for the case of the non-vacuum GW170817-like injection,
the vacuum model analysis leads to an overestimation of the chirp mass and
effective spin, while it underestimates the mass ratio. A recent study [108]
also demonstrated the detectability of environmental effects and the biases on
the inferred parameters of the binary through the study of the (multi-modal)
ringdown phase. The use of environmental waveform models in future parameter
estimation analyses could potentially mitigate this bias.

Figure 3: Posterior distribution of the chirp mass ($\mathcal{M}_{c}$), mass
ratio ($q$), and effective spin ($\chi_{\rm eff}$) using the models
${\mathcal{H}}_{\rm vac}$ and ${\mathcal{H}}_{\rm env}$ with nonzero
$\delta\Phi_{-11}$. We injected GW170817-like waveforms deformed by the effect
of DF for a set of environment densities (shown in bottom $x$-axis), using the
aLIGO sensitivity curve. The top $x$-axis shows the logarithmic Bayes factor
values for ${\mathcal{H}}_{\rm vac}$ versus ${\mathcal{H}}_{\rm env}$. Small
red marks indicate the injection parameter values.

To study the observability of environments with future detectors like ET and
B-DECIGO we use an alternative approach. The Bayes factor between two
competitive models can be approximated as $\ln{\mathcal{B}_{\rm vac}^{\rm
env}}\approx\rho^{2}_{\rm inj}(1-\mathpzc{FF}^{2})/2$ when ${\mathcal{H}}_{\rm
env}$ is true [109], where $\mathpzc{FF}$ refers to the fitting factor of the
injected waveform ($h_{\rm inj}\in{\mathcal{H}}_{\rm env}$) with the vacuum
waveform, such that the match between $h_{\rm vac}$ and $h_{\rm inj}$ is
maximized over all the parameters in ${\mathcal{H}}_{\rm vac}$. We assume
$\mathpzc{FF}\approx 1-\mathpzc{MM}$ ignoring the correlation between
$\delta\Phi_{k}$ and the physical parameters in ${\mathcal{H}}_{\rm vac}$.
Finally, $\mathpzc{MM}$ is a function only of the dephasing parameter
$\delta\Phi_{k}$, which in turn is determined by the environment density and
the binary parameters through the expressions in Tab. 1. Figure 4 shows the
curves of injected SNR versus environment density necessary to achieve
$\log_{10}\mathcal{B}_{\rm vac}^{\rm env}=3$ for events like GW170817 and
GW170608. Our results suggest that ET will be sensitive to the effect of DF in
a GW170817-like event for an environment with $\tilde{\rho}\gtrsim
10^{-3}\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$ and to CA effects in a medium
$\sim 10^{3}$ times denser. As could be already anticipated, due to its better
low frequency coverage, B-DECIGO will be sensitive to much lower densities; DF
effects in a GW170817-like event will be detectable for an environment with
$\tilde{\rho}\gtrsim 10^{-12}\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$ and CA
effects for a density $\sim 10^{4}$ times larger.

A forecast based on the Fisher information matrix was previously performed in
Ref. [34]. Overall, our results are more pessimistic (by a few orders of
magnitude) than the Fisher analysis: we find a threshold density for
detectability larger by roughly 2 orders of magnitude for aLIGO and ET; the
difference is larger for DECIGO, because we consider the 2030s B-DECIGO.

Figure 4: Curves of required SNR for a given density value to achieve
$\log_{10}\mathcal{B}_{\rm vac}^{\rm env}=3$ for a specific environmental
effect, in the configuration of the third-generation detector ET (cyan shade)
and the Japanese space-detector B-DECIGO (gray shade). The dots represent the
expected SNR if we replace the LIGO-Hanford detector with those future
detectors. We omit the BHLA curve since it follows very closely the DF one.

Discussion and Conclusions. We have designed a model-agnostic Bayesian
analysis for detecting the existence of environments surrounding compact
binaries. For that, we modelled the waveform considering environmental
corrections at $-4.5$PN and $-5.5$PN associated to the leading effect of
accretion and DF at linear order in $\epsilon\coloneqq\tilde{\rho}M_{z}^{2}\ll
1$. We checked that higher order terms in $\epsilon$, which affect the phasing
at more negative PN orders, are indeed sub-leading and do not alter our
results for the systems we studied (see Supplemental Material). We analyzed
individually each event in GWTC-1 and found _no_ evidence for the presence of
environments. In addition, we used our Bayesian analysis to derive upper
bounds on the densities of GWTC-1 binary surroundings. Our results for
GW150914, GW151226, GW170104, GW170608, GW170814, and GW170817 decisively
rules out the scenario of dynamical fragmentation of massive stars as a
possible formation channel for these events [24]. Our results for GW150914 are
also in excellent agreement with the findings of Ref. [50], which used
numerical relativity waveforms in gaseous environments.

Even if a compact binary is found evolving in a (low-density) medium, the
detection of environmental effects with LIGO-like detectors is challenging,
since such detectors are not very sensitive to the early inspiral stage, and
the extra dephasing parameter incurs the Occam’s razor penalty, resulting in
less evidence for the environment model. We found that a medium density can
actually be measurable, even when its value is half of the threshold for
detectability ($\log_{10}\mathcal{B}_{\rm vac}^{\rm env}\gtrsim 0$). We also
show that environmental effects can substantially bias the recovered
parameters of the vacuum model, even when they are not detectable. Finally,
our zero-noise injection analyses indicate that LIGO-like detectors will be
capable to see environments (at least) as dilute as $\tilde{\rho}\sim
10\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$—roughly the density of lead at room
temperature on Earth. While they may in the future definitively exclude the
dynamical fragmentation scenario, environments like accretion disk and dark
matter overdensities seem to be out of reach.

We have also analysed the prospects for future (2030s) detectors like ET and
B-DECIGO. We demonstrated that ET will be sensitive to the effects of DF and
BHLA for media as dilute as $\tilde{\rho}\sim
10^{-3}\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$, implying that it may be
capable of detecting environmental effects on a compact binary merging within
dense (thin) accretion disks or superradiant clouds hosted by a SMBH. Our
results are even more promising for B-DECIGO, indicating that it will be
sensitive to DF and BHLA effects from environments as dilute as
$\tilde{\rho}\sim 10^{-12}\>$\mathrm{g}\text{/}{\mathrm{cm}}^{3}$$, which
covers, e.g., most accretion disk densities, superradiant clouds of ultra-
light bosons, or cold dark matter spikes.

In this work we performed a model-agnostic analysis since our goal was to
assess the overall capability of LIGO-like and near-future 2030s detectors in
probing environments and to derive the first (order of magnitude) constraints
on environments from current observations. This also justifies our
(unphysical) choice for considering particular environmental effects
separately. In future work, we plan to focus on specific environments (like
accretion disks and superradiant clouds), considering simultaneously several
environmental effects entering at different PNs, and to study the
distinguishability between different environments with compact binaries (see,
e.g., the study of Ref. [41] for EMRIs). We also plan to extend our modelling
to asymmetric binaries. As the events in GWTC-1 are nearly equal mass
binaries, we analyzed them using a quadrupolar waveform model. To analyze the
events in GWTC-2/3 which contain many asymmetric binaries, we need to include
the contribution of higher multipoles, which become more relevant with
increasing mass ratio. In a follow-up work, we will focus on building an
environmental waveform model including higher harmonics.

###### Acknowledgements.

We gratefully acknowledge comments and feedback from Nathan Johnson-McDaniel
and Stefano Rinaldi. We are grateful to Marc Andrés-Carcasona for useful
discussions throughout the preparation of the paper. This research has made
use of data or software obtained from the Gravitational Wave Open Science
Center (gwosc.org), a service of the LIGO Scientific Collaboration, the Virgo
Collaboration, and KAGRA. This material is based upon work supported by NSF’s
LIGO Laboratory which is a major facility fully funded by the National Science
Foundation, as well as the Science and Technology Facilities Council (STFC) of
the United Kingdom, the Max-Planck-Society (MPS), and the State of
Niedersachsen/Germany for support of the construction of Advanced LIGO and
construction and operation of the GEO600 detector. Additional support for
Advanced LIGO was provided by the Australian Research Council. Virgo is
funded, through the European Gravitational Observatory (EGO), by the French
Centre National de Recherche Scientifique (CNRS), the Italian Istituto
Nazionale di Fisica Nucleare (INFN) and the Dutch Nikhef, with contributions
by institutions from Belgium, Germany, Greece, Hungary, Ireland, Japan,
Monaco, Poland, Portugal, Spain. KAGRA is supported by Ministry of Education,
Culture, Sports, Science and Technology (MEXT), Japan Society for the
Promotion of Science (JSPS) in Japan; National Research Foundation (NRF) and
Ministry of Science and ICT (MSIT) in Korea; Academia Sinica (AS) and National
Science and Technology Council (NSTC) in Taiwan. The authors are grateful for
computational resources provided by the LIGO Laboratory and supported by
National Science Foundation Grants PHY-0757058 and PHY-0823459. S.R. and M.H.
were supported by the research program of the Netherlands Organization for
Scientific Research (NWO). R.V. is supported by grant no. FJC2021-046551-I
funded by MCIN/AEI/10.13039/501100011033 and by the European Union
NextGenerationEU/PRTR. R.V. also acknowledges support by grant no. CERN/FIS-
PAR/0023/2019. This paper has been given LIGO DCC number LIGO-P2300301. This
work is partially supported by the Spanish MCIN/AEI/ 10.13039/501100011033
under the grants SEV-2016-0588, PGC2018-101858-B-I00, and PID2020-113701GB-I00
some of which include ERDF funds from the European Union. IFAE is partially
funded by the CERCA program of the Generalitat de Catalunya. To generate the
waveforms, we have used the LALSimulation package of the LIGO Algorithms
Library (LAL) software suite [110]. To perform the Bayesian parameter
estimation analyses, we have used the nested sampling algorithm [111, 112,
113] implemented in the LALInference [114] and Bilby [115] packages to
evaluate the integral over the model parameter space and calculate the Bayes
factor.

## References

  * Abbott  _et al._ [2019a] B. P. Abbott _et al._ (LIGO Scientific, Virgo), GWTC-1: A Gravitational-Wave Transient Catalog of Compact Binary Mergers Observed by LIGO and Virgo during the First and Second Observing Runs, Phys. Rev. X 9, 031040 (2019a), arXiv:1811.12907 [astro-ph.HE] .
  * Abbott  _et al._ [2021a] R. Abbott _et al._ (LIGO Scientific, Virgo), GWTC-2: Compact Binary Coalescences Observed by LIGO and Virgo During the First Half of the Third Observing Run, Phys. Rev. X 11, 021053 (2021a), arXiv:2010.14527 [gr-qc] .
  * Abbott  _et al._ [2021b] R. Abbott _et al._ (LIGO Scientific, VIRGO, KAGRA), GWTC-3: Compact Binary Coalescences Observed by LIGO and Virgo During the Second Part of the Third Observing Run, arXiv e-prints , arXiv:2111.03606 (2021b), arXiv:2111.03606 [gr-qc] .
  * Abbott  _et al._ [2019b] B. P. Abbott _et al._ (LIGO Scientific, Virgo), Binary Black Hole Population Properties Inferred from the First and Second Observing Runs of Advanced LIGO and Advanced Virgo, Astrophys. J. Lett. 882, L24 (2019b), arXiv:1811.12940 [astro-ph.HE] .
  * Abbott  _et al._ [2021c] R. Abbott _et al._ (LIGO Scientific, Virgo), Population Properties of Compact Objects from the Second LIGO-Virgo Gravitational-Wave Transient Catalog, Astrophys. J. Lett. 913, L7 (2021c), arXiv:2010.14533 [astro-ph.HE] .
  * Abbott  _et al._ [2023a] R. Abbott _et al._ (KAGRA, VIRGO, LIGO Scientific), Population of Merging Compact Binaries Inferred Using Gravitational Waves through GWTC-3, Phys. Rev. X 13, 011048 (2023a), arXiv:2111.03634 [astro-ph.HE] .
  * Abbott  _et al._ [2017] B. P. Abbott _et al._ (LIGO Scientific, Virgo, 1M2H, Dark Energy Camera GW-E, DES, DLT40, Las Cumbres Observatory, VINROUGE, MASTER), A gravitational-wave standard siren measurement of the Hubble constant, Nature 551, 85 (2017), arXiv:1710.05835 [astro-ph.CO] .
  * Abbott  _et al._ [2021d] B. P. Abbott _et al._ (LIGO Scientific, Virgo, VIRGO), A Gravitational-wave Measurement of the Hubble Constant Following the Second Observing Run of Advanced LIGO and Virgo, Astrophys. J. 909, 218 (2021d), arXiv:1908.06060 [astro-ph.CO] .
  * Abbott  _et al._ [2023b] R. Abbott _et al._ (LIGO Scientific, Virgo,, KAGRA, VIRGO), Constraints on the Cosmic Expansion History from GWTC–3, Astrophys. J. 949, 76 (2023b), arXiv:2111.03604 [astro-ph.CO] .
  * Abbott  _et al._ [2021e] R. Abbott _et al._ (LIGO Scientific, VIRGO), GWTC-2.1: Deep Extended Catalog of Compact Binary Coalescences Observed by LIGO and Virgo During the First Half of the Third Observing Run, arXiv e-prints , arXiv:2108.01045 (2021e), arXiv:2108.01045 [gr-qc] .
  * Abbott  _et al._ [2019c] B. P. Abbott _et al._ (LIGO Scientific, Virgo), Tests of General Relativity with the Binary Black Hole Signals from the LIGO-Virgo Catalog GWTC-1, Phys. Rev. D 100, 104036 (2019c), arXiv:1903.04467 [gr-qc] .
  * Abbott  _et al._ [2021f] R. Abbott _et al._ (LIGO Scientific, Virgo), Tests of general relativity with binary black holes from the second LIGO-Virgo gravitational-wave transient catalog, Phys. Rev. D 103, 122002 (2021f), arXiv:2010.14529 [gr-qc] .
  * Abbott  _et al._ [2021g] R. Abbott _et al._ (LIGO Scientific, VIRGO, KAGRA), Tests of General Relativity with GWTC-3, arXiv e-prints , arXiv:2112.06861 (2021g), arXiv:2112.06861 [gr-qc] .
  * Abbott  _et al._ [2018] B. P. Abbott _et al._ (KAGRA, LIGO Scientific, Virgo, VIRGO), Prospects for observing and localizing gravitational-wave transients with Advanced LIGO, Advanced Virgo and KAGRA, Living Rev. Rel. 21, 3 (2018), arXiv:1304.0670 [gr-qc] .
  * Banerjee  _et al._ [2010] S. Banerjee, H. Baumgardt, and P. Kroupa, Stellar-mass black holes in star clusters: implications for gravitational wave radiation, Mon. Not. Roy. Astron. Soc. 402, 371 (2010), arXiv:0910.3954 [astro-ph.SR] .
  * Chatterjee  _et al._ [2017] S. Chatterjee, C. L. Rodriguez, and F. A. Rasio, Binary Black Holes in Dense Star Clusters: Exploring the Theoretical Uncertainties, Astrophys. J. 834, 68 (2017), arXiv:1603.00884 [astro-ph.GA] .
  * Ziosi  _et al._ [2014] B. M. Ziosi, M. Mapelli, M. Branchesi, and G. Tormen, Dynamics of stellar black holes in young star clusters with different metallicities – II. Black hole–black hole binaries, Mon. Not. Roy. Astron. Soc. 441, 3703 (2014), arXiv:1404.7147 [astro-ph.GA] .
  * Rodriguez  _et al._ [2015] C. L. Rodriguez, M. Morscher, B. Pattabiraman, S. Chatterjee, C.-J. Haster, and F. A. Rasio, Binary Black Hole Mergers from Globular Clusters: Implications for Advanced LIGO, Phys. Rev. Lett. 115, 051101 (2015), [Erratum: Phys.Rev.Lett. 116, 029901 (2016)], arXiv:1505.00792 [astro-ph.HE] .
  * Mapelli  _et al._ [2022] M. Mapelli, Y. Bouffanais, F. Santoliquido, M. A. Sedda, and M. C. Artale, The cosmic evolution of binary black holes in young, globular, and nuclear star clusters: rates, masses, spins, and mixing fractions, Mon. Not. Roy. Astron. Soc. 511, 5797 (2022), arXiv:2109.06222 [astro-ph.HE] .
  * Mapelli [2021] M. Mapelli, Formation Channels of Single and Binary Stellar-Mass Black Holes, in  _Handbook of Gravitational Wave Astronomy_ (2021) p. 16, arXiv:2106.00699 [astro-ph.HE] .
  * Sedda  _et al._ [2023] M. A. Sedda, S. Naoz, and B. Kocsis, Quiescent and Active Galactic Nuclei as Factories of Merging Compact Objects in the Era of Gravitational Wave Astronomy, Universe 9, 138 (2023), arXiv:2302.14071 [astro-ph.GA] .
  * Rowan  _et al._ [2023] C. Rowan, T. Boekholt, B. Kocsis, and Z. Haiman, Black hole binary formation in AGN discs: from isolation to merger, Mon. Not. Roy. Astron. Soc. 524, 2770 (2023), arXiv:2212.06133 [astro-ph.GA] .
  * Ishibashi and Gröbner [2020] W. Ishibashi and M. Gröbner, Evolution of binary black holes in AGN accretion discs: Disc-binary interaction and gravitational wave emission, Astron. Astrophys. 639, A108 (2020), arXiv:2006.07407 [astro-ph.GA] .
  * Loeb [2016] A. Loeb, Electromagnetic Counterparts to Black Hole Mergers Detected by LIGO, Astrophys. J. Lett. 819, L21 (2016), arXiv:1602.04735 [astro-ph.HE] .
  * McKernan  _et al._ [2019] B. McKernan, K. E. S. Ford, I. Bartos, M. J. Graham, W. Lyra, S. Marka, Z. Marka, N. P. Ross, D. Stern, and Y. Yang, Ram-pressure stripping of a kicked Hill sphere: Prompt electromagnetic emission from the merger of stellar mass black holes in an AGN accretion disk, Astrophys. J. Lett. 884, L50 (2019), arXiv:1907.03746 [astro-ph.HE] .
  * Perna  _et al._ [2019] R. Perna, D. Lazzati, and W. Farr, Limits on electromagnetic counterparts of gravitational wave-detected binary black hole mergers, Astrophys. J. 875, 49 (2019), arXiv:1901.04522 [astro-ph.HE] .
  * Graham  _et al._ [2020] M. J. Graham _et al._ , Candidate Electromagnetic Counterpart to the Binary Black Hole Merger Gravitational Wave Event S190521g, Phys. Rev. Lett. 124, 251102 (2020), arXiv:2006.14122 [astro-ph.HE] .
  * Kocsis  _et al._ [2011] B. Kocsis, N. Yunes, and A. Loeb, Observable Signatures of EMRI Black Hole Binaries Embedded in Thin Accretion Disks, Phys. Rev. D 84, 024032 (2011), arXiv:1104.2322 [astro-ph.GA] .
  * Yunes  _et al._ [2011] N. Yunes, B. Kocsis, A. Loeb, and Z. Haiman, Imprint of Accretion Disk-Induced Migration on Gravitational Waves from Extreme Mass Ratio Inspirals, Phys. Rev. Lett. 107, 171103 (2011), arXiv:1103.4609 [astro-ph.CO] .
  * Barausse  _et al._ [2014] E. Barausse, V. Cardoso, and P. Pani, Can environmental effects spoil precision gravitational-wave astrophysics?, Phys. Rev. D 89, 104059 (2014), arXiv:1404.7149 [gr-qc] .
  * Barausse  _et al._ [2015] E. Barausse, V. Cardoso, and P. Pani, Environmental Effects for Gravitational-wave Astrophysics, J. Phys. Conf. Ser. 610, 012044 (2015), arXiv:1404.7140 [astro-ph.CO] .
  * Tamanini  _et al._ [2020] N. Tamanini, A. Klein, C. Bonvin, E. Barausse, and C. Caprini, Peculiar acceleration of stellar-origin black hole binaries: Measurement and biases with LISA, Phys. Rev. D 101, 063002 (2020), arXiv:1907.02018 [astro-ph.IM] .
  * Toubiana  _et al._ [2021] A. Toubiana _et al._ , Detectable environmental effects in GW190521-like black-hole binaries with LISA, Phys. Rev. Lett. 126, 101105 (2021), arXiv:2010.06056 [astro-ph.HE] .
  * Cardoso and Maselli [2020] V. Cardoso and A. Maselli, Constraints on the astrophysical environment of binaries with gravitational-wave observations, Astron. Astrophys. 644, A147 (2020), arXiv:1909.05870 [astro-ph.HE] .
  * Derdzinski  _et al._ [2021] A. Derdzinski, D. D’Orazio, P. Duffell, Z. Haiman, and A. MacFadyen, Evolution of gas disc–embedded intermediate mass ratio inspirals in the $LISA$ band, Mon. Not. Roy. Astron. Soc. 501, 3540 (2021), arXiv:2005.11333 [astro-ph.HE] .
  * Zwick  _et al._ [2022] L. Zwick, A. Derdzinski, M. Garg, P. R. Capelo, and L. Mayer, Dirty waveforms: multiband harmonic content of gas-embedded gravitational wave sources, Mon. Not. Roy. Astron. Soc. 511, 6143 (2022), arXiv:2110.09097 [astro-ph.HE] .
  * Kuntz [2022] A. Kuntz, Precession resonances in hierarchical triple systems, Phys. Rev. D 105, 024017 (2022), arXiv:2112.05167 [gr-qc] .
  * Sberna  _et al._ [2022] L. Sberna _et al._ , Observing GW190521-like binary black holes and their environment with LISA, Phys. Rev. D 106, 064056 (2022), arXiv:2205.08550 [gr-qc] .
  * Speri  _et al._ [2023] L. Speri, A. Antonelli, L. Sberna, S. Babak, E. Barausse, J. R. Gair, and M. L. Katz, Probing Accretion Physics with Gravitational Waves, Phys. Rev. X 13, 021035 (2023), arXiv:2207.10086 [gr-qc] .
  * Cardoso  _et al._ [2022] V. Cardoso, K. Destounis, F. Duque, R. Panosso Macedo, and A. Maselli, Gravitational Waves from Extreme-Mass-Ratio Systems in Astrophysical Environments, Phys. Rev. Lett. 129, 241103 (2022), arXiv:2210.01133 [gr-qc] .
  * Cole  _et al._ [2022] P. S. Cole, G. Bertone, A. Coogan, D. Gaggero, T. Karydas, B. J. Kavanagh, T. F. M. Spieksma, and G. M. Tomaselli, Disks, spikes, and clouds: distinguishing environmental effects on BBH gravitational waveforms, arXiv e-prints , arXiv:2211.01362 (2022), arXiv:2211.01362 [gr-qc] .
  * Zwick  _et al._ [2023] L. Zwick, P. R. Capelo, and L. Mayer, Priorities in gravitational waveforms for future space-borne detectors: vacuum accuracy or environment?, Mon. Not. Roy. Astron. Soc. 521, 4645 (2023), arXiv:2209.04060 [gr-qc] .
  * Speeney  _et al._ [2022] N. Speeney, A. Antonelli, V. Baibhav, and E. Berti, Impact of relativistic corrections on the detectability of dark-matter spikes with gravitational waves, Phys. Rev. D 106, 044027 (2022), arXiv:2204.12508 [gr-qc] .
  * Kuntz  _et al._ [2023] A. Kuntz, F. Serra, and E. Trincherini, Effective two-body approach to the hierarchical three-body problem: Quadrupole to 1PN, Phys. Rev. D 107, 044011 (2023), arXiv:2210.13493 [gr-qc] .
  * Kuntz and Leyde [2023] A. Kuntz and K. Leyde, Transverse Doppler effect and parameter estimation of LISA three-body systems, Phys. Rev. D 108, 024002 (2023), arXiv:2212.09753 [gr-qc] .
  * Boudon  _et al._ [2023] A. Boudon, P. Brax, P. Valageas, and L. K. Wong, Gravitational waves from binary black holes in a self-interacting scalar dark matter cloud, arXiv e-prints , arXiv:2305.18540 (2023), arXiv:2305.18540 [astro-ph.CO] .
  * Vijaykumar  _et al._ [2023] A. Vijaykumar, A. Tiwari, S. J. Kapadia, K. G. Arun, and P. Ajith, Waltzing Binaries: Probing the Line-of-sight Acceleration of Merging Compact Objects with Gravitational Waves, Astrophys. J. 954, 105 (2023), arXiv:2302.09651 [astro-ph.HE] .
  * Amaro-Seoane  _et al._ [2017] P. Amaro-Seoane _et al._ , Laser Interferometer Space Antenna, arXiv e-prints (2017), arXiv:1702.00786 [astro-ph.IM] .
  * LIS [2021] LISA—Mission Summary (2021), https://sci.esa.int/s/w5qyMBw [Accessed: 12/05/2023].
  * Fedrow  _et al._ [2017] J. M. Fedrow, C. D. Ott, U. Sperhake, J. Blackman, R. Haas, C. Reisswig, and A. De Felice, Gravitational Waves from Binary Black Hole Mergers inside Stars, Phys. Rev. Lett. 119, 171103 (2017), arXiv:1704.07383 [astro-ph.HE] .
  * Narola  _et al._ [2023] H. Narola, S. Roy, and A. S. Sengupta, Beyond general relativity: Designing a template-based search for exotic gravitational wave signals, Phys. Rev. D 107, 024017 (2023), arXiv:2207.10410 [gr-qc] .
  * Gondolo and Silk [1999] P. Gondolo and J. Silk, Dark matter annihilation at the galactic center, Phys. Rev. Lett. 83, 1719 (1999), arXiv:astro-ph/9906391 .
  * Sadeghian  _et al._ [2013] L. Sadeghian, F. Ferrer, and C. M. Will, Dark matter distributions around massive black holes: A general relativistic analysis, Phys. Rev. D 88, 063522 (2013), arXiv:1305.2619 [astro-ph.GA] .
  * Brito  _et al._ [2015] R. Brito, V. Cardoso, and P. Pani, Superradiance: New Frontiers in Black Hole Physics, Lect. Notes Phys. 906, pp.1 (2015), arXiv:1501.06570 [gr-qc] .
  * Bondi and Hoyle [1944] H. Bondi and F. Hoyle, On the mechanism of accretion by stars, Mon. Not. Roy. Astron. Soc. 104, 273 (1944).
  * Bondi [1952] H. Bondi, On spherically symmetrical accretion, Mon. Not. Roy. Astron. Soc. 112, 195 (1952).
  * Petrich  _et al._ [1989] L. I. Petrich, S. L. Shapiro, R. F. Stark, and S. A. Teukolsky, Accretion onto a Moving Black Hole: A Fully Relativistic Treatment, Astrophys. J. 336, 313 (1989).
  * Chandrasekhar [1943] S. Chandrasekhar, Dynamical Friction. I. General Considerations: the Coefficient of Dynamical Friction, Astrophys. J. 97, 255 (1943).
  * Ruderman and Spiegel [1971] M. A. Ruderman and E. A. Spiegel, Galactic Wakes, Astrophys. J. 165, 1 (1971).
  * Rephaeli and Salpeter [1980] Y. Rephaeli and E. E. Salpeter, Flow past a massive object and the gravitational drag, Astrophys. J. 240, 20 (1980).
  * Tremaine and Weinberg [1984] S. Tremaine and M. D. Weinberg, Dynamical friction in spherical systems., Mon. Not. Roy. Astron. Soc. 209, 729 (1984).
  * Ostriker [1999] E. C. Ostriker, Dynamical friction in a gaseous medium, Astrophys. J. 513, 252 (1999), arXiv:astro-ph/9810324 [astro-ph] .
  * Barausse [2007] E. Barausse, Relativistic dynamical friction in a collisional fluid, Mon. Not. Roy. Astron. Soc. 382, 826 (2007), arXiv:0709.0211 [astro-ph] .
  * Cardoso  _et al._ [2021] V. Cardoso, C. F. B. Macedo, and R. Vicente, Eccentricity evolution of compact binaries and applications to gravitational-wave physics, Phys. Rev. D 103, 023015 (2021), arXiv:2010.15151 [gr-qc] .
  * Peters [1964] P. C. Peters, Gravitational Radiation and the Motion of Two Point Masses, Phys. Rev. 136, B1224 (1964).
  * Desjacques  _et al._ [2022] V. Desjacques, A. Nusser, and R. Bühler, Analytic Solution to the Dynamical Friction Acting on Circularly Moving Perturbers, Astrophys. J. 928, 64 (2022), arXiv:2111.07366 [astro-ph.GA] .
  * Baumann  _et al._ [2022] D. Baumann, G. Bertone, J. Stout, and G. M. Tomaselli, Ionization of gravitational atoms, Phys. Rev. D 105, 115036 (2022), arXiv:2112.14777 [gr-qc] .
  * Buehler and Desjacques [2023] R. Buehler and V. Desjacques, Dynamical friction in fuzzy dark matter: Circular orbits, Phys. Rev. D 107, 023516 (2023), arXiv:2207.13740 [astro-ph.CO] .
  * Vicente and Cardoso [2022] R. Vicente and V. Cardoso, Dynamical friction of black holes in ultralight dark matter, Phys. Rev. D 105, 083008 (2022), arXiv:2201.08854 [gr-qc] .
  * Traykova  _et al._ [2023] D. Traykova, R. Vicente, K. Clough, T. Helfer, E. Berti, P. G. Ferreira, and L. Hui, Relativistic drag forces on black holes from scalar dark matter clouds of all sizes, arXiv e-prints , arXiv:2305.10492 (2023), arXiv:2305.10492 [gr-qc] .
  * Tomaselli  _et al._ [2023] G. M. Tomaselli, T. F. M. Spieksma, and G. Bertone, Dynamical friction in gravitational atoms, JCAP 07, 070, arXiv:2305.15460 [gr-qc] .
  * Kavanagh  _et al._ [2020] B. J. Kavanagh, D. A. Nichols, G. Bertone, and D. Gaggero, Detecting dark matter around black holes with gravitational waves: Effects of dark-matter dynamics on the gravitational waveform, Phys. Rev. D 102, 083006 (2020), arXiv:2002.12811 [gr-qc] .
  * Note [1] The transfer of momentum via accretion is also responsible for an additional (hydrodynamic) drag force which may be included within the factor $\mathcal{I}$ in Eq. (1) [57, Barausse:2007dy, 70].
  * Farris  _et al._ [2010] B. D. Farris, Y. T. Liu, and S. L. Shapiro, Binary Black Hole Mergers in Gaseous Environments: ’Binary Bondi’ and ’Binary Bondi-Hoyle-Lyttleton’ Accretion, Phys. Rev. D 81, 084008 (2010), arXiv:0912.2096 [astro-ph.HE] .
  * Khan  _et al._ [2018] A. Khan, V. Paschalidis, M. Ruiz, and S. L. Shapiro, Disks Around Merging Binary Black Holes: From GW150914 to Supermassive Black Holes, Phys. Rev. D 97, 044036 (2018), arXiv:1801.02624 [astro-ph.HE] .
  * Comerford  _et al._ [2019] T. A. F. Comerford, R. G. Izzard, R. A. Booth, and G. Rosotti, Bondi-Hoyle-Lyttleton accretion by binary stars, Mon. Not. Roy. Astron. Soc. 490, 5196 (2019), arXiv:1910.13353 [astro-ph.SR] .
  * Antoni  _et al._ [2019] A. Antoni, M. MacLeod, and E. Ramirez-Ruiz, The Evolution of Binaries in a Gaseous Medium: Three-Dimensional Simulations of Binary Bondi-Hoyle-Lyttleton Accretion, Astrophys. J. 884, 22 (2019), arXiv:1901.07572 [astro-ph.HE] .
  * Kaaz  _et al._ [2023] N. Kaaz, S. L. Schrøder, J. J. Andrews, A. Antoni, and E. Ramirez-Ruiz, The Hydrodynamic Evolution of Binary Black Holes Embedded within the Vertically Stratified Disks of Active Galactic Nuclei, Astrophys. J. 944, 44 (2023), arXiv:2103.12088 [astro-ph.HE] .
  * Edgar [2004] R. G. Edgar, A Review of Bondi-Hoyle-Lyttleton accretion, New Astron. Rev. 48, 843 (2004), arXiv:astro-ph/0406166 .
  * Eda  _et al._ [2013] K. Eda, Y. Itoh, S. Kuroyanagi, and J. Silk, New Probe of Dark-Matter Properties: Gravitational Waves from an Intermediate-Mass Black Hole Embedded in a Dark-Matter Minispike, Phys. Rev. Lett. 110, 221101 (2013), arXiv:1301.5971 [gr-qc] .
  * Galishnikova  _et al._ [2023] A. Galishnikova, A. Philippov, E. Quataert, F. Bacchini, K. Parfrey, and B. Ripperda, Collisionless Accretion onto Black Holes: Dynamics and Flares, Phys. Rev. Lett. 130, 115201 (2023), arXiv:2212.02583 [astro-ph.HE] .
  * Eddington [1988] A. S. Eddington,  _The Internal Constitution of the Stars_, Cambridge Science Classics (Cambridge University Press, 1988) p. 391.
  * Begelman [1977] M. C. Begelman, Nearly collisionless spherical accretion, Mon. Not. Roy. Astron. Soc. 181, 347 (1977).
  * Barausse  _et al._ [2007] E. Barausse, L. Rezzolla, D. Petroff, and M. Ansorg, Gravitational waves from Extreme Mass Ratio Inspirals in non-pure Kerr spacetimes, Phys. Rev. D 75, 064026 (2007), arXiv:gr-qc/0612123 .
  * Macedo  _et al._ [2013] C. F. B. Macedo, P. Pani, V. Cardoso, and L. C. B. Crispino, Into the lair: gravitational-wave signatures of dark matter, Astrophys. J. 774, 48 (2013), arXiv:1302.2646 [gr-qc] .
  * Cardoso and Macedo [2020] V. Cardoso and C. F. B. Macedo, Drifting through the medium: kicks and self-propulsion of binaries within accretion discs and other environments, Mon. Not. Roy. Astron. Soc. 498, 1963 (2020), arXiv:2008.01091 [astro-ph.HE] .
  * Li  _et al._ [2012] T. G. F. Li, W. Del Pozzo, S. Vitale, C. Van Den Broeck, M. Agathos, J. Veitch, K. Grover, T. Sidery, R. Sturani, and A. Vecchio, Towards a generic test of the strong field dynamics of general relativity using compact binary coalescence, Phys. Rev. D 85, 082003 (2012), arXiv:1110.0530 [gr-qc] .
  * Agathos  _et al._ [2014] M. Agathos, W. Del Pozzo, T. G. F. Li, C. Van Den Broeck, J. Veitch, and S. Vitale, TIGER: A data analysis pipeline for testing the strong-field dynamics of general relativity with gravitational wave signals from coalescing compact binaries, Phys. Rev. D 89, 082001 (2014), arXiv:1311.0420 [gr-qc] .
  * Meidam  _et al._ [2018] J. Meidam _et al._ , Parametrized tests of the strong-field dynamics of general relativity using gravitational wave signals from coalescing binary black holes: Fast likelihood calculations and sensitivity of the method, Phys. Rev. D 97, 044033 (2018), arXiv:1712.08772 [gr-qc] .
  * Husa  _et al._ [2016] S. Husa, S. Khan, M. Hannam, M. Pürrer, F. Ohme, X. Jiménez Forteza, and A. Bohé, Frequency-domain gravitational waves from nonprecessing black-hole binaries. I. New numerical waveforms and anatomy of the signal, Phys. Rev. D 93, 044006 (2016), arXiv:1508.07250 [gr-qc] .
  * Khan  _et al._ [2016] S. Khan, S. Husa, M. Hannam, F. Ohme, M. Pürrer, X. Jiménez Forteza, and A. Bohé, Frequency-domain gravitational waves from nonprecessing black-hole binaries. II. A phenomenological model for the advanced detector era, Phys. Rev. D 93, 044007 (2016), arXiv:1508.07253 [gr-qc] .
  * Hannam  _et al._ [2014] M. Hannam, P. Schmidt, A. Bohé, L. Haegel, S. Husa, F. Ohme, G. Pratten, and M. Pürrer, Simple Model of Complete Precessing Black-Hole-Binary Gravitational Waveforms, Phys. Rev. Lett. 113, 151101 (2014), arXiv:1308.3271 [gr-qc] .
  * Bohe  _et al._ [2016] A. Bohe, M. Hannam, S. Husa, F. Ohme, P. Michael, and P. Schmidt,  _PhenomPv2 – technical notes for the LAL implementation_, LIGO Technical Note T1500602-v4 (LIGO Scientific Collaboration, 2016).
  * Dietrich  _et al._ [2019] T. Dietrich, A. Samajdar, S. Khan, N. K. Johnson-McDaniel, R. Dudi, and W. Tichy, Improving the NRTidal model for binary neutron star systems, Phys. Rev. D 100, 044003 (2019), arXiv:1905.06011 [gr-qc] .
  * Gravitational Wave Open Science Center [2018] Gravitational Wave Open Science Center, Strain data release for gwtc-1: A gravitational-wave transient catalog of compact binary mergers observed by ligo and virgo during the first and second observing runs (2018).
  * Abbott  _et al._ [2021h] R. Abbott _et al._ (LIGO Scientific, Virgo), Open data from the first and second observing runs of Advanced LIGO and Advanced Virgo, SoftwareX 13, 100658 (2021h), arXiv:1912.11716 [gr-qc] .
  * Vallisneri  _et al._ [2015] M. Vallisneri, J. Kanner, R. Williams, A. Weinstein, and B. Stephens, The LIGO Open Science Center, J. Phys. Conf. Ser. 610, 012021 (2015), arXiv:1410.4839 [gr-qc] .
  * Aasi  _et al._ [2015] J. Aasi _et al._ (LIGO Scientific), Advanced LIGO, Class. Quant. Grav. 32, 074001 (2015), arXiv:1411.4547 [gr-qc] .
  * Acernese  _et al._ [2015] F. Acernese _et al._ (VIRGO), Advanced Virgo: a second-generation interferometric gravitational wave detector, Class. Quant. Grav. 32, 024001 (2015), arXiv:1408.3978 [gr-qc] .
  * Abbott  _et al._ [2016] B. P. Abbott _et al._ (LIGO Scientific, Virgo), Binary Black Hole Mergers in the first Advanced LIGO Observing Run, Phys. Rev. X 6, 041015 (2016), [Erratum: Phys.Rev.X 8, 039903 (2018)], arXiv:1606.04856 [gr-qc] .
  * Note [2] This procedure can have a non-trivial prior effect, we re-weight the posterior to enforce a uniform density prior. We do not find any noticeable effect in this process.
  * Barsotti  _et al._ [2018] L. Barsotti, S. Gras, M. Evans, and P. Fritschel,  _The updated Advanced LIGO design curve_, LIGO Technical Note T1800044-v5 (LIGO Scientific Collaboration, 2018) updated from T0900288-v3.
  * Hild  _et al._ [2011] S. Hild _et al._ , Sensitivity Studies for Third-Generation Gravitational Wave Observatories, Class. Quant. Grav. 28, 094013 (2011), arXiv:1012.0908 [gr-qc] .
  * Isoyama  _et al._ [2018] S. Isoyama, H. Nakano, and T. Nakamura, Multiband Gravitational-Wave Astronomy: Observing binary inspirals with a decihertz detector, B-DECIGO, PTEP 2018, 073E01 (2018), arXiv:1802.06977 [gr-qc] .
  * Yagi and Seto [2011] K. Yagi and N. Seto, Detector configuration of DECIGO/BBO and identification of cosmological neutron-star binaries, Phys. Rev. D 83, 044011 (2011), [Erratum: Phys.Rev.D 95, 109901 (2017)], arXiv:1101.3940 [astro-ph.CO] .
  * Romero-Shaw  _et al._ [2020a] I. M. Romero-Shaw _et al._ , Bayesian inference for compact binary coalescences with bilby: validation and application to the first LIGO–Virgo gravitational-wave transient catalogue, Mon. Not. Roy. Astron. Soc. 499, 3295 (2020a), arXiv:2006.00714 [astro-ph.IM] .
  * Romero-Shaw  _et al._ [2020b] I. M. Romero-Shaw _et al._ ,  _Parameter estimation posterior samples for BILBY analysis of GWTC-1_, LIGO Publication Note P2000193-v2 (LIGO Scientific Collaboration, 2020).
  * Leong  _et al._ [2023] S. H. W. Leong, J. Calderón Bustillo, M. Gracia-Linares, and P. Laguna, Impact of ringdown higher-order modes on black-hole mergers in dense environments: the scalar field case, detectability and parameter biases, arXiv e-prints , arXiv:2308.03250 (2023), arXiv:2308.03250 [gr-qc] .
  * Cornish  _et al._ [2011] N. Cornish, L. Sampson, N. Yunes, and F. Pretorius, Gravitational Wave Tests of General Relativity with the Parameterized Post-Einsteinian Framework, Phys. Rev. D 84, 062003 (2011), arXiv:1105.2088 [gr-qc] .
  * LIGO Scientific Collaboration [2023] LIGO Scientific Collaboration, LIGO Algorithm Library - LALSuite, free software (GPL) (2023).
  * Skilling [2006] J. Skilling, Nested sampling for general Bayesian computation, Bayesian Analysis 1, 833 (2006).
  * Veitch and Vecchio [2010] J. Veitch and A. Vecchio, Bayesian coherent analysis of in-spiral gravitational wave signals with a detector network, Phys. Rev. D 81, 062003 (2010), arXiv:0911.3820 [astro-ph.CO] .
  * Speagle [2020] J. S. Speagle, dynesty: a dynamic nested sampling package for estimating bayesian posteriors and evidences, Mon. Not. Roy. Astron. Soc. 493, 3132 (2020).
  * Veitch  _et al._ [2015] J. Veitch _et al._ , Parameter estimation for compact binaries with ground-based gravitational-wave observations using the LALInference software library, Phys. Rev. D 91, 042003 (2015), arXiv:1409.7215 [gr-qc] .
  * Ashton  _et al._ [2019] G. Ashton _et al._ , BILBY: A user-friendly Bayesian inference library for gravitational-wave astronomy, Astrophys. J. Suppl. 241, 27 (2019), arXiv:1811.02042 [astro-ph.IM] .
  * Note [3] These are expected to be the most important environmental effects for LVK compact binaries [30, 50, 86].
  * Note [4] As in the main text, we use units with $G=c=1$.
  * Note [5] In these simulations the BHL radius was associated with the center of mass velocity, since the binaries considered had orbital velocities smaller than (or comparable to) the center of mass or the sound velocity.
  * Note [6] The updated sensitivity curve for ET is available at https://apps.et-gw.eu/tds/?content=3&r=18213, we used the file ET20kmcolumns.txt: HF+LF.

##  I Supplemental Material

###  I.1 Dephasing by dynamical friction (DF) and accretion

Consider a quasi-circular BBH evolving in some non-vacuum environment (e.g.,
an AGN accretion disk) where it may have formed. We focus on a stage when the
binary is sufficiently close to merger that the GW emission drives the
inspiralling, while the environment is responsible only for a subdominant
dephasing in the waveform. This regime is consistent with our assumption of
quasi-circular orbits, since, even though environmental effects may lead to
increase BBH eccentricity [64, 23], the GW radiation reaction circularizes
them [65]. In the inspiral of a LVK binary of total mass $M$, the separation
distance evolves over a length-scale $\sim
10^{3}\textrm{\,km}\,(M/60M_{\odot})^{1/3}$, which is smaller than the typical
length-scale of all environments we are interested in this work (e.g., the
size of AGN accretion disks is $\sim 10^{11}\textrm{\,km}$). Thus, we assume
that the BBH is not able to probe the medium density profile, and effectively
evolves on a _static uniform_ density environment. By assuming a static
density profile, we are neglecting the feedback of the BBH on the environment
(see, e.g., [72]); this is admittedly a strong assumption, since the binary is
expected to dynamically heat and deplete its surroundings, and may even create
a “cavity” devoid of matter (as it happens in circumbinary disks [22]).
Nevertheless, if BBHs are formed close enough to merger and do not describe
too many orbits, the (average) density of the local environment may remain
considerably unchanged (as seen in Ref. [50]).

We now follow the same approach as Ref. [30] to compute the leading PN
contribution of accretion and DF to the dephasing of the GW signal with
respect to the vacuum waveform 333These are expected to be the most important
environmental effects for LVK compact binaries [30, 50, 86].. At leading
order, the rate of energy loss into GWs is given by the quadrupole formula
444As in the main text, we use units with $G=c=1$.

$\dot{E}_{\rm{GW}}\approx\frac{32}{5}\eta^{2}v^{10},$ (3)

using also Kepler’s third law. For close to equal-mass binaries, the rate of
energy loss due to the (tangential) DF is given by Eq. 1 in the main text.
Now, at leading (Newtonian) order, the total orbital energy is

$E_{\mathrm{orb}}=-\frac{1}{2}\eta Mv^{2},$ (4)

so, from conservation of energy
$\dot{E}_{\mathrm{orb}}=-\dot{E}_{\rm{GW}}-\dot{E}_{\rm{DF}}$, we arrive at
the relation

$\dot{f}\approx\frac{96\pi^{\frac{8}{3}}}{5}\eta
M^{\frac{5}{3}}f^{\frac{11}{3}}+12\Big{(}\frac{1-3\eta}{\eta^{2}}\Big{)}\rho\mathcal{I}+\frac{5}{2}\frac{f\dot{M}}{M},$
(5)

where we used that $\dot{\eta}\approx 0$ for close to equal-mass binaries. For
$\dot{M}$ we consider either BHL or collisionless accretion. In the case of
collisional media, LVK binaries are expected to have supersonic orbital
velocities larger than their center of mass velocity with respect to the
medium, so that their individual BHL radii are typically comparable to the
binary separation distance. Led by numerical simulations (e.g., [76, 77], but
note 555In these simulations the BHL radius was associated with the center of
mass velocity, since the binaries considered had orbital velocities smaller
than (or comparable to) the center of mass or the sound velocity.), we assume
that each binary component evolves through BHL accretion, which results in the
total mass growth rate

$\dot{M}_{\rm{BHLA}}\approx\frac{4\pi\rho
M^{2}\lambda}{v^{3}}\Big{(}\frac{1-5\eta(1-\eta)}{\eta^{3}}\Big{)},$ (6)

with the fudge parameter $\lambda\sim O(1)$. In media constituted by particles
with larger mean free path, accretion is less efficient and it is better
described by collisionless accretion [82, 83].

$\dot{M}_{\rm{CA}}\approx\frac{16\pi\rho
M^{2}}{v}\Big{(}\frac{1-3\eta}{\eta}\Big{)}.$ (7)

Figure 5: Mismatch ($\mathpzc{MM}$) as function of the dephasing parameter
($\delta\Phi_{k}$) for a set a non-spinning equal-mass binary systems. The
solid line refers to the effect of CA ($k=-9$), and the dashed line to BHLA or
DF ($k=-11$).

Using the stationary phase approximation, we can write the Fourier transform
of the plus and cross polarization waveforms as

$\tilde{h}_{+,\times}=\mathcal{A}_{+,\times}(f)e^{i\phi_{+,\times}(f)},$ (8)

where at leading order (and for the mode $m=2$),

$\displaystyle\mathcal{A}_{+,\times}\approx\frac{Q_{+,\times}}{D}\eta^{\frac{1}{2}}M^{\frac{5}{6}}f^{-\frac{7}{6}},$
(9) $\displaystyle\phi_{+}\approx 2\pi ft(f)-\varphi(f)-\frac{\pi}{4},$ (10)
$\displaystyle\phi_{\times}=\phi_{+}+\frac{\pi}{2},$ (11)

with $Q_{+,\times}(\iota)$ a real function of the binary inclination $\iota$,
and where $D$ is the binary distance, and $f$ is the GW frequency. The $t(f)$
and $\varphi(f)$ are roughly the instant and the phase when the signal has a
frequency $f$; more precisely, they are defined as

$\displaystyle t(f)$ $\displaystyle\coloneqq
t_{\mathrm{c}}-\int^{+\infty}_{f}\mathrm{d}f^{\prime}\frac{1}{\dot{f}^{\prime}},$
(12) $\displaystyle\varphi(f)$
$\displaystyle\coloneqq\varphi_{\mathrm{c}}-\int^{+\infty}_{f}\mathrm{d}f^{\prime}\,\frac{2\pi
f^{\prime}}{\dot{f}^{\prime}}.$ (13)

Figure 6: Mismatch as function of the medium density for a GW150914-like event
using the actual GW150914’s PSD. The vertical line shows (approximately) the
$90\%$ upper bound from Fig. 1 for DF and BHLA. The $h_{\rm env}^{L}$ and
$h_{\rm env}^{\rm full}$ denote, respectively, the environment waveform
truncated at linear order in $\rho M^{2}$ and the full one. We consider
simultaneously the environmental effects from DF and BHLA.

In our regime of interest, in which the environmental effects are much smaller
than the gravitational radiation-reaction, DF and accretion are responsible
only for a slow dephasing of the waveform with respect to vacuum, which, at
linear order in $\rho M^{2}\ll 1$, can be expressed as

$\displaystyle\phi_{+}$ $\displaystyle\approx 2\pi
ft_{\mathrm{c}}-\varphi_{c}-\frac{\pi}{4}$
$\displaystyle\qquad+\frac{3}{128\eta\,v^{5}}\left(1+\delta_{\mathrm{vac}}^{\mathrm{PN}}+\delta_{\mathrm{DF}}+\delta_{\mathrm{accr}}\right),$
(14)

where $\delta_{\mathrm{vac}}^{\mathrm{PN}}$ contains the PN vacuum GR
corrections, and with the environmental dephasing parameters

$\delta_{\mathrm{DF}}\approx-\frac{25\pi(1-3\eta)}{304\,\eta^{3}}\mathcal{I}\rho
M^{2}v^{-11},$ (15)

and

$\delta_{\mathrm{accr}}\approx\begin{cases}-\frac{125\pi[1-5\eta(1-\eta)]}{1824\,\eta^{4}}\lambda\rho
M^{2}v^{-11}\quad\textrm{(BHLA)},\\\
-\frac{125\pi(1-3\eta)}{357\,\eta^{2}}\rho
M^{2}v^{-9}\quad\textrm{(CA)},\end{cases}$ (16)

Thus, DF and BHLA are both responsible for a $-5.5\textrm{PN}$ correction,
while CA for a $-4.5\mathrm{PN}$ one to the GW phase.

We quantify the difference between a vacuum waveform ($h_{\rm vac}$) and an
environment one ($h_{\rm env}$) by computing the _mismatch_ , which is defined
as

$\mathpzc{MM}\left(h_{\rm vac},h_{\rm env}\right)=1-\max_{t_{\rm
ref},\varphi_{\rm ref}}\tfrac{\langle h_{\rm vac}\,|\,h_{\rm
env}\rangle}{\sqrt{\langle h_{\rm vac}\,|\,h_{\rm vac}\rangle\langle h_{\rm
env}\,|\,h_{\rm env}\rangle}},$ (17)

with the maximization taken over an overall phase $\varphi_{\rm ref}$ and time
$t_{\rm ref}$. The bracket denotes the inner-product weighted by the detector
noise power spectral density (PSD). We considered the aLIGO at design
sensitivity aLIGOZeroDetHighPower [102],

ET 666The updated sensitivity curve for ET is available at https://apps.et-
gw.eu/tds/?content=3&r=18213, we used the file ET20kmcolumns.txt: HF+LF. , and
B-DECIGO [104] (an initial version of DECIGO [105]), with the corresponding
lower cutoff frequencies set to $15$ Hz, $5$ Hz, and $0.1$ Hz, respectively,
while computing the inner product. Figure 5 shows the mismatch curves for a
set of binary systems using the aLIGO sensitivity curve, indicating that with
a BNS system, the aLIGO is sensitive to much smaller dephasing coefficients
($\sim$ 4–6 orders of magnitude smaller) than in the case of BBHs, which is
due to the significantly fewer inspiral cycles in-band for the latter system.
We note that the mismatch curves are symmetric about $\delta\Phi_{k}=0$, but
the waveforms are not: a negative value of $\delta\Phi_{k}$ results in a
shorter waveform than in the vacuum model, while the opposite happens for a
positive value.

The higher order terms in $\rho M^{2}$ contributing to the _full_ environment
waveform are responsible for more negative PN order deformations than the
linear terms that we considered in our analysis, but we have checked that, for
the binary systems studied in this work, they will not change significantly
our results. More precisely, we verified that the mismatch between the
environment waveform truncated at linear order in $\rho M^{2}$ and the full
environment waveform, $\mathpzc{MM}(h_{\rm env}^{L},h_{\rm env}^{\rm full})$,
is always much smaller than the mismatch between the environment waveform
truncated at linear order and the vacuum waveform, $\mathpzc{MM}(h_{\rm
env}^{L},h_{\rm vac})$, for densities smaller than (or comparable to) the
$90\%$ upper bounds derived in Fig. 1. Figure 6 shows one such check for a
GW150914-like event, considering the simultaneous effect of DF and BHLA. The
full waveform deformation (including all powers in $\rho M^{2}$) was expressed
analytically in terms of hypergeometric functions from the evaluation of
integrals (12) and (13).

Figure 7: Marginalized posterior distributions for the dephasing parameters
$\delta\Phi_{-9}$ (in orange) and $\delta\Phi_{-11}$ (in blue) from the GW
observations of compact binaries from GWTC-1. The horizontal lines show the
90% credible intervals.

###  I.2 Posterior distributions

We have analyzed each event in GWTC-1 using the on-source data, PSD, and
calibration uncertainty provided by the Gravitational Wave Open Science Center
(GWOSC) [95, 97, 96]. We have performed the parameter estimation analyses
using the IMRPhenomPv2 and IMRPhenomPv2_NRTidalv2 waveform models for BBH and
BNS systems, respectively. We estimate the magnitude of the environmental
effects by computing the marginalized posterior probability distribution of
the $\delta\Phi_{k}$ parameter within ${\mathcal{H}}_{\rm env}$, integrating
the posterior distribution $p(\vec{\theta}\,|\,d,{\mathcal{H}}_{\rm env})$
over the nuisance parameters,

$p\left(\delta\Phi_{k}\,|\,d\right)=\int\left(\prod_{\vec{\theta}\,\setminus\\{\delta\Phi_{k}\\}}d\theta_{i}\right)\>p(\vec{\theta}\,|\,d,{\mathcal{H}}_{\rm
env}),$ (18)

where $\vec{\theta}\coloneqq\\{\theta_{i}\\}$ is the set of model parameters
for ${\mathcal{H}}_{\rm env}$. Our results for the marginalized posterior
distributions for the dephasing parameters are shown in Fig. 7.

Out of the eleven events analyzed, informative posterior distributions were
successfully obtained for nine events, for both $-4.5$PN and $-5.5$PN
deformations of the (inspiral) waveform. However, for GW170729 and GW170823,
the obtained posterior distributions were uninformative. An uninformative
posterior distribution suggests that the data did not impose stringent
constraints on the parameters. As a result, the posterior distributions
closely resembled the prior distributions. This indicates that the GW data
available for these events did not provide sufficient information to make
precise inferences about the dephasing parameters. Notably, these events
involve binaries with component BHs that stood out as the most massive in the
GWTC-1 catalog. Specifically, GW170729 has progenitor masses of
$50.2M_{\odot}$ and $24M_{\odot}$, resulting in the formation of the most
massive remnant BH in GWTC-1. Similarly, for GW170823, the component masses
were $39.5M_{\odot}$ and $29.0M_{\odot}$.

In contrast, the BNS system—with the smallest total mass in the
catalog—yielded the narrowest distributions for the dephasing parameters: an
impressive $10^{-13}$ and $10^{-15}$ for $k=-9$ and $k=-11$, respectively. The
next highest constraints pertains to the lowest-mass BBH systems: GW151226 and
GW170608. These events produced posterior distributions for the dephasing
parameters of exceptional precision, with widths of $10^{-8}$ and $10^{-9}$
for $k=-9$ and $k=-11$, respectively. Interestingly, GW170608 shows a faint
bimodal tendency. There may be different reasons for this: one factor could
stem from the presence of noise below $30$Hz in the Hanford data, which led us
to choose a lower cutoff frequency of $30$Hz for H1 and $20$Hz for L1 (the
same choice that was made by the LVK).

Despite the large mass of its BBH components, the first observed event
(GW150914) shows an unambiguous posterior distribution for both dephasing
parameters. This is attributed to the event’s high SNR, which compensates for
its brief inspiral stage (in band) due to its large mass. The extent of
information that can be derived from each segment of the signal relies on the
SNR specific to the stage of coalescence being analysed. In analyses conducted
by the LVK collaboration, parameterized inspiral tests are typically employed
solely on signals that exhibit an inspiral SNR exceeding 6 [11, 12, 13].
However, in this work, we refrain from implementing this SNR threshold. If we
were to apply this criterion, it would result in a subset of events comprising
GW150914, GW151226, GW170104, GW170608, GW170814 and GW170817. This indicates
that the remaining events have a lower inspiral SNR, posing challenges in
extracting information from their signals.

The GW151012’s posterior distribution shows an unusual pattern characterized
by a bimodal tendency within a range excluding zero. This event was observed
with a low SNR and one of the least significant events in GWTC-1 [100, 1]. The
combination of these factors, together with the presence of relatively massive
binary components ($23.2M_{\odot}$ and $13.6M_{\odot}$), contributes to the
features observed in its posterior distribution. The same reasoning applies to
GW170818 and GW170809, which also display a particular behavior, with the zero
value positioned at the edge of the $90\%$ confidence interval. These events
have binary components with masses in the range $20\text{--}30M_{\odot}$,
resulting in fewer observable cycles within the frequency range in which the
detector operates. Even though GW170104’s SNR is greater than 6, its posterior
distributions for the dephasing parameters display distinct characteristics,
being notably broad. This large component masses (and associated few cycles in
band) likely contributes to the flatness observed in these distributions.