$\displaystyle+\sum_{j=1}^{I}\sum_{m=1}^{M_{j}}\sum_{g=1}^{G}\bigg{[}B^{j}_{g}-\sum^{S_{j}}_{s=1}x^{jm*}_{s}+\delta^{jm*}_{g}\bigg{]}$
$\displaystyle\times(\mu^{1}_{jmg}-\mu^{1*}_{jmg})\geq 0,$
$\displaystyle\forall(Q^{0},Q^{1},Q^{2},\mathfrak{\Delta}^{0},\mathfrak{\Delta}^{1},d,\lambda^{0},\lambda^{1},\lambda^{2},\mu^{0},\mu^{1})\in\mathcal{K}.$

Clearly, (A) is the sum of the optimality condition (3.2), (3.3), (3.4), and
(30), and thereby, is, according to Definition 1, a cross-sector multi-product
scarce resource SCNE pattern. $\square$


# Impact of Bayesian Priors on the Inferred Masses of Quasi-Circular
Intermediate-Mass Black Hole Binaries

Koustav Chandra<EMAIL_ADDRESS>Department of Physics, Indian
Institute of Technology Bombay, Powai, Mumbai 400 076, India  Archana Pai
<EMAIL_ADDRESS>Department of Physics, Indian Institute of Technology
Bombay, Powai, Mumbai 400 076, India  Samson H. W. Leong
<EMAIL_ADDRESS>Department of Physics, The Chinese University of
Hong Kong, Shatin, N.T., Hong Kong  Juan Calderón Bustillo
<EMAIL_ADDRESS>Instituto Galego de Física de Altas Enerxías,
Universidade de Santiago de Compostela, 15782 Santiago de Compostela, Galicia,
Spain  Department of Physics, The Chinese University of Hong Kong, Shatin,
N.T., Hong Kong

###### Abstract

Observation of gravitational waves from inspiralling binary black holes has
offered a unique opportunity to study the physical parameters of the component
black holes. To infer these parameters, Bayesian methods are employed in
conjunction with general relativistic waveform models that describe the
source’s inspiral, merger, and ringdown. The results depend not only on the
accuracy of the waveform models but also on the underlying fiducial prior
distribution used for the analysis. In particular, when the pre-merger phase
of the signal is barely observable within the detectors’ bandwidth, as is
currently the case with intermediate-mass black hole binary signals in ground-
based gravitational wave detectors, different prior assumptions can lead to
different interpretations. In this study, we utilise the gravitational-wave
inference library, Parallel Bilby, to evaluate the impact of mass prior
choices on the parameter estimation of intermediate-mass black hole binary
signals. While previous studies focused primarily on analysing event data, we
offer a broader, more controlled study by using simulations. Our findings
suggest that the posteriors in total mass, mass ratio and luminosity distance
are contingent on the assumed mass prior distribution used during the
inference process. This is especially true when the signal lacks sufficient
pre-merger information and/or has inadequate power in the higher-order
radiation multipoles. In conclusion, our study underscores the importance of
thoroughly investigating similarly heavy events in current detector
sensitivity using a diverse choice of priors. Absent such an approach,
adopting a flat prior on the binary’s redshifted total mass and mass ratio
emerges as a reasonable choice, preventing biases in the detector-frame mass
posteriors.

††preprint: APS/123-QED

##  I Introduction

The field of gravitational-wave astronomy is experiencing significant
advancements owing to the frequent upgrades made to the International
Gravitational-Wave Observatory Network (IGWN) [1, 2]. Since their inception,
the Advanced LIGO and the Advanced Virgo detectors have observed $\gtrsim 100$
signals, most of which are attributed to binary black hole mergers [3, 4, 5,
6, 7, 8, 9, 10, 11]. This growing catalogue has enabled precision strong-field
tests of general relativity and has shed light on the astrophysical origin and
the nature of the black hole population in the local Universe [12, 13, 14, 15,
16, 17, 18].

Quite a few of these results rely heavily on Bayesian parameter inference. For
instance, the inspiral-merger-ringdown consistency test examines the
consistency of the signal’s low and high-frequency components by independently
inferring the mass and spin of the remnant object from each portion [19, 20].
Similarly, the residual test measures the coherent residual signal-to-noise
ratio in the data after subtracting the best-fit waveform to ensure overall
signal consistency [21]. Likewise, studying secondary parameters, such as the
remnant object’s recoil velocity, requires accurate inference of the binary’s
primary parameters, such as their masses and spins [22, 23, 24, 25, 26, 27,
28, 29, 30]. Also, to determine the shape of the black hole mass or spin
spectrum, we must first accurately obtain the astrophysical properties of the
observed binary black holes. Therefore, addressing the potential issues in
Bayesian parameter inference is essential.

One such issue is the impact of prior choices in signal parameter estimation.
In gravitational-wave astronomy, Bayesian parameter inference involves
analysing signal-containing data segments from a detector network with a prior
distribution to obtain posterior samples of binary parameters. If the analysed
data lacks sufficient information, prior assumptions will significantly
influence the resulting posteriors. For example, GW190521 has had multiple
interpretations due to its barely observable pre-merger phase [31, 32, 33, 34,
35, 36, 37, 38, 39, 40]. Even under the quasi-spherical signal hypothesis, the
binary constituents have a non-negligible probability of occupying the pair-
instability gap or straddling it, depending on the prior choice.

This paper studies the effect of competing mass prior choices on synthetic
intermediate-mass black hole binary signals. This work is motivated by similar
previous studies [39, 40, 41] which highlighted the potential limitations of
certain prior choices in exploring the high-likelihood regions of parameter
space. However, these previous studies primarily focused on analysing data
segments containing the GW190521 event (or other comparatively low-
significance events). In contrast, our study offers a broader and more
controlled perspective as we use simulations.

Our study shows that the posteriors in binary’s redshifted total mass, mass
ratio, and luminosity distance critically depend on the mass prior choice.
This effect is particularly accentuated when the signal has a barely
observable inspiral phase and/or has insufficient power in the higher-order
radiation multipoles. The amplitudes and phases of these higher-order
radiation multipoles relative to the dominant quadrupole mode are determined
by the binary properties, such as mass ratio and inclination and therefore
encode them. We, therefore, advocate mapping the likelihood of similarly high-
mass events with diverse mass prior choices whenever possible, although
reweighting may be enough.

The remainder of this paper is structured as follows. Section II briefly
describes Bayesian inference and summarises the notations and conventions of
the binary parameters used in this article. The following section, Section
III, summarises our simulation setup and our prior choices. The results of our
analysis are in Section IV, and a summary and conclusion of our study are in
Section V.

##  II Preliminaries

###  II.1 Bayesian inference in a nutshell

Bayesian parameter inference libraries, like Bilby, model the calibrated
detector output $d(t)$ as

$d(t)=n(t)+s(t\mid\bm{\theta})$ (1)

where $n(t)$ is the noise, and $s(t|\bm{\theta})$ is the signal strain [42,
43, 44, 45, 46]. These libraries compute the posterior probability
distribution of the signal parameters $\bm{\theta}$ by using the Bayes
theorem:

$p(\bm{\theta}\mid
d,\mathcal{M})=\frac{\mathcal{L}(d\mid\bm{\theta},\mathcal{M})\pi(\bm{\theta}\mid\mathcal{M})}{\mathcal{Z}_{\mathcal{M}}}~{}.$
(2)

Here, $\mathcal{L}(d|\bm{\theta},\mathcal{M})$ is the likelihood function, and
$\pi(\bm{\theta}|\mathcal{M})$ is the prior probability distribution which is
typically chosen to avoid imprinting astrophysical assumptions on the results.
$\mathcal{Z}_{\mathcal{M}}=\int
d\bm{\theta}\mathcal{L}(d|\bm{\theta},\mathcal{M})\pi(\bm{\theta}|\mathcal{M})$
is the probability of observing the data given the model $\mathcal{M}$,
otherwise known as the model evidence. The likelihood function is defined as
follows:

$\mathcal{L}\left(d\mid\bm{\theta},S_{n}\right)\propto\exp\left[-\sum_{i}\frac{2\left|\tilde{d}\left(f_{i}\right)-\tilde{h}\left(f_{i}\mid\bm{\theta}\right)\right|^{2}}{TS_{n}\left(f_{i}\right)}\right]~{}.$
(3)

assuming wide-sense stationary Gaussian noise and involves the data segment
duration $T$, the noise power spectrum $S_{n}(f_{i})$ at each frequency
$f_{i}$, the Fourier representations of the data $\tilde{d}(f_{i})$ and the
template $\tilde{h}(f_{i}\mid\bm{\theta})$. When using data from multiple
detectors, Bilby obtains the joint likelihood by multiplying individual
likelihoods for each detector $j$ while requiring the template to be coherent
across the $N_{\mathrm{ifo}}$ detector network:

$\mathcal{L}(\\{d\\}\mid\bm{\theta},\\{S_{n}\\})=\prod_{j}^{N_{\mathrm{ifo}}}\mathcal{L}\left(d_{j}\mid\bm{\theta},S_{n}\right)~{}.$
(4)

To compute this likelihood function and estimate the posterior probability
distribution described in Eq. (2), Bayesian inference algorithms rely on off-
the-shelf stochastic samplers. Given the dimensionality involved, this can be
computationally expensive for some waveform models and prior choices.
Therefore, it is preferable to use marginalised likelihood over specific
parameters whenever possible. The marginalised parameters can be retrieved
later using post-processing techniques.

Posterior sampling generates a multi-dimensional probability density function
that describes our joint inference about the signal parameters. Understanding
and/or visualising the full posterior probability distribution is tricky. So,
we often focus on just one or two parameters and marginalise the others. For
example, if we are interested in the joint posterior for component masses
$m_{1},~{}m_{2}$ of an event, we integrate over the other parameters as
follows:

$p\left(m_{1},m_{2}\mid\\{d\\}\right)=\sum_{\bm{\theta}_{\left(-m_{1},-m_{2}\right)}}p\left(\bm{\theta}\mid\\{d\\}\right)$
(5)

This marginalised posterior indicates what we can infer about the component
masses when the uncertainties in the omitted parameters are considered. Once
obtained, the marginal posterior samples can be visualised using histograms or
kernel density plots. We can also calculate its summary statistics, such as
the median and estimate the uncertainty in the parameter(s) of interest by
constructing symmetric credible intervals (CI) defined as follows:

$\sum_{\theta_{k}=-\infty}^{a}p\left(\theta_{k}\mid\\{d\\}\right)=\frac{1-X}{2}=\sum_{\theta_{k}=b}^{\infty}p\left(\theta_{k}\mid\\{d\\}\right)~{}.$
(6)

Here, the interval $(a,b)$ is a symmetric 100X% CI for parameter $\theta_{k}$.
In this article, we choose $X=0.68,0.9$ and construct two-dimensional marginal
posteriors to consider the potential influence of one parameter on the other.
The 68% CI corresponds to approximately one standard deviation from the mean
of a normal distribution. Conversely, the 90% CI provides a wider interval,
encompassing a higher degree of uncertainty and has traditionally been
employed in gravitational-wave data analysis.

Finally, when multiple models are available, it is natural to quantify how
much one model, say $A$, is preferred to another, say $B$. Within the Bayesian
framework, the Bayes factor is the primary tool for comparing the preferences
of two competing models. The log Bayes’ factor for $A$ over model $B$ is
calculated by subtracting the log model evidence:
$\ln\mathrm{BF}^{A}_{B}=\ln\mathcal{Z}_{A}-\ln\mathcal{Z}_{B}$. A
$\ln\mathrm{BF}^{A}_{B}>1$ indicates positive support for model $A$ over model
$B$ [47].

###  II.2 Quasi-circular black hole binary parameters: Notations and
Conventions

For quasi-circular or non-precessing black hole binaries 111Precessing
binaries are quasi-spherical as the orbital orientation evolves with time.,
$\bm{\theta}$ represents four intrinsic and seven extrinsic parameters. The
intrinsic parameters include the component’s masses $m_{i}$ and the
dimensionless spins $\chi_{i}$ aligned with each other and the orbital angular
momentum. The extrinsic parameters comprise the source’s two-dimensional sky
location $(\alpha,\delta)$, the luminosity distance $D_{L}$ (or equivalently
the redshift $z$ to the source), the source’s inclination $\iota$, the
polarisation angle $\psi$, the azimuth $\varphi$ and the merger time $t_{c}$.

It is possible to capture the dominant effect of spin on the inspiral rate
using the mass-weighted effective aligned spin parameter
$\chi_{\mathrm{eff}}$. At a given reference frequency, say $f_{\mathrm{ref}}$,
during the binary’s inspiral, it is defined as follows [48, 49, 50]:

$\chi_{\mathrm{eff}}=\frac{\chi_{1}+Q\chi_{2}}{1+Q}~{},$ (7)

using the binary’s spins and the mass ratio $Q=m_{2}/m_{1}\leq 1$. Aligned
spin binaries are characterised by positive values of $\chi_{\mathrm{eff}}$
whereas anti-aligned spin binaries have negative $\chi_{\mathrm{eff}}$. Low-
spin binaries have $\chi_{\mathrm{eff}}\rightarrow 0$.

The redshifted (detector-frame) total mass $M_{T}(1+z)$ is crucial for
detection, as it determines the signal’s overall strength and frequency
content. We express it and the individual component masses in solar mass units
in agreement with other works.

##  III Methodology

###  III.1 Data

Parameter | Symbol | Value
---|---|---
Detector-frame total mass | $M_{T}(1+z)$ | $268.83M_{\odot}$
Azimuth | $\varphi$ | 0.004
Polarization angle | $\psi$ | 2.38
Coalescence GPS time | $t_{c}$ | $1242442967.41\mathrm{~{}s}$
Right ascension | $\alpha$ | 0.16
Declination | $\delta$ | -1.14
Table 1: True parameter values that are shared by both the simulation sets.
They represent the maximum likelihood point of the GW190521 posterior samples
that LIGO-Virgo has publicly released [51]. All angles are in radians.

We aim to explore the differences between posteriors sampled under different
mass priors for systems with detector frame masses similar to GW190521. We,
therefore, create two sets of simulated quasi-circular intermediate-mass black
hole binary signals for our analysis. The first set uses the phenomenological
waveform model IMRPhenomXAS, which models just the dominant waveform harmonic,
$(2,2)$ of a quasi-circular black hole binary [52]. The second set of signals
uses the phenomenological waveform model IMRPhenomXHM, which includes the
modes $(2,1),(3,3),(3,2),(4,4)$ in addition to the $(2,2)$ mode [53].

Parameter | Symbol | Value
---|---|---
Mass ratio | $Q$ | 0.82, 0.25
Inclination | $\iota$ |  $28.6^{\circ}$, $54.4^{\circ}$
Primary Spin | $\chi_{1}$ | 0.054, 0.7, -0.7
Secondary Spin | $\chi_{2}$ | 0.144, 0.4
Table 2: True parameter values for the twelve simulated binaries from the
first simulation set.

These simulated signals share a set of common parameter values, as shown in
Table 1. However, they differ in their $Q$, $\iota$, and $\chi_{\mathrm{eff}}$
values. Even though we perform our study with several combinations of these
varying parameters, we only report the results for simulations with specific
combinations to summarise the features. These combinations are listed in Table
2 and Table 3 for simulation set-1 (quadrupole-only signal) and simulation
set-2 (multipole signal), respectively.

Parameter | Symbol | Value
---|---|---
Mass ratio | $Q$ | 0.82, 0.25, 0.11
Inclination | $\iota$ |  $28.6^{\circ}$, $54.4^{\circ}$, $74.5^{\circ}$
Primary Spin | $\chi_{1}$ | 0.054
Secondary Spin | $\chi_{2}$ | 0.144
Table 3: True parameter values for the nine simulated binaries from the second
simulation set.

We inject these simulated signals to the Gaussian noise of a three-detector
network comprising the Advanced LIGO and Advanced Virgo interferometer
operating at the expected fourth observing run sensitivity [54] 222We present
the results with the analysis performed in Gaussian noise. The qualitative
nature of our analysis does not change even when the analysis is performed in
zero-noise.. The distance to each system is adjusted to fix the network
optimal SNR at 20. Though we show results for a network SNR of 20, our
analysis has produced consistent results for optimal network SNR values of 15
and 30.

Increasing the binary’s mass asymmetry and inclination decreases the system’s
intrinsic luminosity, thus placing it at a relatively closer distance than
mass-symmetric systems. However, the former systems will have significantly
larger contributions from the non-quadrupole modes, particularly for
simulation set 2. This can help reduce the correlation between certain signal
parameters, such as source inclination and distance [55, 56].

Given that the primary objective of our work is to provide a preliminary
understanding of the measurability of mass parameters based on prior choices,
we have intentionally excluded quasi-spherical systems from our analysis.
Relativistic precession introduces amplitude and phase modulations in the
signal inspiral, which is limited for the heavy systems considered in this
study. Consequently, we anticipate that the fundamental qualitative attributes
of our findings would remain unaltered even if we use quasi-spherical signals.
Nonetheless, it is essential to underscore that relativistic precession,
particularly for certain specific spin configurations, can significantly
suppress the inspiral phase before the merger, affecting measurements [57,
58].

###  III.2 Inference Settings

To obtain the posterior samples, we use the parallelised version of the
parameter estimation package Bilby [42, 43] with the nested sampler dynesty
[59]. In our analysis, we fix the number of walks and maxmcmc at 200 and
15000, respectively. Additionally, we set the number of live points, denoted
as nlive, to 1024, and the number of autocorrelation times, denoted as nact,
to 30. Other dynesty attributes, such as naccept, proposals, etc., are set to
their default values.

To expedite convergence, we downsample our simulated data to 1024 Hz and
utilise marginalised likelihood with marginalisation over distance and time,
which is performed numerically following the approach of Farr [60] and Singer
and Price [61]. Further, while using the IMRPhenomXAS waveform model, we use
the likelihood function that is analytically marginalised over the signal
phase [62]. We compute all spin posteriors at $f_{\mathrm{ref}}=11$Hz and use
an identical lower-frequency cutoff of 11Hz for likelihood integration. The
high-frequency cutoff is set to the Nyquist frequency, 512Hz, corresponding to
the data’s sampling rate of 1024Hz.

###  III.3 Priors

Figure 1: The figure compares mass ratios $Q$ resulting from four different
prior choices used in our analysis. Unlike the other three priors, Prior-4
distinctly favours systems with unequal component masses. Also, as evident,
Prior-2 effectively imposes a flat prior on $Q$.

Except for the mass priors, we adopt default prior choices for our quasi-
circular binary black hole analysis, with some adjustments based on our
injections [63, 45, 7]. Table 4 summarises these settings for the signal
parameters. The maximum spin values are set to 0.99 to ensure consistency with
the employed waveform models. We employ the UniformSourceFrame luminosity
distance prior, which implements a uniform prior in comoving volume $V_{c}$
with an initial factor $1/1+z$, assuming the Planck15 cosmology — the factor
of $1/1+z$ accounts for time dilation due to the expansion of the Universe
[45].

Parameter | Shape | Limits | Boundary
---|---|---|---
$\chi_{1},\chi_{2}$ | Uniform | (0, 0.99) | -
$\iota$ | Sinusoidal | $(0,~{}\pi)$ | -
$\psi$ | Uniform | $(0,~{}\pi)$ | Periodic
$\varphi$ | Uniform | $(0,~{}2\pi)$ | Periodic
$\alpha$ | Uniform | $(0,~{}2\pi)$ | Periodic
$\delta$ | Cosinusoidal | $(-\pi/2,~{}\pi/2)$ | -
$D_{L}$ | UniformSourceFrame | (0.1, 10)Gpc | -
Table 4: Default prior settings of signal parameters to study binary black
hole signals.

As for the mass priors, we explore the following four choices for our
analysis. _Prior-1_ adopts a flat prior on the detector frame chirp mass
$M_{c}(1+z)$ in the range of $[30,300]M_{\odot}$ and $Q$ in the range of
$[0.1,1]$ while constraining the individual masses $m_{i}(1+z)$ to lie within
$[5,300]M_{\odot}$. _Prior-2_ aligns with the prior used in LIGO/Virgo’s
GW190521 analysis, setting flat priors on the component masses $m_{i}(1+z)$ in
the range of $[5,300]M_{\odot}$ while ensuring that the mass ratio and the
detector frame total mass range lie within $Q\in[0.1,1]$ and
$M_{T}(1+z)\in[80,400]M_{\odot}$ 333Please note that Abbott  _et al._ [31]
implemented a prior which is $\pi(d_{L})\propto d_{L}^{2}$ which distribute
mergers uniformly throughout Euclidean Universe. Throughout our study, we
adopt the UniformSourceFrame as our distance prior.. _Prior-3_ imposes a flat
prior on $Q\in[0.1,1]$ and $M_{T}(1+z)\in[80,400]M_{\odot}$ while constraining
the individual masses $m_{i}(1+z)$ to lie within $[5,300]M_{\odot}$. Finally,
_Prior-4_ , following Nitz and Capano [33], is flat in $1/Q\in[1,10]$ and
$M_{T}(1+z)\in[80,400]M_{\odot}$. Figure 1 compares the mass ratio
distribution of our prior choices. Prior-4 favours unequal mass ratio systems,
while the first three have no such preference. As meticulously elucidated in
Estellés  _et al._ [39], Olsen  _et al._ [40], Bustillo  _et al._ [41], when
examining the GW190521 event under the quasi-spherical signal hypothesis, such
prior choices can notably impact the signal inference.

Irrespective of our prior choice, we sample in $M_{c}(1+z)$ and $Q$ space
following Bilby’s default setting. Bilby accomplishes this by utilising
appropriate transformations between the chosen parameterisation and the
$M_{c}(1+z)-Q$ representation, thereby folding the chosen prior into
$M_{c}(1+z)-Q$ space.

Finally, to assess the similarity/dissimilarity between posterior probability
distributions obtained using different prior choices, we measure the Jensen-
Shannon divergence (JSD) defined as follows [64]:

$\displaystyle\mathrm{JSD}(p_{1},p_{2})$
$\displaystyle=\frac{1}{2}\sum_{k}p_{1}(\theta_{k}\mid
d)\ln\left(\frac{2p_{1}(\theta_{k}\mid d)}{p_{1}(\theta_{k}\mid
d)+p_{2}(\theta_{k}\mid d)}\right)$
$\displaystyle\quad+\frac{1}{2}\sum_{k}p_{2}(\theta_{k}\mid
d)\ln\left(\frac{2p_{2}(\theta_{k}\mid d)}{p_{1}(\theta_{k}\mid
d)+p_{2}(\theta_{k}\mid d)}\right)~{}.$

Here, $p_{1}(\theta_{i}\mid d)$ and $p_{2}(\theta_{i}\mid d)$ represents the
posterior of parameter $\theta_{i}$ for two different prior choices. A value
of 0 nats 444“nat” is a unit of information based on natural logarithms [65].
indicates that no additional information is gained when transitioning from one
posterior distribution to another, implying that the two distributions are
similar. On the other hand, a value of $\ln(2)=0.69$ nat represents the
maximum divergence between the distributions.

##  IV Results

###  IV.1 Simulation Set-I

(a) For simulated systems at an inclination $\iota=28.6^{\circ}$

(b) For simulated systems at an inclination $\iota=54.4^{\circ}$

Figure 2: Comparison of joint posterior distributions for $M_{T}(1+z)$ and $Q$
for different simulated systems with $Q=0.82$ and varying
$\chi_{\mathrm{eff}}$ as summarised in Table 2. The contours represent the
credible intervals (CI), with the smaller contour indicating the 68% CI and
the larger contour representing the 90% CI. The different colours in the
figure represent the various mass prior choices used in the analysis, while
the dotted line corresponds to the true mass parameters of the injected
signals.

(a) For simulated systems at an inclination $\iota=28.6^{\circ}$

(b) For simulated systems at an inclination $\iota=54.4^{\circ}$

Figure 3: Comparison of joint posterior distributions for $M_{T}(1+z)$ and $Q$
for different simulated systems with $Q=0.25$ and varying
$\chi_{\mathrm{eff}}$ as summarised in Table 2.

We first present the results obtained using quasi-circular dominant waveform
harmonics. As mentioned in Sec. III.1, we simulate and recover set-1, detailed
in Table 2, using the waveform approximant IMRPhenomXAS to prevent any
systematic errors caused by a disagreement between waveform models. Therefore,
the effects discussed in this section are entirely due to different prior
assumptions made for the analysis. We focus on the system’s redshifted total
mass and mass ratio since these parameters are better constrained for the
heavy binaries discussed. We note that most of the signal SNR comes from the
merger and ringdown phases for all our systems. Consequently, during
likelihood evaluation and thus posterior sampling, the sampler prefers
waveforms that better describe the post-inspiral phase.

#### Nearly equal-mass binaries

Figure 2 shows the two-dimensional posterior distribution of $M_{T}(1+z)$ and
$Q$ for systems with injected $Q=0.82$ at varying inclinations and aligned
spin values. We find that the JSD value for $M_{T}(1+z)$ and $Q$ posteriors
lie in the range [0.0008, 0.0100] for Prior-1, 2 and 3, indicating that within
statistical uncertainties, the samples drawn using Prior-1, 2 and 3 have an
identical distribution. Further, the figure illustrates that the posterior’s
68% CI consistently encompasses the injected value across inclinations and
spin values for Prior-1 through Prior-3. However, the behaviour differs for
Prior-4 due to its preference towards unequal mass systems. Notably, at
$\chi_{\mathrm{eff}}=-0.2$ and $\chi_{\mathrm{eff}}=0.56$, the injected value
lies outside the posterior’s 68% CI for this prior choice. Also, this prior
choice leads to a comparatively broader posterior distribution of masses for
all spin values and inclinations, with a notable impact on mass ratios.

The spread of the posterior for all prior choices is least when the systems
have $\chi_{\mathrm{eff}}=0.56$ and most for $\chi_{\mathrm{eff}}=-0.2$ 555For
$\chi_{\mathrm{eff}}=-0.2$ system, $Q=0.43^{+0.44}_{-0.25}$ and for
$\chi_{\mathrm{eff}}=0.09$ system, $Q=0.52^{+0.40}_{-0.31}$ (These
measurements are reported as symmetric 90% CI around the median of the
marginalized posterior distribution).. This is due to the orbital-hangup
effect, which results in a delay (or acceleration) in the merger of aligned
(anti-aligned) spin systems relative to systems with lower spin
($\chi_{\mathrm{eff}}=0.09$) [66]. Consequently, our set of aligned spin
systems spends $\sim 0.05$ seconds more within the detector bandwidth than
anti-aligned systems, thus providing more information about the binary
components and reducing measurement uncertainties. This trend holds
irrespective of the binary’s inclination. Additionally, altering parameters
such as the number of live points (e.g., nlive=2048) or the network optimal
SNR of the system to 30 does not alter the qualitative nature of these
findings.

#### Asymmetric-mass binaries

Similar to $Q=0.82$ systems, we observe that for $Q=0.25$ systems, the JSD
values of $M_{T}(1+z)$ and $Q$ posteriors lie within the range of
$[0.0024,0.0099]$ for Prior-1, 2, and 3, indicating a coherence between their
distributions for these specific prior selections. However, the overall
results for these systems are quite different from $Q=0.82$ systems. As
illustrated in Figure 3, the accuracy of the inference for the first three
priors is unsatisfactory for anti-aligned and low-spin systems, which
contrasts with the previously observed consistent agreement. The inferred
masses’ 90% CI (indicated by lighter contours) fail to enclose the injected
parameter values adequately. This stems from the flat nature of these priors
in the $Q$ parameter and our chosen prior on luminosity distance, which
together favour intrinsically luminous systems — systems that are mass-
symmetric and located further away from us.

Conversely, Prior-4, as depicted in Figure 1, exhibits a propensity for
asymmetric mass systems, resulting in a better recovery of the injected value.
The injected parameter value lies within the 90% CI and the 68% CI for anti-
aligned and low spin systems, respectively, when using Prior-4. As for aligned
spin systems, the maximum a posteriori estimate of $M_{T}(1+z)$ and $Q$ is
close to the injected parameter values for all prior choices, irrespective of
inclination angle choice as evident from Figure 3.

Also, asymmetric systems tend to have a larger total mass posterior than mass-
symmetric binaries. This is expected as an increase in mass asymmetry
corresponds to a decrease in intrinsic loudness of the gravitational wave
signal for a fixed total mass. Consequently, to account for this diminished
signal strength, the resulting posterior distribution for the total mass
naturally extends toward larger values of $M_{T}(1+z)$.

Therefore, in summary, owing to a flat prior on $Q$, the $M_{T}(1+z)$ and $Q$
posteriors obtained using Prior-1, 2 and 3 remain identical for different mass
ratios, inclinations, and spins. However, due to their preference towards
symmetric mass systems, their capability to effectively retrieve higher mass
ratio systems with low spins or anti-aligned configurations is limited for the
quasi-circular quadrupole case.

(a) For simulated systems with $Q=0.82$

(b) For simulated systems with $Q=0.25$

(c) For simulated systems with $Q=0.11$

Figure 4: Joint posterior distributions for $M_{T}(1+z)$ and $Q$ of simulated
systems with parameters summarised in Table 3.

###  IV.2 Simulation Set-2

We now examine the results obtained when analysing the injections described in
Table 3 with IMRPhenomXHM and different prior choices. This set of simulated
binary systems has low spin and is oriented at three different inclination
angles. Given that the posterior distributions obtained using Prior-1,
Prior-2, and Prior-3 yield similar results (with a JSD value of $M_{T}(1+z)$
and $Q$ falling within the range of [0.0016, 0.0099]), we only consider the
results for Prior-2 (LVK-like prior) and Prior-4 presented in Fig 4. The area
inside the $M_{T}(1+z)-Q$ contour decreases compared to the quadrupole
analysis performed with simulation set-1, indicating a lower measurement
uncertainty . This is expected as the presence of higher harmonics in the
signal facilitates a more accurate determination of the binary parameters in
contrast to systems dominated by the quadrupole harmonic alone. Furthermore,
the correlation between the $M_{T}(1+z)$ and $Q$ for a fixed mass ratio
injection decreases as the inclination angles increase from lower to higher
values. Similarly, the correlation becomes less from higher to lower mass
ratios for a fixed inclination angle. This behaviour arises due to the
increasing strength of the higher harmonics relative to the dominant
quadrupole harmonic as the inclination and mass asymmetry rise, reducing the
correlation between these two parameters of the binary system.

However, for Prior-4, the estimate’s uncertainty decreases as we move from low
($\iota=28.6^{\circ}$) to moderate inclination ($\iota=54.4^{\circ}$) but
grows as we move from moderate to high inclination ($\iota=74.5^{\circ}$)
except for $Q=0.82$ system. Also, the posterior probability density of
$M_{T}(1+z)$ is bimodal for $Q=0.25$ and $Q=0.11$ systems. We attribute this
to the prior’s preference towards lower mass ratio systems and the sampler’s
proclivity towards binaries with a lower mass ratio to describe the post-
inspiral part of the signal better.

Lastly, it is essential to note that for the injection with a low mass ratio
and low inclination ($Q=0.11,\iota=28.6^{\circ}$), the actual value of the
parameters, as indicated by the blue lines, lies outside the 90% CI of Prior-2
but inside the 68% CI of Prior-4, indicating that the prior distribution,
rather than the higher harmonic content, is more crucial in such cases. As the
injection’s inclination increases, the actual parameter value lies inside the
68% CI for all prior choices, thus illustrating the opposite.

(a) For simulated systems with $Q=0.82$

(b) For simulated systems with $Q=0.25$

(c) For simulated systems with $Q=0.11$

Figure 5: The 2-dimensional posterior distribution for simulated systems’
redshifted final mass and final spin.

###  IV.3 Impact of mass prior choices on remnant properties

Figure 6: The figure compares the priors on the binary’s redshifted final mass
and final spin. Both the prior choices strongly prefer remnants with final
spin $\sim 0.6-0.7$.

The properties and configuration of an inspiralling black hole binary fully
specify the remnant’s mass, spin and kick. Therefore, in practice, we need
numerical relativity simulations to compute the remnant properties for
arbitrary binary mass and spin sets. However in this article, we use surrogate
fits from Healy and Lousto [67], Hofmann  _et al._ [68], Jiménez-Forteza  _et
al._ [69] to infer the final redshifted mass $M_{f}(1+z)$ and final spin
$\chi_{f}$.

We find that as the binary’s inclination increases and/or the mass ratio
decreases, the uncertainties associated with the final spin estimates
increase, resulting in broader posterior probability distributions. This is
because both priors induce a final spin prior that strongly favours binary
remnants with $\chi_{f}\sim 0.6-0.7$, as shown in Fig. 6. Consequently, the
resulting measurement uncertainty is lower when the true spin is within this
range, making the likelihood peak near the bulk of the prior, as we observe in
the case of the $Q=0.82$ binary system. In contrast, when the true final spin
is within the edge of the prior support, likelihood and prior compete, leading
to broader posterior distributions, as in the case of the $Q=0.25$ and
$Q=0.11$ binaries. The posterior support would progressively move towards the
likelihood peak if the SNR of the injection is increased, similarly to what is
described in Leong  _et al._ [70].

We also observe that the inferred posteriors peak closer to the true values
when Prior-4 is employed. Additionally, for low $Q$ binaries, we observe that
the true value falls outside the 68% CI of the posteriors generated by Prior-2
as the sampler is unable to probe the higher likelihood regions efficiently.

Figure 7: The 2-dimensional posterior distribution for the luminosity distance
(redshift) and inclination of simulated systems summarised in Table 3.

###  IV.4 Impact on inclination and redshift measurement

| Prior-2 | Prior-4
---|---|---
$\iota$ | $m_{1}~{}[M_{\odot}]$ | $m_{2}~{}[M_{\odot}]$ | $\mathbb{P}(m_{1}\in\mathrm{PI})$ | $\mathbb{P}(m_{2}\in\mathrm{PI})$ | $m_{1}~{}[M_{\odot}]$ | $m_{2}~{}[M_{\odot}]$ | $\mathbb{P}(m_{1}\in\mathrm{PI})$ | $\mathbb{P}(m_{2}\in\mathrm{PI})$
$28.6^{\circ}$ | $75.9^{+12.3}_{-9.0}$ | $63.6^{+12.2}_{-14.7}$ | 0.99 | 0.43 | $76.6^{+14.4}_{-9.5}$ | $61.7^{+12.9}_{-18.0}$ | 0.99 | 0.34
$54.4^{\circ}$ | $81.7^{+13.2}_{-11.7}$ | $69.9^{+12.4}_{-15.7}$ | 1.0 | 0.73 | $81.5^{+13.7}_{-11.7}$ | $69.0^{+12.8}_{-18.1}$ | 0.99 | 0.68
$74.5^{\circ}$ | $84.3^{+14.1}_{-12.7}$ | $71.5^{+13.0}_{-14.9}$ | 1.0 | 0.78 | $88.0^{+13.8}_{-13.2}$ | $73.0^{+13.3}_{-16.0}$ | 1.0 | 0.81
Table 5: Component masses and their probability of being within the pair-
instability (PI gap) of multipole simulations with $Q=0.82$. We report the
median values of the masses with their 90% symmetric credible intervals.

To measure the source properties of a binary black hole and study its
astrophysical implications, it is necessary to determine its source frame
masses. Conventionally, we rely on standard cosmology and the gravitational
wave measurement of the luminosity distance to estimate the cosmological
redshift. This redshift estimate is subsequently used to convert detector-
frame masses to source-frame masses. Given this, it is critical to evaluate
the impact of our chosen mass priors on the measurement of luminosity distance
or, equivalently, redshift. This section presents the results obtained using
simulation set-2 as the inclination angle $\iota$ and the distance $D_{L}$
degenerate for quadrupole signals. Further, we only discuss the results
obtained using Prior-2 and Prior-4 as the results obtained using Prior-1,
Prior-2, and Prior-3 exhibit comparable characteristics.

Figure 7 displays the two-dimensional posteriors of redshift and inclination
for each synthetic binary. We find at $Q=0.82$ and $Q=0.25$, the redshift
posterior of low and moderate inclination binaries peaks closer to the true
value of the redshift, irrespective of prior choice. However, when the binary
inclination becomes large, the true redshift value is relatively poorly
recovered for either prior choice.

A subtle disparity emerges for low inclination binaries with $Q=0.11$.
Specifically, the true value lies outside the 90% CI when using Prior-2, yet
it is within the 68% CI when Prior-4 is used. However, as the binary
inclination increases, the measurement bias decreases. Also, the inclination
angle measurement for $Q=0.25$ system is bi-modal, whereas the true
inclination angle for all other binaries lies within the 68% CI for both the
prior choices. Additionally, consistent with other studies [55, 56, 71, 72],
the correlation between distance and inclination angle weakens as the relative
strength of higher harmonics increases with increasing inclination or
decreasing mass ratio. We leave the detailed investigation to future work.

###  IV.5 Pair-instability mass gap

Stellar physics predicts the existence of a black hole mass gap with no first-
generational black holes in the mass range between $65M_{\odot}$ and
$120M_{\odot}$ owing to the development of pair-instability in massive stars
with helium core masses in this range [73, 74, 75, 76]. While the precise
location of the pair-instability mass gap remains theoretically uncertain (See
for e.g. [77, 32] and references therein), we assume that the lower boundary
is at $65M_{\odot}$ and the upper boundary is at $120M_{\odot}$. The
observation of GW190521 by the LVK [31] challenged this prediction, as the
largest merging black hole squarely fits within the mass gap. Given the strong
impact that priors can have on mass estimates, it is relevant to ask whether
prior choices could induce such observation. We test this idea here by
checking if any of our priors falsely places our component masses inside or
outside the mass gap. For instance, Calderón Bustillo  _et al._ [37] showed
that such an effect can happen when highly massive and highly eccentric
sources are mistakenly taken as precessing ones.

Table 5 shows the inferred source-frame component masses ($m_{1},~{}m_{2}$) of
systems with $Q=0.82$. For increasing inclination, these simulations have
$m_{1}$ of $\sim(75,85,105)M_{\odot}$ and $m_{2}$ of
$\sim(63,73,85)M_{\odot}$. Therefore, except for the lowest inclination, both
components lay in the pair-instability mass gap.

We find that for the lowest inclination, $\iota=28.6^{\circ}$, Prior-4 places
the secondary component outside the mass gap with a slightly higher
probability than Prior-2, making it a better choice. Similarly, Prior-4 puts
the secondary black hole within the mass gap at the highest inclination with
somewhat higher confidence, consistent with the true values. However, at a mid
inclination, Prior-2 does a moderately better job. This shows that while prior
choices can mildly affect the astrophysical interpretation of our injected
sources, none of them lead to false conclusive results. For our remaining
injections, for which the true masses straddle the mass gap, both priors place
the masses outside the gap with probabilities $\gtrsim 0.97$.

###  IV.6 Model Selection

| Prior-2 | Prior-4
---|---|---
$\ln\mathrm{BF}^{S}_{N}$ | $\ln\mathcal{L}_{\mathrm{max}}$ | $\ln\mathrm{BF}^{S}_{N}$ | $\ln\mathcal{L}_{\mathrm{max}}$
$Q=0.82~{}\iota=28.6^{\circ}$ | 188.29 | 207.77 | 187.20 | 208.17
$Q=0.82~{}\iota=54.4^{\circ}$ | 186.63 | 207.62 | 185.26 | 207.63
$Q=0.82~{}\iota=74.5^{\circ}$ | 179.08 | 200.93 | 177.15 | 201.75
$Q=0.25~{}\iota=28.6^{\circ}$ | 188.45 | 208.82 | 188.47 | 208.94
$Q=0.25~{}\iota=54.4^{\circ}$ | 187.08 | 208.95 | 187.04 | 208.88
$Q=0.25~{}\iota=74.5^{\circ}$ | 175.30 | 197.98 | 176.16 | 199.31
$Q=0.11~{}\iota=28.6^{\circ}$ | 184.76 | 206.12 | 184.85 | 207.14
$Q=0.11~{}\iota=54.4^{\circ}$ | 201.88 | 225.36 | 203.83 | 226.91
$Q=0.11~{}\iota=74.5^{\circ}$ | 196.52 | 221.15 | 198.94 | 223.15
Table 6: $\ln\mathrm{BF}^{S}_{N}$ and $\ln\mathcal{L}_{\mathrm{max}}$ values
of simulation set-2 for prior choices, 2 and 4.

Given gravitational wave data, the support for a given prior is given by the
Bayes Factor. This is defined as the prior-averaged value of the log-
likelihood $\ln\mathrm{BF}^{S}_{N}$ across the whole parameter space and,
therefore, is bounded by the maximum log-likelihood
$\ln\mathcal{L}_{\mathrm{max}}$. In principle, we expect the latter to peak
near the true injection parameters, provided that the priors do not prevent a
correct sampling of the parameter space.

Table 6 shows the log Bayes factor for the signal versus noise model for our
priors 2 and 4 within the context of the simulation set 2. For $Q=0.82$
systems and the two lowest inclinations, both priors provide similar
$\ln\mathrm{BF}^{S}_{N}$ and $\ln\mathcal{L}_{\mathrm{max}}$. For larger
inclinations and the $Q=0.82$ case, however, we recall that Prior-2 fails to
include the right luminosity distance in its $90\%$ CI while Prior-4 does.
Yet, Prior-2 yields a larger Bayes Factor while finding a worse fit (lower
$\ln\mathcal{L}_{\mathrm{max}}$) to the injection. On the one hand, Prior-2
has very strong support for loud equal-mass sources that emit louder signals,
which are, in turn, strongly favoured by the distance prior. Smaller
inclinations are also preferred for the same reason. On the other hand, since
Prior-2 does not prefer low-$Q$ signals, which generally leads to slightly
weaker sources that are then punished by the luminosity closer to the true
injection values and, therefore, a larger likelihood. However, this difference
in likelihood is not enough to overcome the difference in prior support.

The above situation changes when $Q$ is decreased, and the increasing
contribution of higher-order modes can provide information about the
inclination of the source. On the one hand, we note that for decreasing mass
ratio, the true injection parameters are progressively displaced from the
$90\%$ CI of Prior-2, which prefers distant and louder mass-symmetric systems
at smaller inclinations. This is particularly true for small inclinations, for
which the signal contains weak higher harmonics. In contrast, Prior-4 always
returns $90\%$ CI, including the true values and larger likelihoods. Even in
this situation, the Bayes Factors reveal no clear preference for any of the
models, as the larger likelihood obtained by Prior-4 cannot compensate for the
strong support that the distance prior shows for the mass ratios preferred by
Prior-2. Unlike in the above paragraph, this does not happen for the largest
inclinations, for which the higher-order harmonic modes present in the
injections do provide strong information about the inclination angle. This
way, the maximal likelihood obtained by Prior-4 is enough to overcome the
preference of Prior-2 for high-$Q$ systems, which are further favoured by the
distance prior.

Summarising, we have explicitly shown here how a perfectly reasonable prior, a
standard one in gravitational-wave data analysis, can lead to biased parameter
estimates and still be statistically preferred with respect to a second prior
giving unbiased estimates.

##  V Summary

This study investigates the influence of various mass prior choices on
Bayesian parameter inference of quasi-circular intermediate-mass black hole
binaries. We employ two waveform models, IMRPhenomXAS and IMRPhenomXHM, to
simulate and infer binary parameters, focusing on their detector-frame total
mass and mass ratio. Our results reveal that the choice of mass prior has a
substantial influence on the recovery of parameters for systems with a network
optimal signal-to-noise ratio of 20. Interestingly, changing the network
optimal SNR of the systems to 15 or 30 does not alter the qualitative nature
of our results.

When using the quadrupole-only waveform IMRPhenomXAS, we find that Prior-4,
which is flat in $M_{T}(1+z)$ and $1/Q$, recovers the injected mass parameters
better when the simulated signal corresponds to a short-lived mass-asymmetric
binary. On the other hand, Prior-1, 2 and 3, which are flat in mass ratio,
perform relatively better when recovering nearly equal mass binaries with
component spins aligned relative to each other and the orbital angular
momentum.

We have also shown that when recovering IMRPhenomXHM injections with the same
waveform model, Prior-4 yields better recovery of $M_{T}(1+z)$ and $Q$,
especially when the system is mass-asymmetric and is at a low inclination.
However, as the relative strength of the higher-order modes increases with
inclination or mass asymmetry, the accuracy of inference with conventional
prior choices, such as Prior-2, which is flat in component masses, is better.
Moreover, the correlation between the mass ratio and detector-frame total mass
decreases with increasing higher-order mode content.

Lastly, we highlight a noteworthy observation: the redshift posterior peaks
closer to the true value when using Prior-4, irrespective of the binary
configuration used in the study, implying that this prior provides a more
robust choice. Since accurate luminosity distance and/or redshift inference
can impact efforts to measure source properties of black holes, tests of
general relativity, and the study of secondary parameters, such as the
remnant’s recoil kick, it is recommended to use Prior-4 while analysing such
short-lived signals. While an exhaustive exploration of distance measurability
is beyond the scope of our current study, its importance beckons for future
investigations.

In conclusion, our findings suggest that it is essential to consider multiple
mass prior choices when analysing intermediate-mass black hole binary events
rather than relying on a single prior. Model selection techniques, such as the
Bayes Factor, can then be used with maximum likelihood to evaluate which prior
better explains the observation. While previous works have given indications
of this either when analysing real events [39, 40, 41], or when using a
broader set of injections [70], we explicitly demonstrate this in a controlled
setup.

## Acknowledgements

We thank Simon Stevenson, Thomas Dent and Johann Fernandes for their detailed
comments and valuable suggestions. KC acknowledges the MHRD, the Government of
India, for the fellowship support. AP’s research is supported by SERB-Power
fellowship grant SPF/2021/000036, DST, India. JCB is funded by a fellowship
from the “la Caixa” Foundation (ID100010434) and the European Union’s Horizon
2020 research and innovation programme under the Marie Skłodowska-Curie grant
agreement No 847648. The fellowship code is LCF/BQ/PI20/11760016. JCB is also
supported by the research grant PID2020-118635GB-I00 from the Spain-Ministerio
de Ciencia e Innovación. JAF is supported by the Spanish Agencia Estatal de
Investigación (PGC2018-095984-B-I00) and by the Generalitat Valenciana
(PROMETEO/2019/071). The authors are grateful for the computational resources
provided by the LIGO Laboratory and supported by the National Science
Foundation Grants No. PHY-0757058 and No. PHY-0823459. We are grateful for the
computational resources provided by Cardiff University and funded by an STFC
grant supporting UK Involvement in the Operation of Advanced LIGO. This
document has been assigned the LIGO document number LIGO-P2300217.

## References

  * 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._ [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] .
  * Nitz  _et al._ [2019] A. H. Nitz, C. Capano, A. B. Nielsen, S. Reyes, R. White, D. A. Brown, and B. Krishnan, 1-OGC: The first open gravitational-wave catalog of binary mergers from analysis of public Advanced LIGO data, Astrophys. J. 872, 195 (2019), arXiv:1811.01921 [gr-qc] .
  * Venumadhav  _et al._ [2020] T. Venumadhav, B. Zackay, J. Roulet, L. Dai, and M. Zaldarriaga, New binary black hole mergers in the second observing run of Advanced LIGO and Advanced Virgo, Phys. Rev. D 101, 083030 (2020), arXiv:1904.07214 [astro-ph.HE] .
  * Nitz  _et al._ [2020] A. H. Nitz, T. Dent, G. S. Davies, S. Kumar, C. D. Capano, I. Harry, S. Mozzon, L. Nuttall, A. Lundgren, and M. Tápai, 2-OGC: Open Gravitational-wave Catalog of binary mergers from analysis of public Advanced LIGO and Virgo data, Astrophys. J. 891, 123 (2020), arXiv:1910.05331 [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] .
  * Nitz  _et al._ [2021] A. H. Nitz, C. D. Capano, S. Kumar, Y.-F. Wang, S. Kastha, M. Schäfer, R. Dhurkunde, and M. Cabero, 3-OGC: Catalog of Gravitational Waves from Compact-binary Mergers, Astrophys. J. 922, 76 (2021), arXiv:2105.09151 [astro-ph.HE] .
  * 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 (2021b), arXiv:2111.03606 [gr-qc] .
  * Nitz  _et al._ [2023] A. H. Nitz, S. Kumar, Y.-F. Wang, S. Kastha, S. Wu, M. Schäfer, R. Dhurkunde, and C. D. Capano, 4-OGC: Catalog of Gravitational Waves from Compact Binary Mergers, Astrophys. J. 946, 59 (2023), arXiv:2112.06878 [astro-ph.HE] .
  * Olsen  _et al._ [2022] S. Olsen, T. Venumadhav, J. Mushkin, J. Roulet, B. Zackay, and M. Zaldarriaga, New binary black hole mergers in the LIGO-Virgo O3a data, Phys. Rev. D 106, 043009 (2022), arXiv:2201.02252 [astro-ph.HE] .
  * Abbott  _et al._ [2016a] B. P. Abbott _et al._ (LIGO Scientific, Virgo), Tests of general relativity with GW150914, Phys. Rev. Lett. 116, 221101 (2016a), [Erratum: Phys.Rev.Lett. 121, 129902 (2018)], arXiv:1602.03841 [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._ [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] .
  * Roulet  _et al._ [2020] J. Roulet, T. Venumadhav, B. Zackay, L. Dai, and M. Zaldarriaga, Binary Black Hole Mergers from LIGO/Virgo O1 and O2: Population Inference Combining Confident and Marginal Events, Phys. Rev. D 102, 123022 (2020), arXiv:2008.07014 [astro-ph.HE] .
  * Abbott  _et al._ [2023] 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 (2023), arXiv:2111.03634 [astro-ph.HE] .
  * Abbott  _et al._ [2021c] 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 (2021c), arXiv:2010.14529 [gr-qc] .
  * Abbott  _et al._ [2021d] R. Abbott _et al._ (LIGO Scientific, VIRGO, KAGRA), Tests of General Relativity with GWTC-3 (2021d), arXiv:2112.06861 [gr-qc] .
  * Ghosh  _et al._ [2016] A. Ghosh _et al._ , Testing general relativity using golden black-hole binaries, Phys. Rev. D 94, 021101 (2016), arXiv:1602.02453 [gr-qc] .
  * Ghosh  _et al._ [2018] A. Ghosh, N. K. Johnson-Mcdaniel, A. Ghosh, C. K. Mishra, P. Ajith, W. Del Pozzo, C. P. L. Berry, A. B. Nielsen, and L. London, Testing general relativity using gravitational wave signals from the inspiral, merger and ringdown of binary black holes, Class. Quant. Grav. 35, 014002 (2018), arXiv:1704.06784 [gr-qc] .
  * Abbott  _et al._ [2016b] B. P. Abbott _et al._ (LIGO Scientific, Virgo), Observation of Gravitational Waves from a Binary Black Hole Merger, Phys. Rev. Lett. 116, 061102 (2016b), arXiv:1602.03837 [gr-qc] .
  * Bekenstein [1973] J. D. Bekenstein, Gravitational-Radiation Recoil and Runaway Black Holes, Astrophys. J. 183, 657 (1973).
  * Gonzalez  _et al._ [2007] J. A. Gonzalez, U. Sperhake, B. Bruegmann, M. Hannam, and S. Husa, Total recoil: The Maximum kick from nonspinning black-hole binary inspiral, Phys. Rev. Lett. 98, 091101 (2007), arXiv:gr-qc/0610154 .
  * Lousto and Zlochower [2008] C. O. Lousto and Y. Zlochower, Further insight into gravitational recoil, Phys. Rev. D 77, 044028 (2008), arXiv:0708.4048 [gr-qc] .
  * Lousto  _et al._ [2012] C. O. Lousto, Y. Zlochower, M. Dotti, and M. Volonteri, Gravitational Recoil From Accretion-Aligned Black-Hole Binaries, Phys. Rev. D 85, 084015 (2012), arXiv:1201.1923 [gr-qc] .
  * Lousto and Zlochower [2013] C. O. Lousto and Y. Zlochower, Nonlinear Gravitational Recoil from the Mergers of Precessing Black-Hole Binaries, Phys. Rev. D 87, 084027 (2013), arXiv:1211.7099 [gr-qc] .
  * Calderón Bustillo  _et al._ [2018] J. Calderón Bustillo, J. A. Clark, P. Laguna, and D. Shoemaker, Tracking black hole kicks from gravitational wave observations, Phys. Rev. Lett. 121, 191102 (2018), arXiv:1806.11160 [gr-qc] .
  * Varma  _et al._ [2020] V. Varma, M. Isi, and S. Biscoveanu, Extracting the Gravitational Recoil from Black Hole Merger Signals, Phys. Rev. Lett. 124, 101104 (2020), arXiv:2002.00296 [gr-qc] .
  * Varma  _et al._ [2022] V. Varma, S. Biscoveanu, T. Islam, F. H. Shaik, C.-J. Haster, M. Isi, W. M. Farr, S. E. Field, and S. Vitale, Evidence of Large Recoil Velocity from a Black Hole Merger Signal, Phys. Rev. Lett. 128, 191102 (2022), arXiv:2201.01302 [astro-ph.HE] .
  * Calderón Bustillo  _et al._ [2022] J. Calderón Bustillo, S. H. W. Leong, and K. Chandra, GW190412: measuring a black-hole recoil direction through higher-order gravitational-wave modes (2022), arXiv:2211.03465 [gr-qc] .
  * Abbott  _et al._ [2020a] R. Abbott _et al._ (LIGO Scientific, Virgo), GW190521: A Binary Black Hole Merger with a Total Mass of $150M_{\odot}$, Phys. Rev. Lett. 125, 101102 (2020a), arXiv:2009.01075 [gr-qc] .
  * Fishbach and Holz [2020] M. Fishbach and D. E. Holz, Minding the gap: GW190521 as a straddling binary, Astrophys. J. Lett. 904, L26 (2020), arXiv:2009.05472 [astro-ph.HE] .
  * Nitz and Capano [2021] A. H. Nitz and C. D. Capano, GW190521 may be an intermediate mass ratio inspiral, Astrophys. J. Lett. 907, L9 (2021), arXiv:2010.12558 [astro-ph.HE] .
  * Romero-Shaw  _et al._ [2020a] I. M. Romero-Shaw, P. D. Lasky, E. Thrane, and J. C. Bustillo, GW190521: orbital eccentricity and signatures of dynamical formation in a binary black hole merger signal, Astrophys. J. Lett. 903, L5 (2020a), arXiv:2009.04771 [astro-ph.HE] .
  * Gayathri  _et al._ [2022] V. Gayathri, J. Healy, J. Lange, B. O’Brien, M. Szczepanczyk, I. Bartos, M. Campanelli, S. Klimenko, C. O. Lousto, and R. O’Shaughnessy, Eccentricity estimate for black hole mergers with numerical relativity simulations, Nature Astron. 6, 344 (2022), arXiv:2009.05461 [astro-ph.HE] .
  * Calderón Bustillo  _et al._ [2021a] J. Calderón Bustillo, N. Sanchis-Gual, A. Torres-Forné, J. A. Font, A. Vajpeyi, R. Smith, C. Herdeiro, E. Radu, and S. H. W. Leong, GW190521 as a Merger of Proca Stars: A Potential New Vector Boson of $8.7\times 10^{-13}$ eV, Phys. Rev. Lett. 126, 081101 (2021a), arXiv:2009.05376 [gr-qc] .
  * Calderón Bustillo  _et al._ [2021b] J. Calderón Bustillo, N. Sanchis-Gual, A. Torres-Forné, and J. A. Font, Confusing Head-On Collisions with Precessing Intermediate-Mass Binary Black Hole Mergers, Phys. Rev. Lett. 126, 201101 (2021b), arXiv:2009.01066 [gr-qc] .
  * Gamba  _et al._ [2023] R. Gamba, M. Breschi, G. Carullo, S. Albanesi, P. Rettegno, S. Bernuzzi, and A. Nagar, GW190521 as a dynamical capture of two nonspinning black holes, Nature Astron. 7, 11 (2023), arXiv:2106.05575 [gr-qc] .
  * Estellés  _et al._ [2022] H. Estellés _et al._ , A Detailed Analysis of GW190521 with Phenomenological Waveform Models, Astrophys. J. 924, 79 (2022), arXiv:2105.06360 [gr-qc] .
  * Olsen  _et al._ [2021] S. Olsen, J. Roulet, H. S. Chia, L. Dai, T. Venumadhav, B. Zackay, and M. Zaldarriaga, Mapping the likelihood of GW190521 with diverse mass and spin priors, Phys. Rev. D 104, 083036 (2021), arXiv:2106.13821 [astro-ph.HE] .
  * Bustillo  _et al._ [2022] J. C. Bustillo, N. Sanchis-Gual, S. H. W. Leong, K. Chandra, A. Torres-Forne, J. A. Font, C. Herdeiro, E. Radu, I. C. F. Wong, and T. G. F. Li, Searching for vector boson-star mergers within ligo-virgo intermediate-mass black-hole merger candidates (2022), arXiv:2206.02551 [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] .
  * Smith  _et al._ [2020] R. J. E. Smith, G. Ashton, A. Vajpeyi, and C. Talbot, Massively parallel Bayesian inference for transient gravitational-wave astronomy, Mon. Not. Roy. Astron. Soc. 498, 4492 (2020), arXiv:1909.11873 [gr-qc] .
  * Thrane and Talbot [2019] E. Thrane and C. Talbot, An introduction to bayesian inference in gravitational-wave astronomy: parameter estimation, model selection, and hierarchical models, Publications of the Astronomical Society of Australia 36, e010 (2019).
  * Romero-Shaw  _et al._ [2020b] 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 (2020b), arXiv:2006.00714 [astro-ph.IM] .
  * Christensen and Meyer [2022] N. Christensen and R. Meyer, Parameter estimation with gravitational waves, Rev. Mod. Phys. 94, 025001 (2022), arXiv:2204.04449 [gr-qc] .
  * Kass and Raftery [1995] R. E. Kass and A. E. Raftery, Bayes factors, Journal of the American Statistical Association 90, 773 (1995).
  * Racine [2008] E. Racine, Analysis of spin precession in binary black hole systems including quadrupole-monopole interaction, Phys. Rev. D 78, 044021 (2008), arXiv:0803.1820 [gr-qc] .
  * Ajith  _et al._ [2011] P. Ajith _et al._ , Inspiral-merger-ringdown waveforms for black-hole binaries with non-precessing spins, Phys. Rev. Lett. 106, 241101 (2011), arXiv:0909.2867 [gr-qc] .
  * Santamaria  _et al._ [2010] L. Santamaria _et al._ , Matching post-Newtonian and numerical relativity waveforms: systematic errors and a new phenomenological model for non-precessing black hole binaries, Phys. Rev. D 82, 064016 (2010), arXiv:1005.3306 [gr-qc] .
  * Abbott  _et al._ [2020b] R. Abbott _et al._ (LIGO Scientific, VIRGO, KAGRA), Gw190521 parameter estimation samples and figure data, https://dcc.ligo.org/LIGO-P2000158/public (2020b).
  * Pratten  _et al._ [2020] G. Pratten, S. Husa, C. Garcia-Quiros, M. Colleoni, A. Ramos-Buades, H. Estelles, and R. Jaume, Setting the cornerstone for a family of models for gravitational waves from compact binaries: The dominant harmonic for nonprecessing quasicircular black holes, Phys. Rev. D 102, 064001 (2020), arXiv:2001.11412 [gr-qc] .
  * García-Quirós  _et al._ [2020] C. García-Quirós, M. Colleoni, S. Husa, H. Estellés, G. Pratten, A. Ramos-Buades, M. Mateu-Lucena, and R. Jaume, Multimode frequency-domain model for the gravitational wave signal from nonprecessing black-hole binaries, Phys. Rev. D 102, 064002 (2020), arXiv:2001.10914 [gr-qc] .
  * Abbott  _et al._ [2022] R. Abbott _et al._ (LIGO Scientific, VIRGO, KAGRA), Noise curves used for simulations in the update of the observing scenarios paper, https://dcc.ligo.org/LIGO-T2000012/public (2022).
  * Graff  _et al._ [2015] P. B. Graff, A. Buonanno, and B. S. Sathyaprakash, Missing Link: Bayesian detection and measurement of intermediate-mass black-hole binaries, Phys. Rev. D 92, 022002 (2015), arXiv:1504.04766 [gr-qc] .
  * Usman  _et al._ [2019] S. A. Usman, J. C. Mills, and S. Fairhurst, Constraining the Inclinations of Binary Mergers from Gravitational-wave Observations, Astrophys. J. 877, 82 (2019), arXiv:1809.10727 [gr-qc] .
  * Boyle  _et al._ [2019] M. Boyle _et al._ , The SXS Collaboration catalog of binary black hole simulations, Class. Quant. Grav. 36, 195006 (2019), arXiv:1904.04831 [gr-qc] .
  * Healy  _et al._ [2019] J. Healy, C. O. Lousto, J. Lange, R. O’Shaughnessy, Y. Zlochower, and M. Campanelli, Second RIT binary black hole simulations catalog and its application to gravitational waves parameter estimation, Phys. Rev. D 100, 024021 (2019), arXiv:1901.02553 [gr-qc] .
  * 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), arXiv:1904.02180 [astro-ph.IM] .
  * Farr [2014] W. Farr, Marginalisation of the time parameter in gravitational wave parameter estimation, https://dcc.ligo.org/LIGO-T1400460/public (2014).
  * Singer and Price [2016] L. P. Singer and L. R. Price, Rapid Bayesian position reconstruction for gravitational-wave transients, Phys. Rev. D 93, 024013 (2016), arXiv:1508.03634 [gr-qc] .
  * Veitch and Del Pozzo [2013] J. Veitch and W. Del Pozzo, Analytic marginalisation of phase parameter, https://dcc.ligo.org/LIGO-T1300326/public (2013).
  * 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] .
  * Lin [1991] J. Lin, Divergence measures based on the shannon entropy, IEEE Transactions on Information Theory 37, 145 (1991).
  * MacKay [2005] D. J. MacKay,  _Information Theory, Inference and Learning Algorithms_ (Cambridge University Press, 2005).
  * Campanelli  _et al._ [2006] M. Campanelli, C. O. Lousto, and Y. Zlochower, Spinning-black-hole binaries: The orbital hang up, Phys. Rev. D 74, 041501 (2006), arXiv:gr-qc/0604012 .
  * Healy and Lousto [2017] J. Healy and C. O. Lousto, Remnant of binary black-hole mergers: New simulations and peak luminosity studies, Phys. Rev. D 95, 024037 (2017), arXiv:1610.09713 [gr-qc] .
  * Hofmann  _et al._ [2016] F. Hofmann, E. Barausse, and L. Rezzolla, The final spin from binary black holes in quasi-circular orbits, Astrophys. J. Lett. 825, L19 (2016), arXiv:1605.01938 [gr-qc] .
  * Jiménez-Forteza  _et al._ [2017] X. Jiménez-Forteza, D. Keitel, S. Husa, M. Hannam, S. Khan, and M. Pürrer, Hierarchical data-driven approach to fitting numerical relativity data for nonprecessing binary black holes with an application to final spin and radiated energy, Phys. Rev. D 95, 064024 (2017), arXiv:1611.00332 [gr-qc] .
  * 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 (2023), arXiv:2308.03250 [gr-qc] .
  * Shaik  _et al._ [2020] F. H. Shaik, J. Lange, S. E. Field, R. O’Shaughnessy, V. Varma, L. E. Kidder, H. P. Pfeiffer, and D. Wysocki, Impact of subdominant modes on the interpretation of gravitational-wave signals from heavy binary black hole systems, Phys. Rev. D 101, 124054 (2020), arXiv:1911.02693 [gr-qc] .
  * Kalaghatgi  _et al._ [2020] C. Kalaghatgi, M. Hannam, and V. Raymond, Parameter estimation with a spinning multimode waveform model, Phys. Rev. D 101, 103004 (2020), arXiv:1909.10010 [gr-qc] .
  * Fowler and Hoyle [1964] W. A. Fowler and F. Hoyle, Neutrino Processes and Pair Formation in Massive Stars and Supernovae, Astrophys. J. Suppl. 9, 201 (1964).
  * Barkat  _et al._ [1967] Z. Barkat, G. Rakavy, and N. Sack, Dynamics of Supernova Explosion Resulting from Pair Formation, Phys. Rev. Lett. 18, 379 (1967).
  * Bond  _et al._ [1984] J. R. Bond, W. D. Arnett, and B. J. Carr, The Evolution and fate of Very Massive Objects, Astrophys. J. 280, 825 (1984).
  * Woosley and Heger [2021] S. E. Woosley and A. Heger, The Pair-Instability Mass Gap for Black Holes, Astrophys. J. Lett. 912, L31 (2021), arXiv:2103.07933 [astro-ph.SR] .
  * Farmer  _et al._ [2019] R. Farmer, M. Renzo, S. E. de Mink, P. Marchant, and S. Justham, Mind the gap: The location of the lower edge of the pair instability supernovae black hole mass gap (2019), arXiv:1910.12874 [astro-ph.SR] .

  *[IGWN]: International Gravitational-Wave Observatory Network
  *[CI]: credible intervals