# Optimizing the low-latency localization of gravitational waves

Pierre-Alexandre Duverne  Université Paris-Saclay, CNRS/IN2P3, IJCLab, 91405
Orsay, France  Université Paris Cité, CNRS, Astroparticule et Cosmologie,
F-75013 Paris, France<EMAIL_ADDRESS>Stéphanie Hoang  Université Paris-
Saclay, CNRS/IN2P3, IJCLab, 91405 Orsay, France  Tito Dal Canton  Université
Paris-Saclay, CNRS/IN2P3, IJCLab, 91405 Orsay, France  Sarah Antier  Artemis,
Observatoire de la Côte d’Azur, Université Côte d’Azur, Boulevard de
l’Observatoire, 06304 Nice, France  Nicolas Arnaud  Université Paris-Saclay,
CNRS/IN2P3, IJCLab, 91405 Orsay, France  European Gravitational Observatory
(EGO), I-56021 Cascina, Pisa, Italy  Patrice Hello  Université Paris-Saclay,
CNRS/IN2P3, IJCLab, 91405 Orsay, France  Francesco Pannarale  Dipartimento di
Fisica, Università di Roma “Sapienza”, Piazzale A. Moro 5, I-00185, Roma,
Italy  INFN Sezione di Roma, Piazzale A. Moro 5, I-00185, Roma, Italy

(August 28, 2024)

###### Abstract

Gravitational-wave data from interferometric detectors like LIGO, Virgo and
KAGRA is routinely analyzed by rapid matched-filtering algorithms to detect
compact binary merger events and rapidly infer their spatial position, which
enables the discovery of associated non-GW transients like GRB 170817A and
AT2017gfo. One of the critical requirements for finding such counterparts is
that the rapidly inferred sky location, usually performed by the Bayestar
algorithm, is correct. The reliability of this data product relies on various
assumptions and a tuning parameter in Bayestar, which we investigate in this
paper in the context of PyCBC Live, one of the rapid search algorithms used by
LIGO, Virgo and KAGRA. We perform simulations of compact binary coalescence
signals recovered by PyCBC Live and localized by Bayestar, under various
configurations with different balances between simplicity and realism, and we
test the resulting sky localizations for consistency based on the widely-used
PP plot. We identify some aspects of the search configuration which drive the
optimal setting of Bayestar’s tuning parameter, in particular the properties
of the template bank used for matched filtering. On the other hand, we find
that this parameter does not depend strongly on the nonstationary and non-
Gaussian properties of the detector noise.

##  I Introduction

The rapid detection of gravitational-wave (GW) signals from compact binary
coalescences (CBCs) in data from existing observatories (LIGO [1], Virgo [2]
and KAGRA [3]) is done by running dedicated low-latency searches such as PyCBC
Live [4, 5], GstLAL [6], MBTAOnline [7], SPIIR [8] or cWB online [9]. These
algorithms enable the production of rapid GW alerts usable by other
instruments to search for counterparts. During the O3 observing run [10, 11],
that happened between April 2019 and March 2020, these alerts were publicly
distributed with notices and circulars sent via the General Coordinates
Network (GCN) [12], as it is the case during the current O4 run. Once low-
latency search pipelines identify a candidate CBC signal with a sufficiently
low false-alarm rate, the spatial location of the event is estimated by a
dedicated algorithm called Bayestar [13], which uses a Bayesian framework to
rapidly infer the joint posterior distribution for the source’s sky location
(right ascension $\alpha$ and declination $\delta$) and the luminosity
distance. Other telescopes, or telescope networks, can then use these inferred
quantities to follow up GW alerts. At longer time scales, typically at least
tens of minutes after the merger time, a full parameter inference is performed
on the detected signals [14, 15], which can be used to update the Bayestar
result.

The rapid distribution of this information has played a key role in the
discovery of the GW170817 event, the first detection of a binary neutron star
merger (BNS) associated with the short gamma-ray burst (GRB) GRB 170817A and
the optical counterpart (kilonova) AT2017gfo less than 11 hours after the
merger [16, 17, 18, 19, 20, 21, 22, 23]. For this particular event, the 90%
credible area of the sky location was about $28$ deg2. In general, however,
this area can be much larger, depending on the signal-to-noise ratio (SNR) of
the GW signal and on the status of the detector network at the time of the
signal (one, two or three interferometers observing). For instance, the median
size of the 90% credible region was about $1000$ deg2 for Binary Black Hole
(BBH) alerts during the O3 run, and predictions made before the O4 run give
similar orders of magnitude ($>1000$ deg2) [24, 25].

Since the success of follow-up campaigns depends on the availability of prompt
sky localizations of GW events, confidence in Bayestar results is crucial,
especially for 1) well-localised events and 2) events containing at least one
neutron star à la GW170817. Bayestar’s rapidity (typically $\sim 1$ s) rests
on a simplified treatment of the full parameter-estimation problem, in
particular on the decoupling of the extrinsic parameters describing a CBC
signal (spatial location and orientation) and its intrinsic parameters
(typically masses and spins). Although this simplification has proven to be
extremely useful, in principle, it could be associated with biases in the
spatial localization.

For instance, Bayestar’s self-consistency has been tested for candidates
produced by the GstLAL algorithm [13] and revealed a systematic
underestimation of the credible area of the sky location. This underestimation
was corrected by introducing a correction parameter that we denote $\xi$ for
this work. The optimal value for $\xi$ was then found using a signal
simulation campaign and the GstLAL pipeline, resulting in $0.83$, which is
used to date irrespective of which search pipeline produced the candidate. On
the other hand, [4] pointed out a possible overestimation of the credible area
for candidates reported by the PyCBC Live pipeline, raising the question of
whether the optimal value of $\xi$ might depend on the details of how each GW
candidate is generated. Overestimating the credible area is not a priori
problematic, as the GW source is still inside the sky area covered by the
telescopes. However, it makes the electromagnetic follow-up more demanding,
wasting telescope time. An underestimation of the uncertainty in the sky
location is even less desirable, as the follow-up risks missing the
counterpart if the GW source is outside of the high-probability region.

Two points were not directly addressed in [13]. First, what specific aspect(s)
of the entire detection and inference process introduce the need for such a
correction? [13] lists a number of plausible possibilities, for instance,
certain properties of the template waveforms used for matched-filtering, or
the basic approximation made by Bayestar. Second, does the optimal value
depend on which search pipeline produced the candidate, or is $0.83$
sufficient for all the existing search pipelines? The latter question is
further raised by the PP-plot in Figure 5 of [4], which was constructed using
a PyCBC Live analysis, and shows a curve above the expected diagonal,
consistent with an overestimation of the uncertainty.

Motivated by these considerations, in this paper we present a detailed
exploration of possible sources of bias in the sky localization produced by
Bayestar, focusing on CBC candidates produced by the PyCBC Live pipeline and
on the behavior of the $\xi$ parameter described in [13].

##  II Low latency GW search and localization

The most sensitive search algorithms for GWs from stellar-mass CBCs are based
on a technique called _matched filtering_ that consists in cross-correlating
the data with a CBC template waveform to compute the signal-to-noise ratio.
The template is taken from a pre-computed collection of CBC waveforms called a
template bank. The bank is designed to efficiently cover the parameter space
of the binary’s four intrinsic parameters (the masses and spins of the two
merging objects) under various simplifying assumptions, namely that the binary
orbit is quasi-circular, that the spins are aligned with the orbital angular
momentum, that matter effects do not significantly distort the waveform, and
that the SNR is dominated by the quadrupolar mode of the radiation.

Although there are various implementations of matched filtering, for this
work, we focus on the algorithm called PyCBC Live [4, 5], where matched
filtering is computed in the frequency domain. It is based on the FindChirp
algorithm [26], where the GW data $s$ are correlated with a template waveform
$h$ in the frequency domain to build the complex SNR time series

$z(t)=\frac{4}{\sqrt{\langle
h|h\rangle}}\int^{\infty}_{0}\frac{\tilde{s}(f)\tilde{h}^{*}(f)}{S_{n}(f)}e^{2\pi
ift}\,df$ (1)

where the $\langle.|.\rangle$ operator corresponds to the inner product

$\langle
a|b\rangle=4\int^{\infty}_{0}\frac{\tilde{a}(f)\tilde{b}^{*}(f)}{S_{n}(f)}\,df.$
(2)

The normalisation $S_{n}(f)$ corresponds to the one-sided power spectral
density (PSD) of the noise around time $t$. In PyCBC Live, the PSD is
estimated using the Welch method [27] with a median averaging. A trigger is
identified in the data from a single detector whenever a local maximum in the
magnitude of the SNR time series $\rho(t)=|z(t)|$ crosses a threshold of
$4.5$. The time and template associated with the local maximum are recorded as
part of the trigger. PyCBC Live can then use the set of triggers from multiple
detectors to either identify single-detector candidates (signals detectable in
a single detector only) or coincident candidates (signals detectable in at
least two detectors). A significance (false-alarm rate, FAR) is calculated in
different ways in the two cases. Candidates with a FAR lower than one per two
hours are transmitted to the Gravitational-wave Candidate Event DataBase
(GraceDB) [28].

Bayestar [13] is then used to infer the spatial location of each candidate
uploaded to GraceDB, in terms of the posterior probability for the three-
dimensional position of the source in the universe, under a Bayesian
framework. Bayestar takes as input the point estimates of the masses and spins
from the template reported by PyCBC Live and associated with the candidate, a
portion of the complex SNR time series around the estimated merger time, and
the sensitivities of each involved detector at the time of the candidate.
Based on this information, Bayestar computes the likelihood and infers the
source position. Using the point estimates of masses and spins is an
approximation that avoids having to jointly infer the spatial location and the
intrinsic parameters of the binary (component masses and spins), making the
algorithm rapid. The resulting spatial posterior is produced within seconds,
enabling rapid electromagnetic follow-up of the candidate. The approximation
is justified because the uncertainties in the intrinsic and extrinsic
parameters are typically very weakly correlated; this can be understood
intuitively considering that the spatial location is mainly determined by the
relative arrival times, phases and amplitudes at the different detectors,
while mass and (aligned) spin parameters are mainly determined from the
temporal evolution of the GW waveform.

For this work, we focus on the two-dimensional posterior probability obtained
by marginalizing the three-dimensional one over distance, which is the
quantity of immediate interest for deciding where to point telescopes. We will
refer to this quantity as the _skymap_ for simplicity.

##  III Motivation and Methodology

A common way to test the reliability of parameter inference in GW astronomy is
to use a test based on the Percentile-Percentile plot (PP-plot) [29, 13] and
this is what we use in this study. We simulate a population of CBC signals
(_injections_), add them to noise (either simulated or real, as discussed
case-by-case later), produce a candidate detection for each signal, and use
Bayestar to produce the skymaps for each candidate. We then perform a
cumulative sum of the probability contained in the skymap pixels, starting
from the highest-probability pixel and proceeding in order of descending
probability, and stopping at the pixel that contains the true position. The
cumulative sum corresponds to the cumulative probability in the smallest area
that contains the true location, and is referred to as the _search
probability_ in the following. Repeating this procedure for the whole
population of simulated CBC signals, we then make a cumulative histogram of
the search probability. This histogram is called a percentile-percentile plot
and allows one to evaluate whether the localization is self-consistent. In the
case of a consistent localization, the histogram is expected to be diagonal,
representing the fact that the search probability distribution is uniform. In
other words, one should find $p$% of the simulated signals to have their true
location within the $p$% credible region.

In the following sections, we will perform a sequence of simulations to test
some of the possible origins of this correction, starting with the simplest
controlled simulations and progressing towards full end-to-end simulations
with a realistic configuration of PyCBC Live.

###  III.1 Choices of signals and detectors

The choices of the parameters of the simulated CBC signals are summarised in
Table 1. For BNS signals, we use inspiral-only analytical waveforms based on
the SpinTaylorT4 approximant [30]. For BBH and NSBH signals we use effective-
one-body (EOB) inspiral-merger-ringdown waveforms based on the SEOBNRv4
approximant [31]. For all signals, the spins of the coalescing objects are
assumed to be aligned with the orbital angular momentum.

Parameter | Distribution
---|---
Polarisation angle | $\mathcal{U}([0,2\pi])$
Phase at coalescence | $\mathcal{U}([0,2\pi])$
$\cos{\iota}$ | $\mathcal{U}([-1,1])$
BNS position |  $\mathcal{U}(\mathbb{S}(0,[40,250]$Mpc))
NS mass |  $\mathcal{N}(\mu=1.4$M${}_{\odot},\,\sigma=0.1$M${}_{\odot})$
Aligned comp. of NS spin | $\mathcal{N}(\mu=0,\,\sigma=0.01)$
BBH position |  $\mathcal{U}(\mathbb{S}(0,[100,5000]$Mpc))
BH mass |  $\mathcal{U}([5,100]$M⊙)
Aligned comp. of BH spin | $\mathcal{N}(\mu=0,\,\sigma=0.1)$
Table 1:  Parameters of the simulated CBC signals. $\iota$ is the orbital
inclination of the binary. The notation $\mathcal{U}([a,b])$ indicates a
uniform distribution over the interval $[a,b]$. $\mathcal{N}(\mu,\,\sigma)$
indicates a normal distribution with mean $\mu$ and standard deviation
$\sigma$. The $\mathcal{U}(\mathbb{S}(0,[a,b]$)) notation corresponds to a
uniform distribution in volume in a shell centered on the origin, with inner
radius $a$ and outer radius $b$. Masses are expressed in the detector frame,
i.e. redshifted.

The simulated waveforms are injected in 600 s data segments of Gaussian and
stationary noise, with the exception of Section V where we will use real
noise. The simulated noise has amplitude spectral densities (ASDs) shown in
Figure 1. We use the same ASD for the two LIGO interferometers and a
significantly less sensitive ASD for Virgo, so as to have a network
configuration with relative sensitivities similar to O3 and O4. The simulated
interferometers are overall less sensitive than the O3 instruments,
particularly Virgo in the low-frequency region. However, this is not expected
to significantly impact our conclusions about the localization consistency.

Figure 1: Noise ASD models used to generate the simulated data (dash-dotted
lines). These are specifically aLIGOMidLowSensitivityP1200087 and
AdVEarlyLowSensitivityP1200087 from the LALSimulation library [32]. For
comparison, we also show ASDs representing the three detectors’ strain
sensitivities during O3, taken from the GWTC-2 catalog [11].

In order to simulate the selection effects of CBC search pipelines, we apply a
cut on the optimal SNR of each signal before generating the skymap. The
optimal SNR corresponds to the expected SNR a GW signal would have if the
template was identical to the signal. The optimal SNRs of the individual
detectors involved in a GW detection are summed in quadrature to compute the
network’s optimal SNR. We keep only the signals with a network SNR above 8 and
with at least one individual optimal SNR above $5.5$. These cuts are chosen to
match Table 1 from [10], which summarises the network SNR of the publicly
distributed alerts during O3, and SNR $\sim 8$ is the lower bound for
distributed online GW candidate events. The signals surviving these cuts are
localised with Bayestar to generate the skymaps.

###  III.2 Baseline percentile-percentile test

As a first check, we attempt to reproduce the results in Figure 5 of [4] using
500 BNS and 500 BBH simulated following the distributions in Section III.1.

For our first simulations, we use a tool available in PyCBC, called
pycbc_make_skymap111
https://github.com/gwastro/pycbc/blob/master/bin/pycbc_make_skymap, which
simulates the PyCBC Live search pipeline in a controlled, faster and
simplified way. This tool creates a template waveform from a user-provided set
of parameters, performs the matched filtering on the simulated data,
calculates the SNR time series, formats the results and passes them to
Bayestar to produce the skymap. An important simplification is that
pycbc_make_skymap uses a single template specified by the user, thus removing
the possible effects related to a template bank. After running
pycbc_make_skymap and Bayestar, we calculate the search probability for each
signal and use its distribution to obtain the PP plot.

Apart from showing the PP plot, we also perform a Kolmogorov-Smirnov (KS) test
[33] as a quantitative check of the consistency of the localization. The null
hypothesis of our KS test is that the localization is self-consistent. As
described in Section III, under this hypothesis the search probability is
uniformly distributed in $[0,1]$ and produces a diagonal cumulative
distribution. Correspondingly, we use the KS-test to quantify the agreement of
the search probability distribution with a uniform distribution. In the
following sections, we consider that a KS-test p-value smaller than $3\times
10^{-3}$, corresponding to a $\sim$3-$\sigma$ deviation, is small enough to
reject the null hypothesis.

For the first set of simulations, we use Bayestar in its O3 configuration and
assume that the true intrinsic parameters (masses and spins) of the signals
are exactly known, i.e. we use the true parameters of each signal in
pycbc_make_skymap. The baseline PP plot for this configuration is visible in
Figure 2 with BBHs in orange and BNSs in blue. For both source types, the
curves are clearly inconsistent with the null hypothesis, and we conclude that
the sky-localization uncertainties are overestimated, which appears to confirm
the results in [4].

Figure 2: Baseline PP-plot for the simulated BNSs and BBHs. The blue curve is
the cumulative distribution of the search probability computed with the
skymaps of 500 BNSs. The orange one is for the 500 simulated BBHs. In the case
of a self-consistent analysis, the curve is expected to be diagonal. The grey
shaded areas around the diagonal are the 1, 2 and $3\sigma$ confidence bands,
showing the variability expected from the finite number of simulated signals.

##  IV Origin of the deviation

As introduced in Section III, Bayestar uses a factor $\xi$ to ensure that the
PP-plot is diagonal, and its default value $0.83$ has been found via
simulations recovered by the GstLAL pipeline. Our first variation with respect
to the baseline result is removing the correction introduced via this factor.
We set $\xi=1$ and repeat our simulations, using three separate populations of
$500$ BBH, $500$ BNS and $500$ neutron-star-black-hole (NSBH) mergers to test
whether the effect of $\xi$ could depend on the source type. The results for
these simulations are summarised in Table 2, and the PP plots are visible in
Figure 3. For those three simulated populations, the p-values returned by the
KS-test are not small enough to reject the null hypothesis, i.e. the Bayestar
analysis is self-consistent. Hence, we conclude that the deviation observed in
Figure 2 can be explained by the default value of $\xi$ being different from
$1$, and it appears that our simulations are localized in a perfectly
consistent way without the need for the $\xi$ correction.

CBC Type | p-value | Deviation
---|---|---
BNS | $0.45$ | $0.8\sigma$
BBH | $0.11$ | $1.6\sigma$
NSBH | $0.36$ | $0.9\sigma$
Table 2: Results of the KS-test for our second set of simulations, testing the
effect of setting $\xi=1$. Figure 3: PP plot for three separate populations of
GW sources after setting the $\xi$ parameter to 1. This time the deviations
are within 3-$\sigma$; consequently, we cannot reject the null hypothesis.

After this second test, we are thus left with a number of questions. What is
different between the simulation campaign done in [13] and our simulations? Do
the different values of $\xi$ depend on the fact that we assumed idealized
noise, or on differences between the GstLAL and PyCBC Live analyses, or on the
simplifications made by pycbc_make_skymap? These questions are investigated in
the following sections.

##  V Effect of real detector noise

In the first sets of injections, the signals have been injected in simulated
Gaussian and stationary noise. However, the data produced by a real GW
detector have a non-stationary character and contain instrumental transients,
commonly called glitches, that may have an influence on the localization (see
e.g. GW170817, where a loud glitch in LIGO Livingston’s data had to be removed
before any analysis [34]). Thus, it is necessary to verify if these effects
might be causing a bias in the localization, introducing the need for the
$\xi$ correction. We simulate populations of BBH and BNS in the same way as
described in Section III.1 and inject the signals into publicly available data
from O3a (April to October 2019) and O3b (November 2019 to March 2020). The
injections are made into data segments of O3 by picking a random time and
keeping it if $600$ s of data were available in the three detectors around
that time. As sensitivity and glitch rate significantly improved between the
two halves of O3, we create two separate sets of injections for O3a and O3b.

We inject 1000 BBH signals and 1000 BNS for both O3a and O3b data segments,
and we keep working with the assumption that the intrinsic parameters are
exactly known, as in Section III.2. For each interferometer, an average ASD is
estimated for each week of the O3 run. Then for each signal, the closest ASD
to the picked time is chosen for the injection. Eventually, 908 BBH and 921
BNS are loud enough to be detected and localised with O3a data. For O3b, we
have 912 BBHs and 940 BNS. The results are summarised in Table 3, and Figure 4
shows the PP plots for the injected BNS and BBH in O3b data segments.

For both O3a and O3b, we observe deviations similar to Figure 3 for
$\xi=0.83$, and much smaller deviations for $\xi=1$. The case of the BBH
injected in O3b data shows a deviation of $3.3\sigma$ for $\xi=1$, indicating
a small but statistically significant deviation based on the threshold we
decided in Section III.2. We can see, indeed, that the BBH curve sags below
the diagonal with $\xi=1$, i.e. the uncertainty is slightly underestimated.
This is not the case for BNS signals. This might indicate a small mass
dependence on how real detector noise affects the sky localization. Indeed, a
recent study, which investigates the effect of specific glitch families on sky
localizations produced via PyCBC Live and Bayestar, also draws a similar
conclusion [35].

This mass dependence, however, appears fairly small in our study, and only
detectable with O3b data. We conclude that setting $\xi=0.83$ does not appear
to be required because of the real noise properties like the presence of
glitches in the data. We thus proceed with investigating the remaining
questions using simulated Gaussian noise, which is easier to work with.

Run | CBC | $\xi$ | Deviation | Run | CBC | $\xi$ | Deviation
---|---|---|---|---|---|---|---
| BBH | $0.83$ | $\geq 10\sigma$ |  | BBH | $0.83$ | $7.0\sigma$
O3a |  | $1.0$ | $1.0\sigma$ | O3b |  | $1.0$ | $3.3\sigma$
| BNS | $0.83$ | $\geq 10\sigma$ |  | BNS | $0.83$ | $\geq 10\sigma$
|  | $1.0$ | $0.4\sigma$ |  |  | $1.0$ | $0.5\sigma$
Table 3: Results of the PP-plot and KS tests for the simulations testing
glitches and non-stationarities in O3a data on the left and O3b on the right.
Figure 4: PP-plot for BNS and BBH injection in O3b data. The blue and green
curves are for the default $\xi=0.83$ for BBH and BNS, respectively. The
orange and red curves are produced with a $\xi=1$ for BBH and BNS,
respectively.

##  VI Effect of incorrect template parameters

As stated in [13], one of the possible explanations for why $\xi$ is necessary
might be the inevitable discrepancy between the true signals and the template
waveforms reported with the detection candidates by the matched filtering
algorithms. Such a discrepancy arises mainly from two effects: the detector
noise, and the various approximations made in the algorithms themselves,
including the physics neglected in the models used for the search templates.

In our simulations up to this point, we used the exact same waveform for both
the injection and the localization, which implies assuming that the search
would report exactly the true intrinsic parameters (masses and spins) of the
sources. We now relax this assumption, and perform additional simulations for
which the intrinsic parameters of the template used for Bayestar are
randomized in a way that mimics the effect of an online low-latency search. We
focus on BNS signal only, for they are the most promising GW source for having
EM counterparts, making their low-latency localization critical. In addition,
none of the previous tests show qualitatively different behavior between BNS
and BBH systems.

We inject simulated signals in 600 s segments of Gaussian noise. The aligned-
spin parameters of the templates given to Bayestar are chosen uniformly in
$[-1,1]$, to mimic the fact that search algorithms like PyCBC Live do not
typically provide reliable point estimates of the spins. The redshifted
component masses $m_{1,2}$ of the templates are chosen by expressing them via
the mass ratio

$q=\frac{m_{1}}{m_{2}}\quad m_{1}\geq m_{2},$ (3)

and the chirp mass

$\mathcal{M}_{c}=\frac{(m_{1}m_{2})^{3/5}}{(m_{1}+m_{2})^{1/5}}.$ (4)

Similarly to spins, for the mass ratio $q$ we choose a uniformly distributed
value in $[1,3]$, which is consistent with the range of neutron star masses
commonly assumed within LIGO and Virgo, i.e. $[1,3]M_{\odot}$. Dealing with
the chirp mass $\mathcal{M}_{c}$ is more delicate, for this parameter is very
well reconstructed by an online search, in particular for BNS signals [36]. To
randomize the chirp mass in a realistic way, we use the results obtained with
an end-to-end search performed with PyCBC Live on BNS injections described
later in Section VII.5. The produced results are used to compute the standard
deviation $\sigma_{search}=7.3\times 10^{-4}$ M⊙ of the
$\mathcal{M}_{c,injection}-\mathcal{M}_{c,trigger}$ distribution that is shown
in Figure 5. Based on this, a random number is added to the injection chirp
mass following a Normal distribution with parameters $\mu=0$ and
$\sigma=\sigma_{search}$. Then Bayestar is used to compute the skymaps, and we
test several values of $\xi$, evenly spaced in $[0.7,1]$, plus the default
$0.83$ value.

Figure 5: Difference between the true ($\mathcal{M}_{c,injection}$) and
recovered $(\mathcal{M}_{c,trigger}$) redshifted chirp masses for a population
of BNS signals recovered by an end-to-end search based on PyCBC Live, as a
function of the intrinsic SNR of each signal. These results are used to
estimate the scatter in recovered chirp mass which affects the template used
for Bayestar.

The PP-plots resulting from these simulations are shown in Figure 6. Contrary
to our previous results, here $\xi<1$ is necessary to have a diagonal PP-plot;
with $\xi=1$, the PP-plot sags below the diagonal, similar to the original
observation made in [13]. This result suggests that $\xi$ effectively
compensates for the difference between the true intrinsic parameters of the
component objects, and the “best” values of those parameters identified by the
matched filtering analysis. For a more precise inspection, we plot the
deviation from the null hypothesis of the KS test as a function of $\xi$ in
Figure 7. This plot shows an optimum minimising the KS-test deviation for
$\xi=0.85$, associated with a deviation of $1.9\sigma$. This value is close to
the default setting $\xi=0.83$, showing that the simple test described in this
Section allows us to reproduce almost entirely the effect observed with an
online search.

Figure 6: PP-plots obtained after randomising the template’s intrinsic
parameters compared to their true values before the signal localization.
Figure 7:  Optimisation of the $\xi$ parameter after randomizing the intrinsic
parameters of the CBC template. The arrows show the values of $\xi$ for which
the deviation from the null hypothesis of the KS test is larger than
10$\sigma$.

##  VII End-to-end simulations with PyCBC Live

We now go one step beyond the simplified pycbc_make_skymap tool, and set up an
end-to-end simulation where the injected signals are recovered by an actual
PyCBC Live analysis similar to what has been used during the O3 run. This
gives us a more realistic configuration for optimizing the $\xi$ parameter
than what is possible with pycbc_make_skymap.

As stated earlier, we expect two factors to dominate the discrepancy between
the template and the source’s intrinsic parameters: the noise realisation
around the GW signal, which we cannot control, and the properties of the
template bank used to run the search, which we are entirely free to choose.
The possibility that $\xi$ might be related to the template bank was already
suggested in [13]. Hence we perform analyses with several template banks to
identify the parameters with the most prominent influence on $\xi$. In
particular, we start from the smallest and simplest bank possible and proceed
towards the bank actually used in the O3 run.

###  VII.1 Simulated data and analysis method

All the analyses presented in the following sections use the same dataset,
which comprises 1000 BNS signals injected in simulated Gaussian noise, with
the same distributions and parameters discussed in Section III.1. We choose
simulated Gaussian noise as opposed to real data so as to have as much control
as possible over the data, and avoid the practical complications of data gaps,
data quality issues and actual astrophysical signals already present in real
data. Based on the results of Section V, this choice is not expected to lead
to a different conclusion compared to using injections in real data. The BNS
signals are separated by $250$ s, excluding any overlap between the signals,
and ensuring the absence of a bias in the estimation of the noise ASD
performed by PyCBC Live. The total amount of simulated data represents about
three days of data.

PyCBC Live is run using a configuration equivalent to the one used for the
online search during the O3 run [5], except for the template bank, which
changes for the various tests presented in the following subsections.

The GW candidates found by the search are associated with an injection when
there is less than $40$ ms between the true and estimated coalescence times.
Then, we use Bayestar to localise the detected injections and produce the
associated skymaps for constructing the PP plots. We use various values of
$\xi$ spread in $[0.7,1]$, including the baseline $0.83$ value, and generate
skymaps for each of them. We generate a PP plot in a similar way as in Figure
6 and evaluate the deviation using the KS test for each value of $\xi$. This
process allows us to identify the optimal $\xi$ value for which the deviation
is minimal, in the same way as we did in Figure 7.

###  VII.2 Template bank matching the injections

As we understand $\xi$ as a means to compensate for the difference between the
template and true intrinsic parameters, we first test a “minimal” bank whose
templates have the same intrinsic parameters as the injected signals. Since
the waveforms in this bank are extremely similar to the injected
signals222Apart from minor differences between the SpinTaylorT4 model used for
the injections and the TaylorF2 model used for the templates., we expect the
optimal $\xi$ to be closer to $1$ than the baseline $0.83$. We cannot,
however, exclude that a given injection could be best matched by a template
corresponding to a different, nearby injection, leading to $\xi\neq 1$.
Consequently, we perform a separate test where we explicitly select the PyCBC
Live candidates that exactly match the masses and spins of the injected
signal.

About 500 injections are recovered by PyCBC Live with this bank, which is
similar to the tests presented in the later sections with more conventional
banks. This excludes a pathological behaviour of the search with this bank.
350 triggers have an exact match with the injected signal. Figure 8 presents
the optimisation curves obtained with this bank. The dashed-dotted brown line
corresponds to the results for all the triggers. The one in purple corresponds
to the results for the subset of triggers with an exact match between the
injected and recovered templates. For all the triggers mixed, the optimal
value of $\xi$ is $0.90$ and for the exact match, $\xi$ is optimal at $0.93$.
As expected, both values are indeed closer to $\xi=1$ compared to Section VI,
and the exact-match result has the larger $\xi$. Note, however, that the
optimal $\xi$ in the case of exactly-matching intrinsic parameters is still
noticeably below $1$, which appears different from our previous results based
on pycbc_make_skymap. This effect might originate from the differences between
the SpinTaylorT4 injections and the TaylorF2 templates used by PyCBC Live.

Figure 8: Optimisation of the $\xi$ parameter in Bayestar for PyCBC Live. Each
curve corresponds to a set of triggers produced by a PyCBC Live analysis with
a different template bank. The curves correspond to results presented in
different subsections of Section VII, as detailed in the legend. The dashed-
dotted curves correspond to the results obtained with the template banks based
on the injection set described in Section VII.2. The purple curve corresponds
to the results obtained with all the triggers produced by PyCBC Live and the
brown one to the subset of triggers exactly matching the injection they
recovered. The blue curve corresponds to the template bank containing only BNS
signals without spin and a minimal match of 97%. The red curve corresponds to
a BNS bank including nonzero spins and minimal match of 97%. The green curve
corresponds to the O3 bank. The vertical orange-dashed line corresponds to the
baseline $\xi=0.83$ value implemented in Bayestar.

Together with the results of Section VI, these results suggest that the
template bank does play a role in the underestimation of the skymap
uncertainty.

###  VII.3 Non-spinning banks of varying density

We now extend the search space of our bank to cover the entire range of
component masses of our simulated signals, but set the component spins to zero
for all templates. We also test the hypothesis that using a sparser template
bank may lead to an increase in the discrepancy between the true waveform
parameters and the ones of the template reported by the search, thus driving
the optimal $\xi$ to lower values than found in Section VII.2.

Template banks are designed to have the smallest number of templates for
maintaining computational efficiency without missing signals because of the
mismatch between the model and the true signal. For a given vector of
parameters $\vec{\theta}$, the match between a template $h(\vec{\theta})$ and
a signal $h(\vec{\theta}+\vec{d\theta})$ is defined as [37]:

$M=\mathrm{max}\left(\frac{\langle
h(\vec{\theta})|h(\vec{\theta}+d\vec{\theta})\rangle}{\sqrt{\langle
h(\vec{\theta})|h(\vec{\theta})\rangle\langle
h(\vec{\theta}+d\vec{\theta})|h(\vec{\theta}+d\vec{\theta})\rangle}}\right)$
(5)

where the maximization is carried out over time and phase shifts. A usual
criterion to set the spacing between templates is to introduce a parameter
called _minimal match_ ($MM$) such that $M\geq MM$ for any signal in the
search space. $MM$ fixes the maximum acceptable fractional SNR loss incurred
because of the sparseness of the bank. A commonly adopted value is $MM=97$%,
and this is also the value used for part of the template bank used by PyCBC
Live in O3.

We consider here three banks to test the impact of the sparseness on
localization. These banks are composed of low-mass, non-spinning CBC
waveforms: a sparse bank with $MM=90$%, one with the common $MM=97$% and a
finer bank with $MM=98$%. These banks are generated using a geometrical
placement described in [38], where the templates are placed in the parameter
space on a hexagonal lattice. The generation is made with a PyCBC tool called
pycbc_geom_nonspinbank. The templates of these banks are post-Newtonian
waveforms at the $3.5$ order and have both component masses in the $[1,12]$M⊙
range to safely cover our BNS injections. We then apply a cut on the chirp
mass at $2.7$M⊙, so as to reduce the number of templates as much as possible,
while still covering the BNS injections and avoiding border effects
potentially arising from a cut on component masses.

Figure 9 shows the results for the three banks, and the blue curve of Figure 8
shows the results for $MM=97$% to compare to other template banks of this
Section. Remarkably, the optimal $\xi$ value _decreases_ when the $MM$
increases, although this effect is rather small, as several values of $\xi$
produce deviations below $1\sigma$, indicating acceptable performances. We
also observe that the optimal $\xi$ is further reduced with respect to the
case where the exact injections are used as a template bank. We conclude that
the density of templates appears to have a small but detectable effect on the
optimal choice of $\xi$, as suggested in [13]. The effect goes into a somewhat
counterintuitive direction, as denser banks appear to require a stronger
correction. It also appears that the amount of explored parameter space around
each signal may also have a small effect on the optimal $\xi$.

Figure 9: Optimisation curves for low-mass, zero-spin template banks with
minimal matches of 90% (green-dotted), 97% (blue-solid) and 98% (red-dashed).

###  VII.4 Adding aligned spins to the search space

As the next step, we include spins in the search parameter space. The mass
ranges are the same as in the previous section, but we also include aligned
spins with magnitudes between $0$ and $0.05$ for neutron stars and between $0$
and $0.998$ for black holes. We use a minimal match of 97%, for this is the
usual value for the actual search banks. The curve in red in Figure 8 shows
the results for localising the triggers found with this bank.

The optimal value of $\xi$ is $0.87$, slightly smaller than the zero-spin bank
with the same $MM$ as described in the previous section. This result goes in
the direction consistent with the hypothesis that $\xi$ is somehow related to
the extent (or dimensionality) of the search space around each signal.

###  VII.5 O3 bank

As a last check, we use the same bank used by PyCBC Live during the O3 run on
the simulated data. The bank design is detailed in [5, 39], and the parameter
space explored by this bank includes BNS, BBH and NSBH. For BNS and NSBH
binaries, the neutron stars have masses between 1 and 2 $M_{\odot}$ and
aligned spins with magnitude up to $0.05$. These boundaries account for the
masses and largest spins observed in pulsars [40]. For black holes, the spin
magnitudes are up to $0.998$, accounting for the spins observed in X-ray
binaries [41]. The total mass of the binaries goes up to $100M_{\odot}$ for
equal mass, weakly spinning templates and up to $500M_{\odot}$ for templates
with high mass ratio, and high aligned spins.

The results for the O3 bank are visible as the green curve in Figure 8 and
show an optimal $\xi$ of $0.87$, slightly higher than the baseline $0.83$
value used by default in Bayestar. This optimal value is compatible with the
previous bank with $MM=97$% and nonzero spins. The difference with the default
$\xi$ suggests that the sky localizations for PyCBC Live BNS events are, on
average, slightly broader than they should be. As such, $\xi$ should be
slightly increased when localising BNS candidates produced by PyCBC Live,
particularly considering that the deviation at the baseline $\xi$ approaches
the $3\sigma$ level.

##  VIII Interpretation

The various tests carried out in the previous sections indicate that the
optimal value of $\xi$ is affected in varying amounts by a variety of possible
effects. One first conclusion we can make is that real detector noise, which
is non-stationary and non-Gaussian, does not appear to play a significant role
in terms of $\xi$. This is somewhat consistent with the findings in [35],
especially for BNS mergers, though that study focused on very specific cases.
On the other hand, an important driver of the optimal $\xi$ appears to be
related to the inevitable discrepancy between the true intrinsic parameters of
a CBC signal, and those of the templates identified by the matched filtering
search that produced the candidate event. This was in fact suggested as a
possibility in [13]. In particular, if we pretend that the search will exactly
identify the true intrinsic parameters, then the $\xi$ correction is not
necessary. However, as soon as we allow a deviation, we immediately observe
that a value $\xi\sim 0.87$ becomes necessary to properly compute the
uncertainty in the sky localization. The optimal value depends on how exactly
the deviation takes place: if we allow only a “minimal” deviation,
corresponding to a template bank that only explores the set of simulated
signals, then the optimal $\xi$ is closer to $1$. As we introduce more and
more templates in the vicinity of the true values of each signal (either by
increasing the minimal match of the bank, or by including additional
dimensions in the search space, such as spins) the optimal value tends to
further decrease, albeit by a smaller amount. Perhaps counter-intuitively, we
find that making the bank denser leads to a slightly _smaller_ optimal $\xi$.
Therefore, it appears that the bias is originating from the “freedom” that the
search template has to slightly deviate from the true parameters.

Based on these observations, we propose the following argument. Bayestar makes
the assumption that the uncertainties in the spatial location parameters are
uncorrelated with those in the intrinsic parameters. This assumption is
certainly appropriate, and it introduces no systematic shift of the posterior
distribution of the spatial parameters. However, we argue that it is
implicitly acting as a delta-prior on the intrinsic parameters, effectively
assuming that the search template is exactly matching the true intrinsic
parameters. If this was true, then the spread of Bayestar’s posterior
distribution for the spatial position would match what it would be if Bayestar
carried out a full parameter estimation of the entire vector of unknown signal
parameters, followed by a marginalization over the intrinsic parameters.
However, the matched-filter search returns effectively one sample from the
peak of the likelihood function over the intrinsic parameters, and this sample
will in general not be at the true values of the parameters, and for technical
reasons may not even be at the maximum of the likelihood (see e.g. the
discussion in [5] for PyCBC Live). Therefore, the spread of Bayestar’s
posterior distribution for the spatial position will be slightly
underestimated with respect to a full parameter estimation followed by
marginalization. The $\xi$ factor is then used to compensate this
underestimation by introducing an additional spread in the spatial posterior
distribution. Its optimal value is therefore related to how likely the single
sample of the intrinsic parameters is to be at a certain distance from the
true intrinsic parameters, and this ultimately depends on how many templates
are “available” around the true parameters.

To conclude, we argue that $\xi$ has to be set in a way that depends on the
specific properties of the template bank used by the search, namely its
density (or minimal match), its dimensionality, and the extent of the search
space explored by it. With more work it might be possible, for instance, to
relate $\xi$ to the average match between the true waveforms and the templates
reported by the search. We consider our simulations as an initial step towards
understanding this dependence.

##  IX Conclusion

In this work, we consider GW candidate events produced by the PyCBC Live
detection pipeline and present an analysis for optimizing the low-latency sky
localization of those events with the Bayestar algorithm. The latter provides
skymaps within seconds after detecting a GW signal with online searches. This
information is crucial for multi-messenger astronomy and electromagnetic
follow-up, making the low-latency localization a critical part of the follow-
up. We evaluate the accuracy of Bayestar by simulating a population of CBC
signals, recovering and localizing them under a variety of assumptions, and
checking the consistency of the resulting sky localization using a common test
based on the PP plot, which is expected to be diagonal for correctly estimated
uncertainties.

Bayestar has a tuning parameter $\xi$ that allows one to make the sky
localization self-consistent. It acts as a scaling factor for the SNR reported
by the search pipeline before the Bayesian inference. Originally this factor
had been tuned based on a specific set of assumptions for passing the PP-plot
test, resulting in a default value of $0.83$ [13]. We first show that an
idealized condition exists where such a correction is not necessary
(equivalent to setting $\xi=1$), namely when the true intrinsic parameters of
a signal exactly match those reported by the search pipeline. We also show
that the correction ($\xi<1$) becomes necessary as soon as there is even a
minor mismatch between the true intrinsic parameters of a signal and those
reported by the search pipeline, which of course is always the case in
practice.

Focusing on BNS mergers, which are the most promising sources of EM
counterparts, we use a series of end-to-end tests with PyCBC Live in order to
identify the search parameters influencing $\xi$. Starting with a template
bank for searching for GW signal that closely matches the simulated signal, we
observe that the optimal value for $\xi$ is closer to $1$ compared to the
baseline value. We also use three template banks with only BNS signals but
different minimal matches, a parameter controlling the bank sparseness. Based
on our results, we conclude that the optimal $\xi$ depends to a small extent
both on the sparseness of the template bank and on its dimensionality,
consistently with what was suggested in [13]. Finally, we test the
localization with the template bank used for the O3 analysis, and we find that
the optimal value is $\xi\sim 0.87$, i.e. closer to the default than our
previous tests, but slightly larger. This suggests that a dedicated tuning
campaign of $\xi$ for PyCBC Live, and more generally for any search using
template banks significantly different than what was used in [13] for GstLAL,
might improve the consistency of Bayestar’s sky localizations. We also find
that the consistency of the localization for BNS systems does not appear to
depend significantly on the quality of the data, i.e. on differences between
real detector noise and stationary Gaussian noise. This is somewhat consistent
with the findings of [35].

Tests similar to ours have been carried out recently in [42], using a
dedicated end-to-end mock data challenge in low-latency, but without exploring
the $\xi$ parameter in great detail. The resulting PP plots appear to indicate
consistent sky localizations for all search pipelines currently operating in
the O4 run. This is both reassuring, and perhaps not surprising considering
that the setup of that simulation, and in particular the properties of the
simulated signals, are quite different from our choices here.

In this work, we have only explored the consistency of the posterior
distribution for the two-dimensional position on the sky, focusing on the main
needs of optical telescopes. Similar investigations should be carried out in
the future for the luminosity distance posterior, also provided by Bayestar,
or more generally for the joint three-dimensional spatial distribution.
Another avenue worth investigating is the effect of the post-discovery SNR
maximization done in O3 by PyCBC Live, described in [5]. Since this method
should lead to an improved match between the true source parameters and the
selected template, it should in principle lead to an optimal $\xi$ closer to
$1$. Consequently, a simulation campaign testing the influence of this method
should be made in the future.

###### Acknowledgements.

We thank Leo Singer for useful discussion and help with the Bayestar code, and
Thomas Dent, Michael Coughlin, Ed Porter, Colm Talbot and Jolien Creighton for
useful comments on the study and draft. TD additionally thanks FP for the kind
hospitality at the Sapienza University of Rome during part of this work. SA
thanks the CSI Recherche University Côte d’Azur for their financial support.
Part of our simulations utilized the Virtual Data cloud computing system at
IJCLab. We thank Michel Jouvin and Gerard Marchal-Duval for their prompt
support and advice about this system. The authors are also grateful for
computational resources provided by the LIGO Laboratory and supported by
National Science Foundation Grants PHY-0757058 and PHY-0823459. This
manuscript has LIGO document number LIGO-P2300429.

## References

  * LIGO Scientific Collaboration  _et al._ [2015] LIGO Scientific Collaboration, J. Aasi, B. P. Abbott, R. Abbott, T. Abbott, M. R. Abernathy, K. Ackley, _et al._ , Classical and Quantum Gravity 32, 074001 (2015), arXiv:1411.4547 [gr-qc] .
  * Acernese  _et al._ [2015] F. Acernese, M. Agathos, K. Agatsuma, D. Aisa, N. Allemandou, A. Allocca, J. Amarni, P. Astone, G. Balestri, _et al._ , Classical and Quantum Gravity 32, 024001 (2015), arXiv:1408.3978 [gr-qc] .
  * Akutsu  _et al._ [2021] T. Akutsu, M. Ando, K. Arai, Y. Arai, S. Araki, A. Araya, N. Aritomi, Y. Aso, S. Bae, Y. Bae, L. Baiotti, R. Bajpai, _et al._ , Progress of Theoretical and Experimental Physics 2021, 05A101 (2021), arXiv:2005.05574 [physics.ins-det] .
  * Nitz  _et al._ [2018] A. H. Nitz, T. Dal Canton, D. Davis, and S. Reyes, Phys. Rev. D 98, 024050 (2018), arXiv:1805.11174 [gr-qc] .
  * Dal Canton  _et al._ [2021] T. Dal Canton, A. H. Nitz, B. Gadre, G. S. Cabourn Davies, V. Villa-Ortega, T. Dent, I. Harry, and L. Xiao, APJ 923, 254 (2021), arXiv:2008.07494 [astro-ph.HE] .
  * Cannon  _et al._ [2021] K. Cannon, S. Caudill, C. Chan, B. Cousins, J. D. E. Creighton, B. Ewing, H. Fong, P. Godwin, _et al._ , SoftwareX 14, 100680 (2021), arXiv:2010.05082 [astro-ph.IM] .
  * Adams  _et al._ [2016] T. Adams, D. Buskulic, V. Germain, G. M. Guidi, F. Marion, M. Montani, B. Mours, F. Piergiovanni, and G. Wang, Classical and Quantum Gravity 33, 175012 (2016), arXiv:1512.02864 [gr-qc] .
  * Chu  _et al._ [2020] Q. Chu, M. Kovalam, L. Wen, T. Slaven-Blair, J. Bosveld, Y. Chen, P. Clearwater, A. Codoreanu, Z. Du, X. Guo, X. Guo, K. Kim, T. G. F. Li, V. Oloworaran, F. Panther, J. Powell, A. S. Sengupta, K. Wette, and X. Zhu, arXiv e-prints , arXiv:2011.06787 (2020), arXiv:2011.06787 [gr-qc] .
  * Klimenko  _et al._ [2016] S. Klimenko, G. Vedovato, M. Drago, F. Salemi, V. Tiwari, G. A. Prodi, C. Lazzaro, K. Ackley, S. Tiwari, C. F. Da Silva, and G. Mitselmakher, Phys. Rev. D 93, 042004 (2016), arXiv:1511.05999 [gr-qc] .
  * Abbott  _et al._ [2023] R. Abbott, T. D. Abbott, F. Acernese, K. Ackley, C. Adams, N. Adhikari, _et al._ (LIGO Scientific Collaboration, Virgo Collaboration, and KAGRA Collaboration), Phys. Rev. X 13, 041039 (2023).
  * Abbott  _et al._ [2021] R. Abbott, T. D. Abbott, S. Abraham, F. Acernese, K. Ackley, A. Adams, C. Adams, R. X. Adhikari, V. B. Adya, C. Affeldt, M. Agathos, K. Agatsuma, N. Aggarwal, O. D. Aguiar, L. Aiello, LIGO Scientific Collaboration, Virgo Collaboration, _et al._ , Physical Review X 11, 021053 (2021), arXiv:2010.14527 [gr-qc] .
  * [12] NASA, General coordinates network.
  * Singer and Price [2016] L. P. Singer and L. R. Price, Phys. Rev. D 93, 024013 (2016), arXiv:1508.03634 [gr-qc] .
  * Morisaki  _et al._ [2023] S. Morisaki, R. Smith, L. Tsukada, S. Sachdev, S. Stevenson, C. Talbot, and A. Zimmerman, Rapid localization and inference on compact binary coalescences with the Advanced LIGO-Virgo-KAGRA gravitational-wave detector network (2023), arXiv:2307.13380 [gr-qc] .
  * Ashton  _et al._ [2019] G. Ashton, M. Hübner, P. D. Lasky, C. Talbot, K. Ackley, S. Biscoveanu, Q. Chu, A. Divakarla, P. J. Easter, B. Goncharov, F. Hernandez Vivanco, J. Harms, M. E. Lower, G. D. Meadors, D. Melchor, E. Payne, M. D. Pitkin, J. Powell, N. Sarin, R. J. E. Smith, and E. Thrane, APJS 241, 27 (2019), arXiv:1811.02042 [astro-ph.IM] .
  * Abbott  _et al._ [2017a] B. P. Abbott, R. Abbott, T. D. Abbott, F. Acernese, K. Ackley, C. Adams, T. Adams, P. Addesso, R. X. Adhikari, V. B. Adya, C. Affeldt, M. Afrough, B. Agarwal, M. Agathos, K. Agatsuma, (INTEGRAL, _et al._ , APJL 848, L13 (2017a), arXiv:1710.05834 [astro-ph.HE] .
  * Goldstein  _et al._ [2017] A. Goldstein, P. Veres, E. Burns, M. S. Briggs, R. Hamburg, D. Kocevski, C. A. Wilson-Hodge, R. D. Preece, S. Poolakkil, O. J. Roberts, C. M. Hui, V. Connaughton, J. Racusin, A. von Kienlin, T. Dal Canton, N. Christensen, T. Littenberg, K. Siellez, L. Blackburn, J. Broida, E. Bissaldi, W. H. Cleveland, M. H. Gibby, M. M. Giles, R. M. Kippen, S. McBreen, J. McEnery, C. A. Meegan, W. S. Paciesas, and M. Stanbro, APJL 848, L14 (2017), arXiv:1710.05446 [astro-ph.HE] .
  * LIGO Scientific Collaboration  _et al._ [2017] LIGO Scientific Collaboration, Virgo Collaboration, IceCube Collaboration, AstroSat Team, Insight-HXMT Collaboration, ANTARES Collaboration, Swift Collaboration, AGILE Team, DES Collaboration, GRAWITA Telescope, ATCA Telescope, AGILE Team, ASKAP Group, _et al._ , APJL 848, L12 (2017), arXiv:1710.05833 [astro-ph.HE] .
  * D’Avanzo  _et al._ [2018] P. D’Avanzo, S. Campana, O. S. Salafia, G. Ghirlanda, G. Ghisellini, A. Melandri, M. G. Bernardini, M. Branchesi, E. Chassande-Mottin, S. Covino, V. D’Elia, L. Nava, R. Salvaterra, G. Tagliaferri, and S. D. Vergani, Astronomy and Astrophysics 613, L1 (2018), arXiv:1801.06164 [astro-ph.HE] .
  * Alexander  _et al._ [2017] K. D. Alexander, E. Berger, W. Fong, P. K. G. Williams, C. Guidorzi, R. Margutti, B. D. Metzger, J. Annis, P. K. Blanchard, D. Brout, D. A. Brown, H. Y. Chen, R. Chornock, P. S. Cowperthwaite, M. Drout, T. Eftekhari, J. Frieman, D. E. Holz, M. Nicholl, A. Rest, M. Sako, M. Soares-Santos, and V. A. Villar, APJL 848, L21 (2017), arXiv:1710.05457 [astro-ph.HE] .
  * Hallinan  _et al._ [2017] G. Hallinan, A. Corsi, K. P. Mooley, K. Hotokezaka, E. Nakar, M. M. Kasliwal, D. L. Kaplan, D. A. Frail, _et al._ , Science 358, 1579 (2017), arXiv:1710.05435 [astro-ph.HE] .
  * Troja  _et al._ [2017] E. Troja, L. Piro, H. van Eerten, R. T. Wollaeger, M. Im, O. D. Fox, N. R. Butler, S. B. Cenko, _et al._ , Nature (London) 551, 71 (2017), arXiv:1710.05433 [astro-ph.HE] .
  * Soares-Santos  _et al._ [2017] M. Soares-Santos, D. E. Holz, J. Annis, R. Chornock, K. Herner, E. Berger, D. Brout, H. Y. Chen, _et al._ , APJL 848, L16 (2017), arXiv:1710.05459 [astro-ph.HE] .
  * Petrov  _et al._ [2022] P. Petrov, L. P. Singer, M. W. Coughlin, V. Kumar, M. Almualla, S. Anand, M. Bulla, T. Dietrich, F. Foucart, and N. Guessoum, APJ 924, 54 (2022), arXiv:2108.07277 [astro-ph.HE] .
  * Kiendrebeogo  _et al._ [2023] R. W. Kiendrebeogo, A. M. Farah, E. M. Foley, A. Gray, N. Kunert, A. Puecher, A. Toivonen, R. O. VandenBerg, _et al._ , APJ 958, 158 (2023), arXiv:2306.09234 [astro-ph.HE] .
  * Allen  _et al._ [2012] B. Allen, W. G. Anderson, P. R. Brady, D. A. Brown, and J. D. E. Creighton, Phys. Rev. D 85, 122006 (2012), arXiv:gr-qc/0509116 [gr-qc] .
  * Welch [1967] P. Welch, IEEE Transactions on Audio and Electroacoustics 15, 70 (1967).
  * [28] LIGO Scientific Collaboration, Virgo Collaboration, and KAGRA Collaboration, GraceDB.
  * Sidery  _et al._ [2014] T. Sidery, B. Aylott, N. Christensen, B. Farr, W. Farr, F. Feroz, J. Gair, K. Grover, _et al._ , Phys. Rev. D 89, 084060 (2014), arXiv:1312.6013 [astro-ph.IM] .
  * Buonanno  _et al._ [2009] A. Buonanno, B. R. Iyer, E. Ochsner, Y. Pan, and B. S. Sathyaprakash, Phys. Rev. D 80, 084043 (2009), arXiv:0907.0700 [gr-qc] .
  * Bohé  _et al._ [2017] A. Bohé, L. Shao, A. Taracchini, A. Buonanno, S. Babak, I. W. Harry, I. Hinder, S. Ossokine, M. Pürrer, V. Raymond, T. Chu, H. Fong, P. Kumar, H. P. Pfeiffer, M. Boyle, D. A. Hemberger, L. E. Kidder, G. Lovelace, M. A. Scheel, and B. Szilágyi, Phys. Rev. D 95, 044028 (2017), arXiv:1611.03703 [gr-qc] .
  * LIGO Scientific Collaboration [2018] LIGO Scientific Collaboration, LIGO Algorithm Library - LALSuite, free software (GPL) (2018).
  * Dodge [2008] Y. Dodge, Kolmogorov–smirnov test, in  _The Concise Encyclopedia of Statistics_ (Springer New York, NY, 2008) pp. 283–287.
  * Abbott  _et al._ [2017b] B. P. Abbott, R. Abbott, T. D. Abbott, F. Acernese, K. Ackley, C. Adams, T. Adams, P. Addesso, R. X. Adhikari, V. B. Adya, C. Affeldt, _et al._ , Phys. Rev. Lett. 119, 161101 (2017b), arXiv:1710.05832 [gr-qc] .
  * Macas  _et al._ [2022] R. Macas, J. Pooley, L. K. Nuttall, D. Davis, M. J. Dyer, Y. Lecoeuche, J. D. Lyman, J. McIver, and K. Rink, Phys. Rev. D 105, 103021 (2022), arXiv:2202.00344 [astro-ph.HE] .
  * Biscoveanu  _et al._ [2019] S. Biscoveanu, S. Vitale, and C.-J. Haster, APJL 884, L32 (2019), arXiv:1908.03592 [astro-ph.HE] .
  * Owen [1996] B. J. Owen, Phys. Rev. D 53, 6749 (1996), arXiv:gr-qc/9511032 [gr-qc] .
  * Cokelaer [2007] T. Cokelaer, Phys. Rev. D 76, 102004 (2007), arXiv:0706.4437 [gr-qc] .
  * Dal Canton and Harry [2017] T. Dal Canton and I. W. Harry, arXiv e-prints , arXiv:1705.01845 (2017).
  * Lorimer [2008] D. R. Lorimer, Living Reviews in Relativity 11, 8 (2008), arXiv:0811.0762 [astro-ph] .
  * McClintock  _et al._ [2014] J. E. McClintock, R. Narayan, and J. F. Steiner, Space Science Reviews 183, 295 (2014), arXiv:1303.1583 [astro-ph.HE] .
  * Sharma Chaudhary  _et al._ [2023] S. Sharma Chaudhary, A. Toivonen, G. Waratkar, G. Mo, D. Chatterjee, S. Antier, P. Brockill, _et al._ , arXiv e-prints , arXiv:2308.04545 (2023), arXiv:2308.04545 [astro-ph.HE] .