# The unresolved stochastic background from compact binary mergers detectable
by next-generation ground-based gravitational-wave observatories

Darsan S. Bellie<EMAIL_ADDRESS>Center for Interdisciplinary
Exploration and Research in Astrophysics (CIERA), Northwestern University,
1800 Sherman Ave, Evanston, IL 60201, USA  Department of Physics and
Astronomy, Northwestern University, 2145 Sheridan Road, Evanston, IL 60201,
USA  Sharan Banagiri  Center for Interdisciplinary Exploration and Research in
Astrophysics (CIERA), Northwestern University, 1800 Sherman Ave, Evanston, IL
60201, USA  Zoheyr Doctor  Center for Interdisciplinary Exploration and
Research in Astrophysics (CIERA), Northwestern University, 1800 Sherman Ave,
Evanston, IL 60201, USA  Vicky Kalogera  Center for Interdisciplinary
Exploration and Research in Astrophysics (CIERA), Northwestern University,
1800 Sherman Ave, Evanston, IL 60201, USA  Department of Physics and
Astronomy, Northwestern University, 2145 Sheridan Road, Evanston, IL 60201,
USA

###### Abstract

The next generation of ground-based gravitational-wave detectors will look
much deeper into the Universe and have unprecedented sensitivities and low-
frequency capabilities. Especially alluring is the possibility of detecting an
early-Universe cosmological stochastic background that could provide important
insights into the beginnings of our Universe and fundamental physics at
extremely high energies. However, even if next-generation detectors are
sensitive to cosmological stochastic backgrounds, they will be masked by more
dominant astrophysical backgrounds, namely the residual background from the
imperfect subtraction of resolvable compact binary coalescences (CBCs) as well
as the CBC background from individually unresolvable CBCs. Using our latest
knowledge of masses, rates, and delay time distributions, we present a data-
driven estimate of the unresolvable CBC background that will be seen by next-
generation detectors. Accounting for statistical and systematic errors, this
estimate quantifies an important piece in the CBC noise budget for next-
generation detectors and can help inform detector design and subtraction
algorithms. We compare our results with predictions for backgrounds from
several cosmological sources in the literature, finding that the unresolvable
background will likely be a significant impediment for many models. This
motivates the need for simultaneous inference methods or other statistical
techniques to detect early-Universe cosmological backgrounds.

##  I Introduction

In 2015, the Laser Interferometer Gravitational-Wave Observatory (LIGO) Aasi
_et al._ (2015) directly detected gravitational waves for the first time by
capturing the coalescence of two stellar-mass binary black holes (BBHs) Abbott
_et al._ (2016). This was eventually followed by the detection of GWs from
binary neutron stars (BNSs) Abbott  _et al._ (2017a) and neutron-star black-
hole binaries (NSBHs) Abbott  _et al._ (2021a) by the LIGO and Virgo
observatories Acernese  _et al._ (2015). Constituting the category of compact
binary coalescences (CBCs), these sources have heralded the arrival of GW
astronomy with nearly a hundred detected candidates already Abbott  _et al._
(2019a, 2021b, 2021c).

While we have only detected GWs from CBCs so far, several other sources have
been theorized and expected. In particular, several models postulate the
presence of “cosmological” GWs, predominantly showing up as a cosmological
stochastic gravitational wave background (SGWB) Caprini and Figueroa (2018);
Christensen (2019); Sakellariadou (2022). These might be accessible by next-
generation ground-based GW detectors (henceforth known as XG detectors) like
Cosmic Explorer (CE) and Einstein Telescope (ET) Evans  _et al._ (2021);
Punturo  _et al._ (2010); Regimbau  _et al._ (2017); Maggiore  _et al._
(2020). Cosmological SGWBs originating from the early Universe might carry
imprints of the physics of the earliest epochs in the Universe, inaccessible
through any other channels. Detection of such a background would therefore
drastically expand our understanding of the first instants in the Universe and
unlock fundamental physics at very high energies well beyond the reach of
particle accelerators Caprini and Figueroa (2018).

However, such cosmological SGWBs can be shielded by the GWs arising from
various astrophysical sources Regimbau (2011); Regimbau and Mandic (2008);
Buonanno  _et al._ (2005) and in particular CBCs of stellar origin since they
are expected to be the largest contributor to the GW power in the frequency
band of XG detectors. The GWs from CBCs will pose a major challenge for the
capability of XG detectors to detect cosmological backgrounds. There are two
components of such astrophysical shielding, coming from “resolved” and
“unresolved” CBCs. Resolved sources, here, refer to those sources with GW
signals that can be individually and confidently detected within the detector
noise.

The proposed XG detectors are expected to have a high redshift reach to
astrophysical CBC sources owing to their better sensitivities and low-
frequency capabilities of up to $\sim$ 5 Hz Kalogera  _et al._ (2019); Evans
_et al._ (2023, 2021); Kalogera  _et al._ (2019). It is expected that XG
detectors would be able to resolve nearly all BBHs of stellar origin in the
Universe as well as many BNSs and NSBHs up to redshifts of a few Kalogera  _et
al._ (2019); Evans  _et al._ (2021, 2023). To minimize their impact on
searches for cosmological SGWBs, significant ongoing work is devoted to
developing techniques to subtract resolvable CBCs from GW data Zhou  _et al._
(2023, 2022); Zhong  _et al._ (2023); Sachdev  _et al._ (2020); Regimbau
(2022); Sharma and Harms (2020). However, the subtraction residue from
imperfect subtraction could still be a challenge for XG detectors Zhou  _et
al._ (2022, 2023).

In addition to the resolvable CBC signals, we also expect an astrophysical
SGWB originating from the incoherent superposition of unresolvable CBC
sources. Such an astrophysical SGWB is a key detection target for current-
generation ground-based detectors Abbott  _et al._ (2017b, 2019b, 2021d), but
it can also shield cosmological SGWBs. Crucially, this unresolved CBC
background represents a noise source that is independent of the fidelity and
efficacy of subtraction. Accurately characterizing this background will enable
us to understand the accessibility of cosmological SGWBs to XG detectors, will
inform us of the necessary levels of subtraction of resolvable signals, and
will provide an important data point in comparing detector designs and
networks. Therefore, in this paper, we provide a data-driven estimate of the
levels of the unresolvable CBC SGWB expected for different configurations of
XG detector networks. We draw upon several models used in the population
inference with LIGO-Virgo-KAGRA (LVK) observations and Galactic binary neutron
stars; by utilizing multiple models, we provide robust estimates that account
for both statistical and systematic uncertainty.

The rest of the paper is organized as follows. In Sec. II we briefly motivate
the paper and place our study in the context of previous literature for XG
detectors. In Sec. III, we mathematically define the SGWB and describe how an
astrophysical SGWB can be calculated from a set of sources. Sec. IV describes
our assumptions on the astrophysical merger rates, mass distributions, and
redshift distributions which we use to compute the unresolved background, and
Sec. V describes our choices of detector network configurations. In Sec. VI,
we detail our formalism to numerically compute the unresolved background from
CBC sources and in Sec. VII, we discuss the cosmological SGWB models that we
use. Sec. VIII presents our results, and Sec. IX discusses their implications
and comments on future work needed.

Throughout this paper, $G$ is the gravitational constant, $c$ the speed of
light, $H_{0}$ the Hubble constant, and we define our cosmology using the
$\Lambda$CDM model with cosmological parameters taken from Planck 2018 Aghanim
_et al._ (2020).

##  II The unresolved SGWB from CBCs

The total GW energy-density spectrum $\Omega_{\text{GW}}$ as seen by a GW
detector is the sum of the GW energy-density spectrum from cosmological
sources $\Omega_{\text{cosmo}}$ and from astrophysical sources 111i.e. of
stellar origin. $\Omega_{\text{astro}}$ Sachdev  _et al._ (2020),

$\Omega_{\text{GW}}=\Omega_{\text{astro}}+\Omega_{\text{cosmo}}.$ (1)

In general, the astrophysical component can mask the cosmological component if
$\Omega_{\text{astro}}>\Omega_{\text{cosmo}}$. For the rest of this paper, we
assume that $\Omega_{\text{astro}}=\Omega_{\text{cbc}}$ and that the GW energy
from other sources can be neglected so that

$\Omega_{\text{GW}}\approx\Omega_{\text{cbc}}+\Omega_{\text{cosmo}}.$ (2)

Various methods have been developed to subtract resolvable CBC signals in
order to minimize their impact on searches for $\Omega_{\text{cosmo}}$ Zhou
_et al._ (2023, 2022); Zhong  _et al._ (2023); Sachdev  _et al._ (2020);
Regimbau (2022); Sharma and Harms (2020). Such techniques, however, are always
imperfect, as they depend on parameter estimation results and can also suffer
from waveform systematics. The imperfect subtraction leaves behind a residue
$\Omega_{\text{cbc, residue}}$ that contributes to the effective SGWB from
CBCs Sachdev  _et al._ (2020); Zhou  _et al._ (2023).

In addition to the resolved CBC signals, the incoherent superposition of the
unresolvable CBC signals gives rise to an astrophysical SGWB that we call
$\Omega_{\text{cbc, unres}}$. Therefore, effectively,

$\Omega_{\text{cbc}}=\Omega_{\text{cbc, unres}}+\Omega_{\text{cbc, residue}}.$
(3)

While some studies find that $\Omega_{\text{cbc, residue}}$ remains a major
challenge for detecting cosmological SGWBs Zhou  _et al._ (2023, 2022), others
find that $\Omega_{\text{cbc, unres}}$ will provide the noise floor Sharma and
Harms (2020); Zhong  _et al._ (2023). Since $\Omega_{\text{cbc, unres}}$ is a
purely stochastic noise, it will have to be fit simultaneously with any
cosmological SGWB.

While some constraints were previously placed on $\Omega_{\text{cbc, unres}}$
for XG detectors by Refs. Zhou  _et al._ (2023); Sachdev  _et al._ (2020);
Regimbau  _et al._ (2017); Wu  _et al._ (2012); Zhong  _et al._ (2023); Zhu
_et al._ (2013) they were either limited by the dearth of data before O3, or
by simple population models or both. Ref.Zhou  _et al._ (2022) in particular
accounts for the uncertainties in the latest LVK BBH and BNS rate inferences
in their calculation of the residue but not for $\Omega_{\text{cbc, unres}}$.
Furthermore, papers in the literature generally do not include the
contribution from NSBHs, whose unresolved SGWB, we find, is comparable with
that of BNSs 222Although Ref. Zhu  _et al._ (2013) includes NSBHs in their
estimate, this was not data-informed..

In our study, we implement a data-driven estimate of $\Omega_{\text{cbc,
unres}}$ for XG detectors, using the latest inferences on the CBC population
from LVK data. Wherever possible, we consistently incorporate the uncertainty
from the rate and from the population modeling. In addition, we consider
several population models to propagate the systematic uncertainty stemming
from parameterized GW analyses to present uncertainty envelopes for
$\Omega_{\text{cbc, unres}}$.

##  III Calculating the SGWB from CBCs

The SGWB from any type of GW source is generally characterized by the
dimensionless GW energy-density spectrum Romano and Cornish (2017)

$\Omega_{\text{GW}}(f)=\frac{f}{\rho_{c}}\frac{d\rho_{\text{GW}}}{df}(f),$ (4)

where $f$ is the observed GW frequency, $\rho_{c}=3H_{0}^{2}c^{2}/8\pi G$ is
the critical density needed to close the Universe, and $\rho_{\text{GW}}$ is
the GW energy density at $f$.

The GW energy density $\rho_{\text{GW}}$ arising from astrophysical sources is
given by

$\frac{d\rho_{\text{GW}}}{df}(f)=\frac{1}{c}F(f),$ (5)

where $F$ is the GW energy flux in the observer frame Ferrari  _et al._
(1999). We define the energy flux by summing up the individual fluxes in a
population of astrophysical sources between redshifts $z_{\rm low}$ to $z_{\rm
up}$ and with an associated set of source parameters $\theta$ as Regimbau
(2011, 2022)

$F(f)=\int_{\theta}p(\theta)d\theta\int_{z_{\text{low}}}^{z_{\text{up}}}\frac{dR_{o}(z)}{dz}\frac{\frac{dE_{\text{gw}}(f_{\text{s}},\theta)}{df_{\text{s}}}(1+z)^{2}}{4\pi
d_{L}^{2}}dz,$ (6)

where $p(\theta)$ is the distribution of the source parameters and $d_{L}$ is
the luminosity distance. The source-frame GW energy spectrum emitted by each
astrophysical source is given by
${dE_{\text{gw}}(f_{\text{s}},\theta)}/{df_{\text{s}}}$, where
$f_{\text{s}}=f(1+z)$ is the source-frame frequency.

For astrophysical sources distributed isotropically across the Universe, the
rate of observed signals $R_{o}$ in some redshift slice $dz$ is

$\begin{multlined}\frac{dR_{o}(z)}{dz}=\frac{dN(z)}{dz\>dt_{o}}=\frac{dN(z)}{dV_{c}\>dt_{s}}\frac{dt_{s}}{dt_{o}}\frac{dV_{c}}{dz}(z)\\\
=R_{\nu}(z)\frac{1}{1+z}\frac{dV_{c}}{dz}(z),\end{multlined}\frac{dR_{o}(z)}{dz}=\frac{dN(z)}{dz\>dt_{o}}=\frac{dN(z)}{dV_{c}\>dt_{s}}\frac{dt_{s}}{dt_{o}}\frac{dV_{c}}{dz}(z)\\\
=R_{\nu}(z)\frac{1}{1+z}\frac{dV_{c}}{dz}(z),$ (7)

where $N$ specifies the number of events occurring in a given cosmic slice and
$dt_{s}/dt_{o}=1/(1+z)$ accounts for time dilation between the source and
observer. $R_{\nu}(z)=dN(z)/dV_{c}\>dt_{s}$ is the source-frame rate density
per comoving volume $V_{c}$ of a specific type of astrophysical source $\nu$.
The differential comoving volume element is given by Dominik  _et al._ (2015)

$\frac{dV_{c}}{dz}(z)=\frac{4\pi d_{L}^{2}}{(1+z)^{2}}\frac{c}{H_{0}E(z)},$
(8)

where, for a flat $\Lambda$CDM cosmology and ignoring radiation density,
$E(z)=\sqrt{\Omega_{M}(1+z)^{3}+\Omega_{\Lambda}}$ Regimbau (2011, 2022).
Bringing this all together, we arrive at the energy-density spectrum for a
SGWB arising from astrophysical sources distributed isotropically across the
Universe between some redshifts $z_{\rm low}$ to $z_{\rm up}$:

$\Omega_{\text{GW}}(f)=\frac{f}{\rho_{c}H_{0}}\int_{\theta}\int_{z_{\text{low}}}^{z_{\text{up}}}\frac{R_{\nu}(z)\frac{dE_{\text{gw}}(f_{\text{s}},\theta)}{df_{\text{s}}}p(\theta)}{(1+z)E(z)}d\theta
dz.$ (9)

Finally, the spectral energy $E_{\text{gw}}$ emitted by any astrophysical
source can be related to the amplitudes of the plus ($+$) and cross ($\times$)
GW polarizations $\tilde{h}_{+}$ and $\tilde{h}_{\times}$ via Phinney (2001)

$\frac{dE_{\rm
gw}}{df_{s}}=\frac{2\pi^{2}c^{3}d^{2}_{L}f^{2}}{G(1+z)^{2}}\biggr{<}|\tilde{h}_{+}(f,\theta)|^{2}+|\tilde{h}_{\times}(f,\theta)|^{2}\biggr{>}_{\Omega_{s}},$
(10)

where the right-hand side is averaged over the source orientations
${\Omega_{s}}$.

##  IV The CBC population models

In the previous sections, we provided the motivation and the methods to
calculate $\Omega_{\text{cbc, unres}}$, the unresolved SGWB from CBCs. Now we
describe the models we use for the astrophysical population of CBCs, focusing
in particular on detectability by XG GW detectors. Several studies have shown
that almost all of the BBHs will be individually resolvable by these detectors
Evans  _et al._ (2021, 2023); Zhou  _et al._ (2023) and that the unresolvable
background from BBHs will be several orders of magnitude smaller than the
corresponding background from BNS. Therefore, we assume that the unresolvable
SGWB arising from the BBH population is negligible and limit ourselves to the
NSBH and BNS populations.

Since XG detectors will likely be able to detect all compact binaries that
merge at low redshifts, such events will not contribute to the unresolved CBC
SGWB. On the other hand, mergers at high redshifts will often not be
individually resolvable, making them important contributors to the unresolved
CBC SGWB. Therefore, we incorporate realistic star-formation-based redshift
distribution models paired with a model of merger delay-time distributions.
This is described further in Sec. IV.1.

For all considered astrophysical populations, we assume that tidal effects and
any eccentricity effects are negligible as they will be subdominant. The GW
signal of such binaries can then be described by a set of 15 source parameters
$\theta$, of which 8 are intrinsic and 7 are extrinsic. Intrinsic parameters
include the component masses $m_{1}$ and $m_{2}$ and the three-dimensional
spin vectors $\vec{\chi}_{1}$ and $\vec{\chi}_{2}$. Since spins are expected
to have a subdominant effect on the SGWB, we also set them to zero Zhou  _et
al._ (2023); Safarzadeh  _et al._ (2020). Extrinsic parameters include the
redshift $z$, the right ascension and declination ($\alpha,\delta$), the
polarization angle $\psi$, the inclination angle $\iota$, the coalescence
phase $\phi_{c}$, and the coalescence time $t_{c}$ 333Any choice of $t_{c}$
does not affect the calculation of $\Omega_{\text{cbc, unres}}$, since it
depends only on the GW power.. We draw
$\text{cos}\;\iota,\>\text{cos}\;\delta,\>\alpha,\>\psi,\>\text{and}\>\phi_{c}$
from uniform distributions.

The rate and mass models are described in the following subsections and the
full distributions are summarized in Tab. 1.

###  IV.1 The redshift distributions

Compact binaries experience a delay time $t_{d}$ between their formation at
$z_{f}$ and merger at $z$, where $z_{f}$ is the zero-age main-sequence (ZAMS)
“formation” redshift. We calculate the delay time as the difference between
the cosmological lookback time Hogg (1999) at $z_{f}$ and $z$,

$t_{d}=t_{L}(z_{f})-t_{L}(z),$ (11)

We assume a delay time distribution of $p(t_{d})\propto t_{d}^{-1}$ as
suggested by stellar evolution and population synthesis models Piran (1992);
Sana  _et al._ (2012). The GW and Galactic pulsar observations Fishbach and
Kalogera (2021); Beniamini and Piran (2019) are also consistent with such a
distribution, although Galactic populations likely harbor an excess of sources
that will merge rapidly. An analysis using localized short gamma-ray bursts
potentially finds steeper time-delay distributions, but it does not include
selection effects Zevin  _et al._ (2022). We set the maximum delay time
$t_{d}^{max}$ to the Hubble time to limit ourselves to binaries that will
merge in the age of the Universe. We set a fiducial minimum delay time
$t_{d}^{min}$ of 20 Myr, which is approximately how long massive binaries take
to evolve into two neutron stars Regimbau  _et al._ (2012); Belczynski and
Kalogera (2001); Belczynski  _et al._ (2006). It is also consistent with
observations of binary pulsar merger times and of short gamma-ray bursts in
both late- and early-type galaxies Regimbau  _et al._ (2012); Berger  _et al._
(2007).

We convolve $t_{d}$ with a star formation rate (SFR) model $R_{f}(z_{f})$ to
calculate the source-frame CBC merger rate density per comoving volume
Fishbach and Kalogera (2021)

$R_{m}(z)\propto\int_{t_{d}^{min}}^{t_{d}^{max}}R_{f}(\tilde{z}[t_{L}(z)+t_{d}])p(t_{d})dt_{d},$
(12)

where $\tilde{z}$ is the formation redshift. We normalize Eq. 12 as

$R_{\nu,cbc}(z)=\frac{\mathcal{R}_{0}}{R_{m}(z=0)}R_{m}(z),$ (13)

such that $R_{\nu,cbc}(z=0)=\mathcal{R}_{0}$ is the inferred local source-
frame CBC merger rate density per comoving volume obtained from population
analysis (see Sec. IV.2) for CBCs of type $\nu$. Using Eq. 7, we calculate the
observed differential CBC merger rate

$\frac{dR_{o,cbc}(z)}{dz}=R_{\nu,cbc}(z)\frac{1}{1+z}\frac{dV_{c}}{dz}(z).$
(14)

We simulate CBCs up to redshift $z=10$ as we expect minimal star formation
beyond this redshift and therefore no CBCs Callister  _et al._ (2020). Our
fiducial SFR model is the Madau-Fragos SFR Madau and Fragos (2017),

$R_{f}(z_{f})\propto\frac{(1+z_{f})^{2.6}}{1+(\frac{1+z_{f}}{3.2})^{6.2}}.$
(15)

###  IV.2 The mass distributions and rates

We now describe the mass distributions and local merger rate densities
$\mathcal{R}_{0}$ we use to simulate the CBC population. We consider several
population models, fit to LVK data up until the end of the third observing run
via Bayesian inference. We also consider models of mass distribution fit to
Galactic double neutron stars. Considering several models enables us to
characterize systematic modeling uncertainties that could come from using
strongly parameterized models. Each model has a statistical uncertainty on the
population that we marginalize over. Additionally in each model, a single set
of hyperparameters $\Lambda_{\text{GW}}$ characterizes a single CBC
population. We assume that only the overall rate, not the mass distribution,
evolves with redshift according to Eq. 12.

For each model, we self-consistently include the rates inferred by it when
available May . We refer to such models as “rate-inclusive” models. For models
where such information is not available – which we refer to as “rate-
marginalized” models – we randomly pair each hyperparameter sample with a
$\mathcal{R}_{0}$ posterior from an LVK O3b population model Abbott  _et al._
(2023) that estimates rates. For clarity, we rewrite Eq. 14 for cases where
$\mathcal{R}_{0}$ depends on $\Lambda_{\text{GW}}$ (i.e.,
$\mathcal{R}_{0}(\Lambda_{\text{GW}})$) as

$\frac{dR_{o,cbc}(z,\Lambda_{\text{GW}})}{dz}=R_{v,cbc}(z,\Lambda_{\text{GW}})\frac{1}{1+z}\frac{dV_{c}}{dz}(z).$
(16)

The various mass-rate model configurations are described below and summarized
in Tab. 1.

| Model
---|---
Parameter |  pdb-pl |  nsbh-pl |  nsbh-g |  bns-g |  bns-pl
$\mathcal{R}_{0}(\Lambda_{\text{GW}})$ |  Rate-inclusive |  Rate-inclusive |  Rate-inclusive |  Drawn from pdb-ind analysis Farah  _et al._ (2022); Fishbach  _et al._ (2020a); Abbott  _et al._ (2023, 2021c) |  Drawn from pdb-ind analysis Farah  _et al._ (2022); Fishbach  _et al._ (2020a); Abbott  _et al._ (2023, 2021c)
${m}_{1}$ |  PDB with power-law pairing Farah  _et al._ (2022); Abbott  _et al._ (2023, 2021c) |  Power Law Biscoveanu  _et al._ (2022) |  Power Law Biscoveanu  _et al._ (2022) |  Double Gaussian Farrow  _et al._ (2019) |  Power Law Abbott  _et al._ (2023); Landry and Read (2021); Abbott  _et al._ (2021c)
${m}_{2}$ |  PDB with power-law pairing Farah  _et al._ (2022); Abbott  _et al._ (2023, 2021c) |  Power Law Pairing Biscoveanu  _et al._ (2022) |  Gaussian Pairing Biscoveanu  _et al._ (2022) |  Uniform Farrow  _et al._ (2019) |  Random Pairing Abbott  _et al._ (2023); Landry and Read (2021); Abbott  _et al._ (2021c)
$z$ | Madau-Fragos SFR Madau and Fragos (2017)
$t^{min}_{d}$ | 20 Myr
$\vec{\chi}_{1}$, $\vec{\chi}_{2}$ | 0
cos $\iota$, cos $\delta$ | Uniform in $[-1,1]$
$\alpha$, $\psi$, $\phi_{c}$ | Uniform in $[0,\boldmath 2\pi]$
Table 1: Different populations of the various source parameters used in
calculating the SGWB from CBCs. See Section IV for further details.

####  IV.2.1 Neutron-star black-hole models

We model the NSBH mass distribution using the Bayesian population analyses
from Ref. Biscoveanu  _et al._ (2022). We assume that the distribution of the
primary mass ${m}_{1}$ follows the truncated power-law Fishbach and Holz
(2017)

$\pi(m_{1}|\gamma,m_{1,\text{min}},m_{1,\text{max}})\propto$

$\begin{cases}m_{1}^{-\gamma},\;\text{if}\;m_{1,\text{min}}\leq m_{1}\leq
m_{1,\text{max}}\\\ 0,\;\text{otherwise}\end{cases},$ (17)

with a power-law index $\gamma$, minimum $m_{1}$ cutoff $m_{1,\text{min}}$,
and maximum $m_{1}$ cutoff $m_{1,\text{max}}$.

We consider two pairing functions to get the distribution of the secondary
mass. The first is a power-law pairing function Fishbach and Holz (2020),
which we refer to as the nsbh-pl model

$\pi(q|\beta,m_{1},m_{2,\text{max}})\propto$

$\begin{cases}q^{\beta},\;\text{if}\;q_{\text{min}}(m_{1})\leq q\leq
q_{\text{max}}(m_{1},m_{2,\text{max}}))\\\ 0,\;\text{otherwise}\end{cases},$
(18)

where $q=m_{2}/m_{1}$ is the mass ratio and $\beta$ is a power-law index. The
second is a truncated-Gaussian pairing function Biscoveanu  _et al._ (2022),
which we refer to as the nsbh-g model

$\pi(q|\mu,\sigma,m_{1},m_{2,\text{max}})\propto$

$\begin{cases}\mathcal{N}(q|\mu,\sigma),\;\text{if}\;q_{\text{min}}(m_{1})\leq
q\leq q_{\text{max}}(m_{1},m_{2,\text{max}}))\\\
0,\;\text{otherwise}\end{cases},$ (19)

where $\mathcal{N}(q|\mu,\sigma)$ is a Gaussian with mean $\mu$ and standard
deviation $\sigma$. For both models, we set the minimum NS mass
$m_{2,\text{min}}=1\ M_{\odot}$ so that the minimum mass ratio cutoff
$q_{\text{min}}=1/m_{1}$. The maximum mass ratio cutoff is set as
$q_{\text{max}}=\text{min}(m_{2,\text{max}}/m_{1},1)$, where the maximum NS
mass $m_{2,\text{max}}$ is a free parameter drawn uniformly between $1.97\
M_{\odot}$ and $2.7\ M_{\odot}$. While the ranges of $m_{1}$ and $m_{2}$
differ based on the particular hyperparameter values, the broadest possible
range is $m_{1}\in[2,20]\ M_{\odot}$ and $m_{2}\in[1,2.7]\ M_{\odot}$. We
refer the reader to Ref. Biscoveanu  _et al._ (2022) for a more detailed
overview of the models.

Both the nsbh-pl and nsbh-g models are rate-inclusive models that self-
consistently calculate $\mathcal{R}_{0}$. Hence, for each hyperparameter
sample in both analyses, we set $\mathcal{R}_{0}$ to the associated posterior
$\mathcal{R}_{0}(\Lambda_{\text{GW}})$.

####  IV.2.2 Power Law + Dip + Break model

We next consider the Power Law + Dip + Break (PDB) model used in the LVK
GWTC-3 analysis Farah  _et al._ (2022); Fishbach  _et al._ (2020a); Abbott
_et al._ (2023, 2021c). This model fits the entire mass spectrum of CBCs
eschewing the difference between black holes and neutron stars. The PDB model
fits both the primary and secondary masses using a common distribution
$\pi_{\rm pdb}(m|\vec{\Lambda})$ that is based on a broken-power law and
includes a “notch filter” to model a potential mass gap between $3\text{ to
}5M_{\odot}$. We refer the reader to Ref. Farah  _et al._ (2022) for a
detailed overview of the model.

We then use a power-law pairing function to pair the primary and secondary
masses Fishbach and Holz (2020); Abbott  _et al._ (2023):

$g(m_{1},m_{2},\beta)=\begin{cases}q^{\beta},\;\text{if}\;m_{2}\leq m_{1}\\\
0,\;\text{otherwise}\end{cases}.$ (20)

This gives the following joint mass distribution, which we refer to as the
pdb-pl model:

$\pi_{\textsc{pdb-pl}}(m_{1},m_{2}|\vec{\Lambda},\beta)\propto\pi_{\rm
pdb}(m_{1}|\vec{\Lambda})\pi_{\rm
pdb}(m_{2}|\vec{\Lambda})g(m_{1},m_{2},\beta)\;.$ (21)

While pdb-pl is a rate-inclusive model, we need to account for the fact that
the model natively incorporates BBH mergers as well. In order to calculate
rates that only include the NSBH and BNS parts of the spectrum, we correct the
rates by the fraction of NSBH and BNS mergers that this model predicts (at
each point in the hyperparameter space). Since we need systems to have at
least one neutron star, we limit $m_{2}$ to span $[1,3]\>M_{\odot}$ while
allowing $m_{1}$ to span $[1,100]\>M_{\odot}$ as defined in the original model
Abbott  _et al._ (2023).

####  IV.2.3 Binary neutron star models

We model the BNS mass distribution using two different rate-marginalized
models. The first model is motivated through observations of Galactic double
neutron stars, which provide the largest and most well-constrained sample of
such systems. This model, which we refer to as bns-g, is based on the most
favored model from the population analysis of the 17 observed Galactic BNSs by
Ref. Farrow  _et al._ (2019). In particular, the primary neutron-star mass
follows a double-Gaussian mixture model

$\begin{split}\pi(m_{1}|\gamma,\mu_{1},\sigma_{1},\mu_{2},\sigma_{2})=&\frac{\gamma}{\sigma_{1}\sqrt{2\pi}}{\rm
e}^{-\frac{(m_{1}-\mu_{1})^{2}}{2\sigma_{1}}}\\\
&+\frac{1-\gamma}{\sigma_{2}\sqrt{2\pi}}{\rm
e}^{-\frac{(m_{2}-\mu_{2})^{2}}{2\sigma_{2}}},\end{split}$ (22)

where $\gamma$ ($1-\gamma$), $\mu_{1}$ ($\mu_{2}$), and $\sigma_{1}$
($\sigma_{2}$) give the weight, mean, and standard deviation of the first
(second) peak, respectively. The secondary neutron-star mass is drawn
uniformly within the range

$\pi(m_{2}|m_{2,\text{min}},m_{2,\text{max}})=U[m_{2,\text{min}},m_{2,\text{max}}],$
(23)

with minimum and maximum cutoff $m_{2,min}$ and $m_{2,max}$, respectively. As
usual, we assume that $m_{2}\leq m_{1}$. We note that both the double-Gaussian
and uniform distributions are truncated at 0.8 $M_{\odot}$ and 2 $M_{\odot}$
to reflect the prior ranges of the analysis from Ref. Farrow  _et al._ (2019).

While the analysis in Ref. Farrow  _et al._ (2019) relies on observations of
only Galactic BNSs, we assume in our usage of the bns-g model that it is
extendable to all redshifts and metallicities. Note, however, that the BNS
mass distribution inferred with Galactic observations is potentially
inconsistent with GW observations. This is especially highlighted by the BNS
merger GW190425 Abbott  _et al._ (2020), which has a total mass heavier than
that of typical Galactic double neutron stars and could potentially have
formed through different formation channels Romero-Shaw  _et al._ (2020).

Hence, we include a second model, which we call bns-pl, based on the BNS power
mass distribution inferred from all the neutron star systems detected through
GWs Abbott  _et al._ (2021c, 2023); Landry and Read (2021). In the bns-pl
model, both the primary and secondary masses are paired randomly after each
being drawn from a power-law distribution

$\pi(m|\gamma,m_{\text{min}},m_{\text{max}})\propto m^{\gamma},$ (24)

with a power-law index $\gamma$ and minimum (maximum) mass cutoff
$m_{\text{min}}$ ($m_{\text{max}}$). While the neutron star power model is
used in Ref. Abbott  _et al._ (2023) to infer the mass distribution from both
BNSs and NSBHs, we only consider the BNS case here.

Since our BNS models are all rate-marginalized, we model their rate
distributions by drawing $\mathcal{R}_{0}$ posteriors from the rate-inclusive
PDB random-pairing model (which we refer to as pdb-ind) Farah  _et al._
(2022); Abbott  _et al._ (2023) fit to LVK GWTC-3 data. Similar to Sec.
IV.2.2, we correct the rates inferred by pdb-ind to include only the fraction
of systems that are BNSs. This means restricting the masses to
$m_{1},m_{2}\in[1,3]\>M_{\odot}$ for the bns-pl model, but to
$m_{1},m_{2}\in[1,2]\>M_{\odot}$ for the bns-g model to be consistent with the
mass ranges from Refs. Farrow  _et al._ (2019); Abbott  _et al._ (2023).

Figure 1:  Amplitude spectral densities $\sqrt{S_{n}}$ for the detectors in
our networks.

##  V Detector Networks

We consider four different possible networks of XG detectors, using the latest
detector designs described in Refs. Evans  _et al._ (2023); Gupta  _et al._
(2023). Using the GW simulation package Gwbench Borhanian (2021), we consider
the CE-40 and CE-20 options for the proposed 40-km arm and 20-km arm CE
sensitivities respectively Evans  _et al._ (2023); Gupta  _et al._ (2023), the
A# option for the proposed 4-km arm $\text{A}^{\\#}$ sensitivity Fritschel
_et al._ (2022), and the ET-10-XYL option for the proposed triangular
“xylophone” 10-km arm ET sensitivity ET_ . Figure 1 shows the amplitude
spectral density $\sqrt{S_{n}(f)}$ of each detector considered.

Since an SGWB is detected by cross-correlating data observed separately by at
least two detectors, we consider only multiple-detector networks Allen and
Romano (1999). We consider the CE-A (coast of Washington, USA), CE-B (coast of
Texas, USA), ETS (slightly South of Virgo’s current location at Cascina,
Italy), and LLO (current location of the LIGO-Livingston Observatory at
Livingston, Louisiana, USA) facility locations specified in Tab. II of Ref.
Gupta  _et al._ (2023), which we point to for the details.

We consider four different networks of detectors:

  * •

A fiducial three-detector network including one CE-40 at CE-A, one CE-20 at
CE-B, and one ET-10-XYL at ETS.

  * •

An alternate three-detector network including one CE-40 at CE-A, one ET-10-XYL
at ETS, and one A#-upgraded detector at LLO.

  * •

An alternate two-detector network including one CE-40 at CE-A, and one CE-20
at CE-B.

  * •

An alternate two-detector network including one CE-40 at CE-A, and one
ET-10-XYL at ETS.

We choose the same three-detector networks as those proposed in Refs. Evans
_et al._ (2023); Gupta  _et al._ (2023) to make it easy to compare with
existing and future literature. While the future development of ET has been
officially confirmed, we specifically explore a two-detector network
configuration that does not include ET in order evaluate the usefulness of and
aid in the detector design of CE facilities in detecting astrophysically and
cosmologically-arising SGWBs even in the absence of ET.

Throughout the paper, we set a minimum frequency of 5 Hz, corresponding to the
proposed lower limit of CE and $\text{A}^{\\#}$ design sensitivities Evans
_et al._ (2023); Gupta  _et al._ (2023); Fritschel  _et al._ (2022) as shown
in Fig. 1. In addition, we set a maximum frequency of 2000 Hz. For each
network, we use the same $3\sigma$ power-law integrated (PI) curves Thrane and
Romano (2013); Zhou  _et al._ (2023) as Ref. Gupta  _et al._ (2023)
444Obtained through private communication with the authors. in order to
measure the ability of the network to detect SGWB signals. We refer the reader
to Ref. Gupta  _et al._ (2023) for the overlap reduction functions Thrane and
Romano (2013); Zhou  _et al._ (2023) of the various detector pairs comprising
the networks described in this section.

##  VI Simulating the SGWB : A Monte-Carlo approach

In the previous sections, we have described the population models and the
detector networks that we use. In this section, we describe the simulation and
the calculation of the unresolved SGWB using Eq. 9 and Eq. 10 for which we
adopt a Monte-Carlo approach Regimbau (2022); Regimbau  _et al._ (2017, 2012,
2014); Meacher  _et al._ (2015). We use the Gwbench simulation platform
Borhanian (2021) to generate GW waveforms for non-spinning quasi-circular
binaries neglecting tidal effects. Gwbench naturally accounts for the Earth’s
rotation and its impact on the antenna patterns which is important since BNSs
can last for several hours in the observing band of XG detectors Borhanian
(2021). We use IMRPhenomD Husa  _et al._ (2016); Khan  _et al._ (2016) as our
waveform approximant 555Note that the choice of waveform approximant is
expected to be subdominant in calculating the SGWB Zhou  _et al._ (2023)..

For each population configuration described in Sec. IV and summarized in Tab.
1, we are interested in calculating $\Omega_{\text{cbc, unres}}$. To account
for astrophysical uncertainties, we draw 2000 hyperparameter samples from the
inferred population distributions and estimate $\Omega_{\text{cbc, unres}}$
for each sample.

To do this in a computationally tractable way, we first simulate $10^{5}$
waveforms each for NSBHs and BNSs drawn from broad fiducial population
distributions 666In this fiducial population, for NSBHs we uniformly draw
$m_{1}\sim[2,50]$ and $m_{2}\sim[1,3]$ and for BNSs we uniformly draw both
$m_{1},m_{2}\sim[0.5,3]$.. We then apply rejection sampling to draw a
population corresponding to any particular hyperparameter draw. We estimate
the mean number of sources needed for a CBC model $\nu$ with a set of
hyperparameters $\Lambda_{\text{GW}}$ and observation time $T$ as

$\langle
N_{\nu}(\Lambda_{\text{GW}})\rangle=\int^{z_{\text{max}}}_{0}R_{\nu}(z,\Lambda_{\text{GW}})\frac{dV_{c}}{dz}\frac{dz}{1+z}\;T,$
(25)

and draw the actual number of sources through a Poisson draw

$N_{\nu}(\Lambda_{\text{GW}})\sim\texttt{Poiss}(\lambda=\langle
N_{\nu}(\Lambda_{\text{GW}})\rangle).$ (26)

In order to estimate $\Omega_{\text{cbc, unres}}$, we first need to extract
only the unresolvable signals from this population. To do this, we first
compute the matched-filter signal-to-noise ratio (SNR) $\rho_{j}^{\rm mf}$ of
the signal, defined for detector $j$ with noise PSD $S_{n,j}(f)$ as Moore  _et
al._ (2015):

$\rho_{j}^{\rm
mf}=\left[4\int_{0}^{\infty}\frac{|\tilde{h}_{j}(f)|^{2}}{S_{n,j}(f)}df\right]^{1/2}.$
(27)

This is, however, the optimal matched-filter SNR; in order to account for the
measurement uncertainty due to detector noise, we correct the SNR as Fishbach
_et al._ (2020b); Callister  _et al._ (2020); Fishbach  _et al._ (2018)

$\rho_{j}^{\rm obs}=\rho_{j}^{\rm mf}+\mathcal{N}(0,1).$ (28)

We then define the optimal network SNR $\rho_{\text{net}}^{\rm obs}$ for a
network of $D$ detectors by summing the individual observed SNRs in quadrature
as

$\rho_{\text{net}}^{\rm obs}=\sqrt{\sum^{D}_{j=1}(\rho^{\rm obs}_{j})^{2}}.$
(29)

A signal is then labeled as resolved if $\rho_{\text{net}}^{\rm obs}$ is
greater than some frequency-independent threshold $\rho_{\text{thresh}}$.

Once we have a population $N_{\nu}^{\rm unres}(\Lambda_{\text{GW}})$ of
unresolvable sources with respect to $\rho_{\text{thresh}}$, we calculate the
GW energy flux (Eq. 6) using a Monte-Carlo sum over our simulated population
Zhou  _et al._ (2023); Regimbau  _et al._ (2017); Meacher  _et al._ (2015):

$F_{\nu}(f;\Lambda_{\text{GW}})\approx\frac{\pi
c^{3}f^{2}}{2GT}\sum_{i=1}^{N_{\nu}^{\text{unres}}(\Lambda_{\text{GW}})}\biggr{[}|\tilde{h}_{+}^{i}(f,\theta^{i})|^{2}+|\tilde{h}_{\times}^{i}(f,\theta^{i})|^{2}\biggr{]}.$
(30)

(a) $\Omega_{\text{cbc, unres}}$ estimates for the various BNS models in Tab.
1. The bands show the central $90\%$ credible levels.

(b) $\Omega_{\text{cbc, unres}}$ estimates for the various NSBH models in Tab.
1. The bands show the central $90\%$ credible levels.

Figure 2: This plot shows the $\Omega_{\text{cbc, unres}}$ estimate for the
various BNS and NSBH models for the fiducial 3-detector network described in
Sec. V. The purple band in both plots is the pdb-pl model that simultaneously
incorporates both the NSBH and BNS systems. Also overlaid is the $3\sigma$ PI
curve that shows that the unresolved CBC backgrounds will likely be very loud.
Figure 3: A comparison of various cosmological SGWB models along with
$\Omega_{\text{cbc, unres}}$ estimates for the fiducial network described in
Sec. V. The blue band shows the pdb-pl model while the pink band shows the
joint-envelope constructed from the various BNS-only and NSBH-only models in
Tab. 1. Both bands represent $90\%$ credible intervals. The solid black curve
is the $3\sigma$ PI curve.

##  VII Cosmological SGWB Models

To gauge the impact of the $\Omega_{\text{cbc, unres}}$, we consider several
models of early Universe cosmological SGWBs that could be accessible with XG
detectors in principle. These models include SGWBs from domain walls Saikawa
(2017), first-order phase transition (FOPT) sound waves Romero  _et al._
(2021); Espinosa  _et al._ (2010), FOPT bubble collisions Romero  _et al._
(2021), stiff equation of state Figueroa and Tanin (2019), and preheating
Easther  _et al._ (2007); Gupta  _et al._ (2023). We also use the model
spectrum for Nambu-Goto oscillating cosmic string loops (the model C1 from
Ref.Abbott  _et al._ (2021e)).

We note that XG detectors will not be sensitive to standard slow-roll
inflation, for which $\Omega_{\text{GW}}\sim 10^{-15}-10^{-17}$ depending on
the model Caprini and Figueroa (2018); Zhou  _et al._ (2022). In general,
these various SGWBs are highly sensitive to the choice of model parameters,
but in our plots, we chose illustrative curves for $\Omega_{\rm cosmo}$.

##  VIII Results

As described in Sec. VI, for each population model in Tab. 1, we have 2000
different draws for the unresolved CBC SGWB. These draws incorporate both the
uncertainty in the rate of mergers and the astrophysical uncertainty from the
population model. In Fig. 2, we show the $90\%$ credible bands for the
unresolved SGWB for the various BNS-only and NSBH-only population models, for
the fiducial 3-detector network from Sec. V using a threshold SNR of 12. We
overlay them on top of the $3\sigma$ PI curve for this network. We assume an
observation time of $T=1$ year throughout the paper. The purple band in both
plots is from the pdb-pl model and is therefore higher than either of the BNS-
only and NSBH-only models.

The unresolved SGWB from CBCs is above the $3\sigma$ PI curve, implying that
it will very likely be detectable by XG detectors with one year of observing.
In general, the width of each band in Fig. 2 – which represents the total
uncertainty in our estimate of the unresolved SGWB for that model – spans
about an order of magnitude and the bands fall out of detectability at $\sim
80$ Hz even with XG detectors. The combined pdb-pl background, however, might
be significantly detectable up to $\sim 200$ Hz. Nevertheless, it would seem
that we will only observe the SGWB due to the inspiral phase of the mergers.
Similarly, while the simulations show a clear turnover for the NSBH SGWB – the
morphology of which probably depends on the waveform used – it clearly happens
at high frequencies and XG detectors will probably not be able to probe this.

Figure 4: The SGWB for different BNS models evaluated at three different rates
using a threshold SNR of 12. We use the following NS-system rates (corrected
BNS rates for bns-pl and bns-g): Rate 1 = 45.5 (32.2)
$\text{Gpc}^{-3}\>\text{Yr}^{-1}$, Rate 2 = 209.3 (186.6)
$\text{Gpc}^{-3}\>\text{Yr}^{-1}$, Rate 3 = 377.0 (362.6)
$\text{Gpc}^{-3}\>\text{Yr}^{-1}$. While there is significant difference
between the curves for various rates, the differences between the models is
minimal.

In addition to the pdb-pl model which naturally incorporates both BNS and NBSH
mergers, we create a maximum-uncertainty envelope, which we refer to as the
joint-envelope as a means of accounting for model systematics. This is
generated by combining the BNS-only and NSBH-only models (i.e., nsbh-pl,
nsbh-g, bns-pl, bns-g) in Tab. 1 in all possible combinations and choosing the
widest band at each frequency bin.

Figure 3 shows the $90\%$ credible band of the joint-envelope along with the
pdb-pl estimate of $\Omega_{\text{cbc, unres}}$ for our fiducial 3-detector
network. We contrast these bands with several cosmological SGWB models,
demonstrating that the unresolved background will likely be a major source of
noise for accessing the said cosmological signals. Since the effective CBC
background is the sum of unresolved and resolved components (see Eq. 3), it
also depends on the efficacy of subtraction of resolved signals. Since the
subtraction residue could contribute significantly to the effective CBC
background (see for e.g. Zhong  _et al._ (2023); Sharma and Harms (2020)), the
unresolved background we show in Fig. 3 therefore represents the floor of the
noise from CBC sources. Either way, a joint simultaneous analysis of the CBC
background and the cosmological background will be necessary for the detection
of the latter. We note again that the cosmological SGWB curves themselves can
change depending on the choice of parameters. In Appendix. A, we present
similar results for the three alternate detector networks described in Sec. V.

Further exploring the uncertainty in the $\Omega_{\text{cbc, unres}}$, we see
remarkable consistency between the two NSBH-only models, and similarly between
the two BNS-only models. This suggests that the uncertainty in the rate is
dominant, as opposed to the astrophysical uncertainty on the mass
distribution. This is further confirmed by Fig. 4, where we plot the SGWB from
different BNS models for three different rates. The consistency in SGWB
between models at any given frequency shows that the astrophysical uncertainty
is a much smaller contributor to the SGWB uncertainty than the uncertainty in
CBC rate. We find similar results for NSBHs as well. Since the number of BNSs
and NSBHs observed thus far is fairly small Abbott  _et al._ (2019a, 2021b,
2021c), this means that any future detections in O4 and O5 can substantially
lower the uncertainty of these bands. This is also consistent with studies
such as Ref. Jenkins  _et al._ (2019), which show that the monopole of the CBC
SGWB is particularly sensitive to the local merger rate.

Figure 5: 90% credible intervals of $\Omega_{\text{cbc, unres}}$ for the pdb-
pl model with three different SFR assumptions. While the green and the saffron
bands are for a gamma-ray burst based SFR and the Madau-Fragos SFR
respectively, the pink band is for a conservative 1 Gyr minimum time delay.

To explore how SFR models impact our estimates of $\Omega_{\text{cbc,
unres}}$, we also used an SFR extracted from long gamma-ray burst rates
following Zhou  _et al._ (2023); Vangioni  _et al._ (2015):

$R_{f}(z_{f})\propto\nu\frac{ae^{b(z_{f}-z_{p})}}{a-b+be^{a(z_{f}-z_{p})}},$
(31)

where $\nu=0.146\>M_{\odot}\>{\text{yr}}^{-1}\>{\text{Mpc}}^{3}$,
$z_{p}=1.72$, $a=2.80$, and $b=2.46$. As the pdb-pl bands in Fig. 5 show, we
find that the impact of the SFR model is minimal as long as the minimum delay
times are small. A minimum delay time of $1$ Gyr gives a significantly
smaller, albeit still loud, SGWB. Note, however, that studies of Galactic
systems have actually found an excess of systems with small delay times
Beniamini and Piran (2019), meaning that the value of $1$ Gyr is very
conservative. We therefore conclude that the key observation of this paper,
that the unresolved background will significantly impact searches for
cosmological backgrounds, is likely robust to uncertainties in SFR. We have
also tested our assumption of $z_{\rm max}=10$ and find that the GW power from
higher redshifts is very small.

The LVK population inference at the end of O3 tested multiple versions of the
power model Abbott  _et al._ (2023) including only confident events, or
confident and marginal events, or confident events and GW190814 777which is a
population outlier Essick  _et al._ (2022). While we use the last version in
all our plots as the bns-pl model, we find that the impact of this choice is
negligible. In particular, our estimates of $\Omega_{\text{bns, unres}}$ for
our implemented version differs from those for the other two versions at most
by $3\%$ at 25 Hz, demonstrating that our estimates are robust with respect to
this uncertainty.

Finally, we also estimate the CBC SGWB at multiple $\rho_{\rm thresh}$, which
determines what we define as a resolvable signal. We note that our fiducial
threshold of $\rho_{\rm thresh}=12$ is somewhat conservative for the purposes
of the unresolved SGWB. Figure 6 shows the pdb-pl and the joint envelopes for
two different SNR choices. The threshold $\rho_{\rm thresh}=20$, in
particular, was shown to be the frequency-independent optimal threshold for
minimizing the effective background from BNS systems, including both the
subtraction residue and the unresolved background Zhou  _et al._ (2023).

Figure 6: 90% credible intervals of $\Omega_{\text{cbc, unres}}$ for the pdb-
pl and the joint-envelope models with two different $\rho_{\rm thresh}$
values. The filled bands show the credible intervals with $\rho_{\rm
thresh}=12$ and the dashed-dotted lines bound the intervals for $\rho_{\rm
thresh}=20$.

##  IX Implications and Conclusions

In this paper, we have presented a data-driven estimate of the strength of the
unresolvable CBC background seen by XG detectors, accounting for both
statistical and systematic uncertainties. A robust implication from our study
is that $\Omega_{\rm cbc,unres}$ will be an impediment for the direct
detection of many cosmological SGWB models, independent of the fidelity and
efficacy of the subtraction of resolved CBC signals. Some manner of
simultaneous inference of the astrophysical and cosmological background (e.g.,
see Refs. Martinovic  _et al._ (2021); Yu and Wang (2022); Kaiser  _et al._
(2022)) will be required at a minimum.

Given the large number of resolvable CBC signals that we expect, it is also
possible that their population inference can be used as a strong prior for
$\Omega_{\rm cbc,unres}$. While this could make inferring multiple backgrounds
easier, this effort could be hindered by the efficacy of subtraction.
Moreover, using the resolvable signals to form a prior could create a
systematic bias because the resolvable and unresolvable sources come from
different redshifts Callister  _et al._ (2016).

However, this also makes $\Omega_{\text{cbc, unres}}$ astrophysically
interesting in its own right by allowing us to study compact binary
populations at high redshifts Callister  _et al._ (2020); Mukherjee  _et al._
(2022). Since we will confidently detect $\Omega_{\rm cbc,unres}$ with XG
detectors, this represents a clear science case that can inform detector
design. For this purpose, in Tab. 2 we provide figures-of-merit of the level
of $\Omega_{\text{cbc, unres}}$ for various detector networks at two reference
frequencies.

Our results also help inform the target efficacy of subtraction techniques for
XG detectors. The optimal level of subtraction is one that leaves a residue
not significantly higher than $\Omega_{\rm cbc,unres}$. At the same time, a
subtraction technique that leaves a residue substantially lower than
$\Omega_{\rm cbc,unres}$ will not be ideal and likely computationally
expensive from a statistical point of view. Our figures-of-merit for the
various networks in Tab. 2 will again be useful in informing the optimal
efficacy of subtraction that algorithms should target.

| $\log_{10}\Omega_{\rm cbc,unres}$
---|---
Model | $f$ | $\rho_{\rm thresh}=12$ | $\rho_{\rm thresh}=20$ | $\rho_{\rm thresh}=40$
| 40 km CE A, 20 km CE B, 10 km ET |
joint-envelope | 25 Hz | $-10.8^{+0.4}_{-0.4}$ | $-10.3^{+0.4}_{-0.4}$ | $-10.0^{+0.3}_{-0.4}$
| 250 Hz | $-10.3^{+0.5}_{-0.4}$ | $-9.9^{+0.4}_{-0.4}$ | $-9.6^{+0.4}_{-0.4}$
pdb-pl | 25 Hz | $-10.4^{+0.4}_{-0.5}$ | $-10.0^{+0.4}_{-0.4}$ | $-9.8^{+0.4}_{-0.4}$
| 250 Hz | $-9.8^{+0.4}_{-0.5}$ | $-9.5^{+0.4}_{-0.5}$ | $-9.3^{+0.4}_{-0.4}$
| 40 km CE A, 10 km ET, 4 km A# LLO |
joint-envelope | 25 Hz | $-10.6^{+0.4}_{-0.4}$ | $-10.3^{+0.4}_{-0.4}$ | $-10.0^{+0.4}_{-0.4}$
| 250 Hz | $-10.2^{+0.4}_{-0.4}$ | $-9.8^{+0.4}_{-0.4}$ | $-9.6^{+0.4}_{-0.4}$
pdb-pl | 25 Hz | $-10.3^{+0.4}_{-0.5}$ | $-10.0^{+0.4}_{-0.5}$ | $-9.8^{+0.4}_{-0.5}$
| 250 Hz | $-9.8^{+0.4}_{-0.5}$ | $-9.4^{+0.4}_{-0.5}$ | $-9.3^{+0.4}_{-0.5}$
| 40 km CE A, 20 km CE B |
joint-envelope | 25 Hz | $-10.6^{+0.4}_{-0.4}$ | $-10.2^{+0.4}_{-0.4}$ | $-10.0^{+0.4}_{-0.3}$
| 250 Hz | $-10.1^{+0.4}_{-0.4}$ | $-9.8^{+0.4}_{-0.4}$ | $-9.6^{+0.4}_{-0.3}$
pdb-pl | 25 Hz | $-10.2^{+0.4}_{-0.5}$ | $-9.9^{+0.4}_{-0.5}$ | $-9.8^{+0.4}_{-0.4}$
| 250 Hz | $-9.7^{+0.4}_{-0.5}$ | $-9.4^{+0.4}_{-0.5}$ | $-9.2^{+0.4}_{-0.4}$
| 40 km CE A, 10 km ET |
joint-envelope | 25 Hz | $-10.6^{+0.4}_{-0.4}$ | $-10.3^{+0.4}_{-0.4}$ | $-10.0^{+0.4}_{-0.4}$
| 250 Hz | $-10.2^{+0.4}_{-0.4}$ | $-9.8^{+0.4}_{-0.4}$ | $-9.6^{+0.4}_{-0.4}$
pdb-pl | 25 Hz | $-10.3^{+0.4}_{-0.5}$ | $-10.0^{+0.4}_{-0.5}$ | $-9.8^{+0.4}_{-0.5}$
| 250 Hz | $-9.8^{+0.4}_{-0.5}$ | $-9.4^{+0.4}_{-0.5}$ | $-9.3^{+0.4}_{-0.5}$
Table 2: $\Omega_{\text{cbc, unres}}$ for different networks at two reference
frequencies for the pdb-pl model and the joint-envelope. The quoted numbers
correspond to median values with $90\%$ uncertainties. We present these
estimates as figures of merit for both detector design and subtraction
algorithms.

Instead of subtraction, one can also infer the parameters of CBC signals
without any threshold. This Bayesian “global fit” technique was first
developed in Ref. Smith and Thrane (2018), with various extensions studied in
Refs. Smith  _et al._ (2020); Talbot  _et al._ (2021); Banagiri  _et al._
(2020); Biscoveanu  _et al._ (2020). Ref. Biscoveanu  _et al._ (2020) in
particular showed how a cosmological (non-CBC) background can be incorporated
into this formalism for simultaneous inference. While this statistical
formalism holds much promise because it removes both the need for subtracting
the resolved CBC signals and for separately inferring the unresolved CBC
background, it is also computationally very expensive, which can make running
it on long stretches of data untenable. Moreover, the method is susceptible to
subtle systematics Talbot  _et al._ (2021) and can struggle to deal with
overlapping signals, both of which can make an application on real XG data
challenging.

Finally, we discuss some of the caveats of this study. One of the most
important assumptions we made is the specific form of the SFR distribution.
While we have tried two different SFR models as we show in Fig. 5,
uncertainties on the rate naturally increase at higher redshift. For instance,
observations with the James Webb Space Telescope show heightened star
formation at high redshifts, which might imply a louder CBC background, and
more observation in the near future might make the picture clearer Finkelstein
_et al._ (2023); Harikane  _et al._ (2023). Another effect that is somewhat
related, which we did not consider here is the effect of metallicity Lehoucq
_et al._ (2023); Périgois  _et al._ (2021); Kowalska-Leszczynska  _et al._
(2015); Martinovic  _et al._ (2022); Nakazato  _et al._ (2016); Dvorkin  _et
al._ (2016). This would again predominantly impact the rate at higher
redshifts.

On the GW waveform side, we have assumed zero spins and tidal effects. While
these are all subdominant compared to the mass distributions, it is possible
that the spins in particular might have a noticeable impact on $\Omega_{\rm
cbc,unres}$ given the lengths of the signals in XG detectors.

##  X Acknowledgements

We thank Salvatore Vitale, Emanuele Berti, Vuk Mandic, Haowen Zhong, and
Vishal Baibhav for helpful discussions. We thank Amanda Farah and Philippe
Landry for their help with the LVK population analyses and Sylvia Biscoveanu
for her help with the NSBH population analyses and for providing PI curves for
XG detectors. D.S.B. is supported by the National Science Foundation (NSF)
Graduate Research Fellowship Program under grant DGE-2234667. Any opinion,
findings, and conclusions or recommendations expressed in this material are
those of the authors and do not necessarily reflect the views of the NSF.
D.S.B. was supported by a CIERA Post-Baccalaureate Research Fellowship during
the majority of this research. D.S.B. acknowledges partial support for this
project from NASA through a NASA Illinois Space Grant. S.B. acknowledges
support from the NSF grant PHY-2207945. Z.D. is supported by the CIERA Board
of Visitors Research Professorship. V.K. is partially supported through a
CIFAR Senior Fellowship, a Guggenheim Fellowship, the Gordon and Betty Moore
Foundation (grant award GBMF8477), and from Northwestern University, including
the Daniel I. Linzer Distinguished University Professorship fund.

This material is based upon work supported by NSF’s LIGO Laboratory which is a
major facility fully funded by the National Science Foundation. This research
has made use of data obtained from the Gravitational Wave Open Science Center
(gwosc.org), a service of LIGO Laboratory, the LIGO Scientific Collaboration,
the Virgo Collaboration, and KAGRA. LIGO Laboratory and Advanced LIGO are
funded by the United States NSF as well as the STFC of the United Kingdom, the
Max-Planck-Society (MPS), and the State of Niedersachsen/Germany for support
of the construction of Advanced LIGO and construction and operation of the GEO
600 detector. Additional support for Advanced LIGO was provided by the
Australian Research Council. Virgo is funded, through the European
Gravitational Observatory (EGO), by the French Centre National de Recherche
Scientifique (CNRS), the Italian Istituto Nazionale di Fisica Nucleare (INFN)
and the Dutch Nikhef, with contributions by institutions from Belgium,
Germany, Greece, Hungary, Ireland, Japan, Monaco, Poland, Portugal, Spain.
KAGRA is supported by Ministry of Education, Culture, Sports, Science and
Technology (MEXT), Japan Society for the Promotion of Science (JSPS) in Japan;
National Research Foundation (NRF) and Ministry of Science and ICT (MSIT) in
Korea; Academia Sinica (AS) and National Science and Technology Council (NSTC)
in Taiwan. The authors are grateful for computational resources provided by
the LIGO Laboratory and supported by NSF Grants PHY-0757058 and PHY-0823459.
Lastly, we acknowledge the efforts of the CE Consortium and ET Collaboration
in the planning and development of XG GW detectors. This manuscript carries an
LSC DCC number LIGO-P2300334.

## References

  * Aasi  _et al._ (2015) J. Aasi _et al._ (LIGO Scientific), “Advanced LIGO,” Class. Quant. Grav. 32, 074001 (2015), arXiv:1411.4547 [gr-qc] .
  * Abbott  _et al._ (2016) B. P. Abbott _et al._ (LIGO Scientific, Virgo), “Observation of Gravitational Waves from a Binary Black Hole Merger,” Phys. Rev. Lett. 116, 061102 (2016), arXiv:1602.03837 [gr-qc] .
  * Abbott  _et al._ (2017a) B. P. Abbott _et al._ (LIGO Scientific, Virgo), “GW170817: Observation of Gravitational Waves from a Binary Neutron Star Inspiral,” Phys. Rev. Lett. 119, 161101 (2017a), arXiv:1710.05832 [gr-qc] .
  * Abbott  _et al._ (2021a) R. Abbott _et al._ (LIGO Scientific, KAGRA, VIRGO), “Observation of Gravitational Waves from Two Neutron Star–Black Hole Coalescences,” Astrophys. J. Lett. 915, L5 (2021a), arXiv:2106.15163 [astro-ph.HE] .
  * 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] .
  * Abbott  _et al._ (2021b) 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 (2021b), arXiv:2010.14527 [gr-qc] .
  * Abbott  _et al._ (2021c) 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,” (2021c), arXiv:2111.03606 [gr-qc] .
  * Caprini and Figueroa (2018) Chiara Caprini and Daniel G. Figueroa, “Cosmological Backgrounds of Gravitational Waves,” Class. Quant. Grav. 35, 163001 (2018), arXiv:1801.04268 [astro-ph.CO] .
  * Christensen (2019) Nelson Christensen, “Stochastic Gravitational Wave Backgrounds,” Rept. Prog. Phys. 82, 016903 (2019), arXiv:1811.08797 [gr-qc] .
  * Sakellariadou (2022) Mairi Sakellariadou, “Gravitational Waves: The Theorist’s Swiss Knife,” Universe 8, 132 (2022), arXiv:2202.00735 [astro-ph.CO] .
  * Evans  _et al._ (2021) Matthew Evans _et al._ , “A Horizon Study for Cosmic Explorer: Science, Observatories, and Community,” (2021), arXiv:2109.09882 [astro-ph.IM] .
  * Punturo  _et al._ (2010) M. Punturo _et al._ , “The Einstein Telescope: A third-generation gravitational wave observatory,” Class. Quant. Grav. 27, 194002 (2010).
  * Regimbau  _et al._ (2017) T. Regimbau, M. Evans, N. Christensen, E. Katsavounidis, B. Sathyaprakash, and S. Vitale, “Digging deeper: Observing primordial gravitational waves below the binary black hole produced stochastic background,” Phys. Rev. Lett. 118, 151105 (2017), arXiv:1611.08943 [astro-ph.CO] .
  * Maggiore  _et al._ (2020) Michele Maggiore _et al._ , “Science Case for the Einstein Telescope,” JCAP 03, 050 (2020), arXiv:1912.02622 [astro-ph.CO] .
  * Regimbau (2011) Tania Regimbau, “The astrophysical gravitational wave stochastic background,” Res. Astron. Astrophys. 11, 369–390 (2011), arXiv:1101.2762 [astro-ph.CO] .
  * Regimbau and Mandic (2008) T. Regimbau and V. Mandic, “Astrophysical Sources of Stochastic Gravitational-Wave Background,” Class. Quant. Grav. 25, 184018 (2008), arXiv:0806.2794 [astro-ph] .
  * Buonanno  _et al._ (2005) Alessandra Buonanno, Gunter Sigl, Georg G. Raffelt, Hans-Thomas Janka, and Ewald Muller, “Stochastic gravitational wave background from cosmological supernovae,” Phys. Rev. D 72, 084001 (2005), arXiv:astro-ph/0412277 .
  * Kalogera  _et al._ (2019) Vassiliki Kalogera _et al._ , “Deeper, Wider, Sharper: Next-Generation Ground-Based Gravitational-Wave Observations of Binary Black Holes,” (2019), arXiv:1903.09220 [astro-ph.HE] .
  * Evans  _et al._ (2023) Matthew Evans _et al._ , “Cosmic Explorer: A Submission to the NSF MPSAC ngGW Subcommittee,” (2023), arXiv:2306.13745 [astro-ph.IM] .
  * Zhou  _et al._ (2023) Bei Zhou, Luca Reali, Emanuele Berti, Mesut Çalışkan, Cyril Creque-Sarbinowski, Marc Kamionkowski, and B. S. Sathyaprakash, “Subtracting compact binary foregrounds to search for subdominant gravitational-wave backgrounds in next-generation ground-based observatories,” Phys. Rev. D 108, 064040 (2023), arXiv:2209.01310 [gr-qc] .
  * Zhou  _et al._ (2022) Bei Zhou, Luca Reali, Emanuele Berti, Mesut Çalışkan, Cyril Creque-Sarbinowski, Marc Kamionkowski, and B. S. Sathyaprakash, “Compact Binary Foreground Subtraction in Next-Generation Ground-Based Observatories,” (2022), arXiv:2209.01221 [gr-qc] .
  * Zhong  _et al._ (2023) Haowen Zhong, Rich Ormiston, and Vuk Mandic, “Detecting cosmological gravitational wave background after removal of compact binary coalescences in future gravitational wave detectors,” Phys. Rev. D 107, 064048 (2023), arXiv:2209.11877 [gr-qc] .
  * Sachdev  _et al._ (2020) Surabhi Sachdev, Tania Regimbau, and B. S. Sathyaprakash, “Subtracting compact binary foreground sources to reveal primordial gravitational-wave backgrounds,” Phys. Rev. D 102, 024051 (2020), arXiv:2002.05365 [gr-qc] .
  * Regimbau (2022) Tania Regimbau, “The Quest for the Astrophysical Gravitational-Wave Background with Terrestrial Detectors,” Symmetry 14, 270 (2022).
  * Sharma and Harms (2020) Ashish Sharma and Jan Harms, “Searching for cosmological gravitational-wave backgrounds with third-generation detectors in the presence of an astrophysical foreground,” Phys. Rev. D 102, 063009 (2020), arXiv:2006.16116 [gr-qc] .
  * Abbott  _et al._ (2017b) Benjamin P. Abbott _et al._ (LIGO Scientific, Virgo), “Upper Limits on the Stochastic Gravitational-Wave Background from Advanced LIGO’s First Observing Run,” Phys. Rev. Lett. 118, 121101 (2017b), [Erratum: Phys.Rev.Lett. 119, 029901 (2017)], arXiv:1612.02029 [gr-qc] .
  * Abbott  _et al._ (2019b) B. P. Abbott _et al._ (LIGO Scientific, Virgo), “Search for the isotropic stochastic background using data from Advanced LIGO’s second observing run,” Phys. Rev. D 100, 061101 (2019b), arXiv:1903.02886 [gr-qc] .
  * Abbott  _et al._ (2021d) R. Abbott _et al._ (KAGRA, Virgo, LIGO Scientific), “Upper limits on the isotropic gravitational-wave background from Advanced LIGO and Advanced Virgo’s third observing run,” Phys. Rev. D 104, 022004 (2021d), arXiv:2101.12130 [gr-qc] .
  * Aghanim  _et al._ (2020) N. Aghanim _et al._ (Planck), “Planck 2018 results. VI. Cosmological parameters,” Astron. Astrophys. 641, A6 (2020), [Erratum: Astron.Astrophys. 652, C4 (2021)], arXiv:1807.06209 [astro-ph.CO] .
  * Wu  _et al._ (2012) C. Wu, V. Mandic, and T. Regimbau, “Accessibility of the Gravitational-Wave Background due to Binary Coalescences to Second and Third Generation Gravitational-Wave Detectors,” Phys. Rev. D 85, 104024 (2012), arXiv:1112.1898 [gr-qc] .
  * Zhu  _et al._ (2013) Xing-Jiang Zhu, Eric J. Howell, David G. Blair, and Zong-Hong Zhu, “On the gravitational wave background from compact binary coalescences in the band of ground-based interferometers,” Mon. Not. Roy. Astron. Soc. 431, 882–899 (2013), arXiv:1209.0595 [gr-qc] .
  * Romano and Cornish (2017) Joseph D. Romano and Neil J. Cornish, “Detection methods for stochastic gravitational-wave backgrounds: a unified treatment,” Living Rev. Rel. 20, 2 (2017), arXiv:1608.06889 [gr-qc] .
  * Ferrari  _et al._ (1999) Valeria Ferrari, Sabino Matarrese, and Raffaella Schneider, “Stochastic background of gravitational waves generated by a cosmological population of young, rapidly rotating neutron stars,” Mon. Not. Roy. Astron. Soc. 303, 258 (1999), arXiv:astro-ph/9806357 .
  * Dominik  _et al._ (2015) Michal Dominik, Emanuele Berti, Richard O’Shaughnessy, Ilya Mandel, Krzysztof Belczynski, Christopher Fryer, Daniel E. Holz, Tomasz Bulik, and Francesco Pannarale, “Double Compact Objects III: Gravitational Wave Detection Rates,” Astrophys. J. 806, 263 (2015), arXiv:1405.7016 [astro-ph.HE] .
  * Phinney (2001) E. S. Phinney, “A Practical theorem on gravitational wave backgrounds,” (2001), arXiv:astro-ph/0108028 .
  * Safarzadeh  _et al._ (2020) Mohammadtaher Safarzadeh, Sylvia Biscoveanu, and Abraham Loeb, “Constraining the delay time distribution of compact binary objects from the stochastic gravitational wave background searches,” Astrophys. J. 901, 137 (2020), arXiv:2004.12999 [astro-ph.HE] .
  * Hogg (1999) David W. Hogg, “Distance measures in cosmology,” (1999), arXiv:astro-ph/9905116 .
  * Piran (1992) Tsvi Piran, “The Implications of the Compton (GRO) Observations for Cosmological Gamma-Ray Bursts,” Astrophys. J. Lett. 389, L45 (1992).
  * Sana  _et al._ (2012) H. Sana, S. E. de Mink, A. de Koter, N. Langer, C. J. Evans, M. Gieles, E. Gosset, R. G. Izzard, J. B. Le Bouquin, and F. R. N. Schneider, “Binary interaction dominates the evolution of massive stars,” Science 337, 444 (2012), arXiv:1207.6397 [astro-ph.SR] .
  * Fishbach and Kalogera (2021) Maya Fishbach and Vicky Kalogera, “The Time Delay Distribution and Formation Metallicity of LIGO-Virgo’s Binary Black Holes,” Astrophys. J. Lett. 914, L30 (2021), arXiv:2105.06491 [astro-ph.HE] .
  * Beniamini and Piran (2019) Paz Beniamini and Tsvi Piran, “The Gravitational waves merger time distribution of binary neutron star systems,” Mon. Not. Roy. Astron. Soc. 487, 4847–4854 (2019), arXiv:1903.11614 [astro-ph.HE] .
  * Zevin  _et al._ (2022) Michael Zevin, Anya E. Nugent, Susmita Adhikari, Wen-fai Fong, Daniel E. Holz, and Luke Zoltan Kelley, “Observational Inference on the Delay Time Distribution of Short Gamma-Ray Bursts,” Astrophys. J. Lett. 940, L18 (2022), arXiv:2206.02814 [astro-ph.HE] .
  * Regimbau  _et al._ (2012) Tania Regimbau _et al._ , “A Mock Data Challenge for the Einstein Gravitational-Wave Telescope,” Phys. Rev. D 86, 122001 (2012), arXiv:1201.3563 [gr-qc] .
  * Belczynski and Kalogera (2001) Krzysztof Belczynski and Vassiliki Kalogera, “A New formation channel for double neutron stars without recycling: Implications for gravitational wave detection,” Astrophys. J. Lett. 550, L183–L187 (2001), arXiv:astro-ph/0012172 .
  * Belczynski  _et al._ (2006) Krzysztof Belczynski, Rosalba Perna, Tomasz Bulik, Vassiliki Kalogera, Natalia Ivanova, and Donald Q. Lamb, “A study of compact object mergers as short gamma-ray burst progenitors,” Astrophys. J. 648, 1110–1116 (2006), arXiv:astro-ph/0601458 .
  * Berger  _et al._ (2007) Edo Berger _et al._ , “A New Population of High Redshift Short-Duration Gamma-Ray Bursts,” Astrophys. J. 664, 1000–1010 (2007), arXiv:astro-ph/0611128 .
  * Callister  _et al._ (2020) Thomas Callister, Maya Fishbach, Daniel Holz, and Will Farr, “Shouts and Murmurs: Combining Individual Gravitational-Wave Sources with the Stochastic Background to Measure the History of Binary Black Hole Mergers,” Astrophys. J. Lett. 896, L32 (2020), arXiv:2003.12152 [astro-ph.HE] .
  * Madau and Fragos (2017) Piero Madau and Tassos Fragos, “Radiation Backgrounds at Cosmic Dawn: X-Rays from Compact Binaries,” Astrophys. J. 840, 39 (2017), arXiv:1606.07887 [astro-ph.GA] .
  * (50) _LIGO Document T2000100-v6_ , Tech. Rep.
  * 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] .
  * Farah  _et al._ (2022) Amanda M. Farah, Maya Fishbach, Reed Essick, Daniel E. Holz, and Shanika Galaudage, “Bridging the Gap: Categorizing Gravitational-wave Events at the Transition between Neutron Stars and Black Holes,” Astrophys. J. 931, 108 (2022), arXiv:2111.03498 [astro-ph.HE] .
  * Fishbach  _et al._ (2020a) Maya Fishbach, Reed Essick, and Daniel E. Holz, “Does Matter Matter? Using the mass distribution to distinguish neutron stars and black holes,” Astrophys. J. Lett. 899, L8 (2020a), arXiv:2006.13178 [astro-ph.HE] .
  * Biscoveanu  _et al._ (2022) Sylvia Biscoveanu, Philippe Landry, and Salvatore Vitale, “Population properties and multimessenger prospects of neutron star–black hole mergers following GWTC-3,” Mon. Not. Roy. Astron. Soc. 518, 5298–5312 (2022), arXiv:2207.01568 [astro-ph.HE] .
  * Farrow  _et al._ (2019) Nicholas Farrow, Xing-Jiang Zhu, and Eric Thrane, “The mass distribution of Galactic double neutron stars,” Astrophys. J. 876, 18 (2019), arXiv:1902.03300 [astro-ph.HE] .
  * Landry and Read (2021) Philippe Landry and Jocelyn S. Read, “The Mass Distribution of Neutron Stars in Gravitational-wave Binaries,” Astrophys. J. Lett. 921, L25 (2021), arXiv:2107.04559 [astro-ph.HE] .
  * Fishbach and Holz (2017) Maya Fishbach and Daniel E. Holz, “Where Are LIGO’s Big Black Holes?” Astrophys. J. Lett. 851, L25 (2017), arXiv:1709.08584 [astro-ph.HE] .
  * Fishbach and Holz (2020) Maya Fishbach and Daniel E. Holz, “Picky Partners: The Pairing of Component Masses in Binary Black Hole Mergers,” Astrophys. J. Lett. 891, L27 (2020), arXiv:1905.12669 [astro-ph.HE] .
  * Abbott  _et al._ (2020) B. P. Abbott _et al._ (LIGO Scientific, Virgo), “GW190425: Observation of a Compact Binary Coalescence with Total Mass $\sim 3.4M_{\odot}$,” Astrophys. J. Lett. 892, L3 (2020), arXiv:2001.01761 [astro-ph.HE] .
  * Romero-Shaw  _et al._ (2020) Isobel M. Romero-Shaw, Nicholas Farrow, Simon Stevenson, Eric Thrane, and Xing-Jiang Zhu, “On the origin of GW190425,” Mon. Not. Roy. Astron. Soc. 496, L64–L69 (2020), arXiv:2001.06492 [astro-ph.HE] .
  * Gupta  _et al._ (2023) Ish Gupta _et al._ , “Characterizing Gravitational Wave Detector Networks: From A♯ to Cosmic Explorer,” (2023), arXiv:2307.10421 [gr-qc] .
  * Borhanian (2021) Ssohrab Borhanian, “GWBENCH: a novel Fisher information package for gravitational-wave benchmarking,” Class. Quant. Grav. 38, 175014 (2021), arXiv:2010.15202 [gr-qc] .
  * Fritschel  _et al._ (2022) P. Fritschel, K. Kuns, J. Driggers, A. Effler, B. Lantz, D. Ottaway, S. Ballmer, K. Dooley, R. Adhikari, M. Evans, B. Farr, G. Gonzalez, P. Schmidt, and S. Raja,  _Report from the LSC Post-O5 Study Group_ , Tech. Rep. T2200287 (LIGO, 2022).
  * (64) “Et-gw tds,” .
  * Allen and Romano (1999) Bruce Allen and Joseph D. Romano, “Detecting a stochastic background of gravitational radiation: Signal processing strategies and sensitivities,” Phys. Rev. D 59, 102001 (1999), arXiv:gr-qc/9710117 .
  * Thrane and Romano (2013) Eric Thrane and Joseph D. Romano, “Sensitivity curves for searches for gravitational-wave backgrounds,” Phys. Rev. D 88, 124032 (2013), arXiv:1310.5300 [astro-ph.IM] .
  * Regimbau  _et al._ (2014) Tania Regimbau, Duncan Meacher, and Michael Coughlin, “Second Einstein Telescope mock science challenge: Detection of the gravitational-wave stochastic background from compact binary coalescences,” Phys. Rev. D 89, 084046 (2014), arXiv:1404.1134 [astro-ph.CO] .
  * Meacher  _et al._ (2015) Duncan Meacher, Michael Coughlin, Sean Morris, Tania Regimbau, Nelson Christensen, Shivaraj Kandhasamy, Vuk Mandic, Joseph D. Romano, and Eric Thrane, “Mock data and science challenge for detecting an astrophysical stochastic gravitational-wave background with Advanced LIGO and Advanced Virgo,” Phys. Rev. D 92, 063002 (2015), arXiv:1506.06744 [astro-ph.HE] .
  * Husa  _et al._ (2016) Sascha Husa, Sebastian Khan, Mark Hannam, Michael Pürrer, Frank Ohme, Xisco Jiménez Forteza, and Alejandro Bohé, “Frequency-domain gravitational waves from nonprecessing black-hole binaries. I. New numerical waveforms and anatomy of the signal,” Phys. Rev. D 93, 044006 (2016), arXiv:1508.07250 [gr-qc] .
  * Khan  _et al._ (2016) Sebastian Khan, Sascha Husa, Mark Hannam, Frank Ohme, Michael Pürrer, Xisco Jiménez Forteza, and Alejandro Bohé, “Frequency-domain gravitational waves from nonprecessing black-hole binaries. II. A phenomenological model for the advanced detector era,” Phys. Rev. D 93, 044007 (2016), arXiv:1508.07253 [gr-qc] .
  * Moore  _et al._ (2015) C. J. Moore, R. H. Cole, and C. P. L. Berry, “Gravitational-wave sensitivity curves,” Class. Quant. Grav. 32, 015014 (2015), arXiv:1408.0740 [gr-qc] .
  * Fishbach  _et al._ (2020b) Maya Fishbach, Will M. Farr, and Daniel E. Holz, “The Most Massive Binary Black Hole Detections and the Identification of Population Outliers,” Astrophys. J. Lett. 891, L31 (2020b), arXiv:1911.05882 [astro-ph.HE] .
  * Fishbach  _et al._ (2018) Maya Fishbach, Daniel E. Holz, and Will M. Farr, “Does the Black Hole Merger Rate Evolve with Redshift?” Astrophys. J. Lett. 863, L41 (2018), arXiv:1805.10270 [astro-ph.HE] .
  * Saikawa (2017) Ken’ichi Saikawa, “A review of gravitational waves from cosmic domain walls,” Universe 3, 40 (2017), arXiv:1703.02576 [hep-ph] .
  * Romero  _et al._ (2021) Alba Romero, Katarina Martinovic, Thomas A. Callister, Huai-Ke Guo, Mario Martínez, Mairi Sakellariadou, Feng-Wei Yang, and Yue Zhao, “Implications for First-Order Cosmological Phase Transitions from the Third LIGO-Virgo Observing Run,” Phys. Rev. Lett. 126, 151301 (2021), arXiv:2102.01714 [hep-ph] .
  * Espinosa  _et al._ (2010) Jose R. Espinosa, Thomas Konstandin, Jose M. No, and Geraldine Servant, “Energy Budget of Cosmological First-order Phase Transitions,” JCAP 06, 028 (2010), arXiv:1004.4187 [hep-ph] .
  * Figueroa and Tanin (2019) Daniel G. Figueroa and Erwin H. Tanin, “Ability of LIGO and LISA to probe the equation of state of the early Universe,” JCAP 08, 011 (2019), arXiv:1905.11960 [astro-ph.CO] .
  * Easther  _et al._ (2007) Richard Easther, John T. Giblin, Jr., and Eugene A. Lim, “Gravitational Wave Production At The End Of Inflation,” Phys. Rev. Lett. 99, 221301 (2007), arXiv:astro-ph/0612294 .
  * Abbott  _et al._ (2021e) R. Abbott _et al._ (LIGO Scientific, Virgo, KAGRA), “Constraints on Cosmic Strings Using Data from the Third Advanced LIGO–Virgo Observing Run,” Phys. Rev. Lett. 126, 241102 (2021e), arXiv:2101.12248 [gr-qc] .
  * Jenkins  _et al._ (2019) Alexander C. Jenkins, Richard O’Shaughnessy, Mairi Sakellariadou, and Daniel Wysocki, “Anisotropies in the astrophysical gravitational-wave background: The impact of black hole distributions,” Phys. Rev. Lett. 122, 111101 (2019), arXiv:1810.13435 [astro-ph.CO] .
  * Vangioni  _et al._ (2015) E. Vangioni, K. A. Olive, T. Prestegard, J. Silk, P. Petitjean, and V. Mandic, “The Impact of Star Formation and Gamma-Ray Burst Rates at High Redshift on Cosmic Chemical Evolution and Reionization,” Mon. Not. Roy. Astron. Soc. 447, 2575 (2015), arXiv:1409.2462 [astro-ph.GA] .
  * Essick  _et al._ (2022) Reed Essick, Amanda Farah, Shanika Galaudage, Colm Talbot, Maya Fishbach, Eric Thrane, and Daniel E. Holz, “Probing Extremal Gravitational-wave Events with Coarse-grained Likelihoods,” Astrophys. J. 926, 34 (2022), arXiv:2109.00418 [astro-ph.HE] .
  * Martinovic  _et al._ (2021) Katarina Martinovic, Patrick M. Meyers, Mairi Sakellariadou, and Nelson Christensen, “Simultaneous estimation of astrophysical and cosmological stochastic gravitational-wave backgrounds with terrestrial detectors,” Phys. Rev. D 103, 043023 (2021), arXiv:2011.05697 [gr-qc] .
  * Yu and Wang (2022) Wang-Wei Yu and Shao-Jiang Wang, “Searching for double-peak and doubly-broken gravitational-wave spectra from Advanced LIGO-Virgo’s first three observing runs,” (2022), arXiv:2211.13111 [gr-qc] .
  * Kaiser  _et al._ (2022) Andrew R. Kaiser, Nihan S. Pol, Maura A. McLaughlin, Siyuan Chen, Jeffrey S. Hazboun, Luke Zoltan Kelley, Joseph Simon, Stephen R. Taylor, Sarah J. Vigeland, and Caitlin A. Witt, “Disentangling Multiple Stochastic Gravitational Wave Background Sources in PTA Data Sets,” Astrophys. J. 938, 115 (2022), arXiv:2208.02307 [astro-ph.CO] .
  * Callister  _et al._ (2016) Thomas Callister, Letizia Sammut, Shi Qiu, Ilya Mandel, and Eric Thrane, “The limits of astrophysics with gravitational-wave backgrounds,” Phys. Rev. X 6, 031018 (2016), arXiv:1604.02513 [gr-qc] .
  * Mukherjee  _et al._ (2022) Suvodip Mukherjee, Matthew S. P. Meinema, and Joseph Silk, “Prospects of discovering subsolar primordial black holes using the stochastic gravitational wave background from third-generation detectors,” Mon. Not. Roy. Astron. Soc. 510, 6218–6224 (2022), arXiv:2107.02181 [astro-ph.CO] .
  * Smith and Thrane (2018) Rory Smith and Eric Thrane, “Optimal Search for an Astrophysical Gravitational-Wave Background,” Phys. Rev. X 8, 021019 (2018), arXiv:1712.00688 [gr-qc] .
  * Smith  _et al._ (2020) Rory J. E. Smith, Colm Talbot, Francisco Hernandez Vivanco, and Eric Thrane, “Inferring the population properties of binary black holes from unresolved gravitational waves,” Mon. Not. Roy. Astron. Soc. 496, 3281–3290 (2020), arXiv:2004.09700 [astro-ph.HE] .
  * Talbot  _et al._ (2021) Colm Talbot, Eric Thrane, Sylvia Biscoveanu, and Rory Smith, “Inference with finite time series: Observing the gravitational Universe through windows,” Phys. Rev. Res. 3, 043049 (2021), arXiv:2106.13785 [astro-ph.IM] .
  * Banagiri  _et al._ (2020) Sharan Banagiri, Vuk Mandic, Claudia Scarlata, and Kate Z. Yang, “Measuring angular N-point correlations of binary black hole merger gravitational-wave events with hierarchical Bayesian inference,” Phys. Rev. D 102, 063007 (2020), arXiv:2006.00633 [astro-ph.CO] .
  * Biscoveanu  _et al._ (2020) Sylvia Biscoveanu, Colm Talbot, Eric Thrane, and Rory Smith, “Measuring the primordial gravitational-wave background in the presence of astrophysical foregrounds,” Phys. Rev. Lett. 125, 241101 (2020), arXiv:2009.04418 [astro-ph.HE] .
  * Finkelstein  _et al._ (2023) Steven L. Finkelstein, Micaela B. Bagley, Henry C. Ferguson, Stephen M. Wilkins, Jeyhan S. Kartaltepe, Casey Papovich, L. Y. Aaron Yung, Pablo Arrabal Haro, Peter Behroozi, Mark Dickinson, Dale D. Kocevski, Anton M. Koekemoer, Rebecca L. Larson, Auré lien Le Bail, Alexa M. Morales, Pablo G. Pérez-González, Denis Burgarella, Romeel Davé, Michaela Hirschmann, Rachel S. Somerville, Stijn Wuyts, Volker Bromm, Caitlin M. Casey, Adriano Fontana, Seiji Fujimoto, Jonathan P. Gardner, Mauro Giavalisco, Andrea Grazian, Norman A. Grogin, Nimish P. Hathi, Taylor A. Hutchison, Saurabh W. Jha, Shardha Jogee, Lisa J. Kewley, Allison Kirkpatrick, Arianna S. Long, Jennifer M. Lotz, Laura Pentericci, Justin D. R. Pierel, Nor Pirzkal, Swara Ravindranath, Russell E. Ryan, Jonathan R. Trump, Guang Yang, Rachana Bhatawdekar, Laura Bisigello, Véronique Buat, Antonello Calabrò, Marco Castellano, Nikko J. Cleri, M. C. Cooper, Darren Croton, Emanuele Daddi, Avishai Dekel, David Elbaz, Maximilien Franco, Eric Gawiser, Benne W. Holwerda, Marc Huertas-Company, Anne E. Jaskot, Gene C. K. Leung, Ray A. Lucas, Bahram Mobasher, Viraj Pandya, Sandro Tacchella, Benjamin J. Weiner, and Jorge A. Zavala, “CEERS key paper. i. an early look into the first 500 myr of galaxy formation with JWST,” The Astrophysical Journal Letters 946, L13 (2023).
  * Harikane  _et al._ (2023) Yuichi Harikane, Masami Ouchi, Masamune Oguri, Yoshiaki Ono, Kimihiko Nakajima, Yuki Isobe, Hiroya Umeda, Ken Mawatari, and Yechi Zhang, “A comprehensive study of galaxies at z $\sim$ 9–16 found in the early JWST data: Ultraviolet luminosity functions and cosmic star formation history at the pre-reionization epoch,” The Astrophysical Journal Supplement Series 265, 5 (2023).
  * Lehoucq  _et al._ (2023) Leonard Lehoucq, Irina Dvorkin, Rahul Srinivasan, Clement Pellouin, and Astrid Lamberts, “Astrophysical Uncertainties in the Gravitational-Wave Background from Stellar-Mass Compact Binary Mergers,” (2023), arXiv:2306.09861 [astro-ph.HE] .
  * Périgois  _et al._ (2021) C. Périgois, C. Belczynski, T. Bulik, and T. Regimbau, “StarTrack predictions of the stochastic gravitational-wave background from compact binary mergers,” Phys. Rev. D 103, 043002 (2021), arXiv:2008.04890 [astro-ph.CO] .
  * Kowalska-Leszczynska  _et al._ (2015) I. Kowalska-Leszczynska, T. Regimbau, T. Bulik, M. Dominik, and K. Belczynski, “Effect of metallicity on the gravitational-wave signal from the cosmological population of compact binary coalescences,” Astron. Astrophys. 574, A58 (2015), arXiv:1205.4621 [astro-ph.CO] .
  * Martinovic  _et al._ (2022) Katarina Martinovic, Carole Perigois, Tania Regimbau, and Mairi Sakellariadou, “Footprints of Population III Stars in the Gravitational-wave Background,” Astrophys. J. 940, 29 (2022), arXiv:2109.09779 [astro-ph.SR] .
  * Nakazato  _et al._ (2016) Ken Nakazato, Yuu Niino, and Norichika Sago, “Gravitational-Wave Background from Binary Mergers and Metallicity Evolution of Galaxies,” Astrophys. J. 832, 146 (2016), arXiv:1605.02146 [astro-ph.HE] .
  * Dvorkin  _et al._ (2016) Irina Dvorkin, Elisabeth Vangioni, Joseph Silk, Jean-Philippe Uzan, and Keith A. Olive, “Metallicity-constrained merger rates of binary black holes and the stochastic gravitational wave background,” Mon. Not. Roy. Astron. Soc. 461, 3877–3885 (2016), arXiv:1604.04288 [astro-ph.HE] .

##  Appendix A Stochastic backgrounds from alternate detector networks

In this appendix, we show the results for the alternate detector networks that
we consider in this study. Figure 7 shows our $\Omega_{\text{cbc, unres}}$
estimates for the pdb-pl and joint-envelope models for the various alternate
detector networks we consider in Sec. V. Also overlaid are the different
cosmological models discussed in Sec. VII. The equivalent plot for our
fiducial network is Fig. 3. Finally, Fig. 8 overlays the pdb-pl model
$\Omega_{\text{cbc, unres}}$ estimates for all our considered networks as a
means of direct comparison.

(a) $\Omega_{\text{cbc, unres}}$ estimates for a network consisting of a 40-km
CE A and a 20-km CE B.

(b) $\Omega_{\text{cbc, unres}}$ estimates for a network consisting of a 40-km
CE A and an ET.

(c) $\Omega_{\text{cbc, unres}}$ estimates for a network consisting of a 40-km
CE A, an $\text{A}^{\\#}$ at LLO, and an ET.

Figure 7: A comparison of various cosmological SGWB models along with
$\Omega_{\text{cbc, unres}}$ estimates for the various alternate detector
networks we consider in Sec. V. The blue band shows the pdb-pl model while the
pink band shows the joint-envelope constructed from the various BNS-only and
NSBH-only models in Tab. 1. Both bands represent $90\%$ credible intervals.
The solid black curve is the $3\sigma$ PI curve. Figure 8: A comparision of
the 90% credible intervals of $\Omega_{\text{cbc, unres}}$ for the pdb-pl
model for the four different detector networks defined in Sec. V.

  *[GW]: gravitational wave
  *[GWs]: gravitational wave
  *[SGWB]: stochastic gravitational wave background
  *[CE]: Cosmic Explorer
  *[ET]: Einstein Telescope
  *[SGWBs]: stochastic gravitational wave background