New binary black hole mergers in the second observing run of Advanced LIGO and Advanced Virgo

TL;DR

This paper reports the discovery of six new binary black hole merger events in LIGO and Virgo O2 data, nearly doubling the known sample and providing new insights into black hole populations and their properties.

Contribution

The paper introduces a new search pipeline that identified six additional black hole mergers in existing data, expanding the gravitational wave event catalog.

Findings

Six new black hole merger events detected.

One event (GW170121) has a negative effective spin parameter.

Detection volume of the pipeline may be twice that of previous methods.

Abstract

We report the detection of new binary black hole merger events in the publicly available data from the second observing run of advanced LIGO and advanced Virgo (O2). The mergers were discovered using the new search pipeline described in Venumadhav et al. [Phys. Rev. D 100, 023011 (2019)], and are above the detection thresholds as defined in Abbott et al. (LIGO Scientific and Virgo Collaborations) [Phys. Rev. X 9, 031040 (2019).] Three of the mergers (GW170121, GW170304, GW170727) have inferred probabilities of being of astrophysical origin $p_{\rm astro} > 0.98$. The remaining three (GW170425, GW170202, GW170403) are less certain, with $p_{\rm astro}$ ranging from 0.5 to 0.8. The newly found mergers largely share the statistical properties of previously reported events, with the exception of GW170403, the least secure event, which has a highly negative effective spin parameter $\chi_{\rm eff}$ . The most secure new event, GW170121 ($p_{\rm astro} > 0.99$), is also notable due to its inferred negative value of $\chi_{\rm eff}$, which is inconsistent with being positive at the ~95.8% confidence level. The new mergers nearly double the sample of gravitational wave events reported from O2, and present a substantial opportunity to explore the statistics of the binary black hole population in the Universe. The number of detected events is not surprising since we estimate that the detection volume of our pipeline may be larger than that of other pipelines by as much as a factor of two (with significant uncertainties in the estimate). The increase in volume is larger when the constituent detectors of the network have very different sensitivities, as is likely to be the case in current and future runs.

Figures

Tables

Table 1. Table 1 : Events already reported by the LIGO–Virgo Collaboration Abbott  et al. ( 2018 ) as detected with our pipeline. The rate distributions used to compute p astro subscript 𝑝 astro p_{\rm astro} are shown in Fig.  4 . The maximum likelihood rates are ℛ max = 8 / O2 subscript ℛ max 8 O2 \mathcal{R}_{\rm max}=8/{\rm O2} and 5 / O2 5 O2 5/{\rm O2} in banks BBH 3 and BBH 4 , respectively.

Name	Bank	GPS time111The times given are the ‘linear-free’ times of the best fit templates in our bank in the Hanford detector; with this time as the origin, the phase of the template is orthogonal to shifts in time, given the fiducial PSD.	${\rho }_H^2$	${\rho }_L^2$	IFAR (O2)222The inverse false alarm rates, or IFARs, are computed within each bank, and do not include any additional trials factors; our BBH analysis has 5 chirp-mass banks. The IFAR is given in terms of “O2” instead of physical time, since the ranking statistic includes a time-dependent volumetric correction factor to account for the significant and systematic changes in the network’s sensitivity over the run. If the network’s sensitivity were constant during the observing run, the unit $``\mathrm{O2}"\approx 118$ days.	$\frac{W\left(\text{event}\right)}{ℛ\left(\text{event}\mid 𝒩\right)}$ (O2)	$p_{\mathrm{astro}}$
GW170104	BBH (3,0)	$1167559936.582$	$85.1$	$104.3$	$>2\text{×}{10}^4$	$>100$	$>0.99$
GW170809	BBH (3,0)	$1186302519.740$	$40.5$	$113$	$>2\text{×}{10}^4$	$>100$	$>0.99$
GW170814	BBH (3,0)	$1186741861.519$	$90.2$	$170$	$>2\text{×}{10}^4$	$>100$	$>0.99$
GW170818	BBH (3,0)	$1187058327.075$	$19.4$	$95.1$	1.7333See discussion in Sec. IV.	—	—3
GW170729	BBH (3,1)	$1185389807.311$	$62.1$	$53.6$	$>2\text{×}{10}^4$	$>100$	$>0.99$
GW170823	BBH (3,1)	$1187529256.500$	$46.0$	$90.7$	$>2\text{×}{10}^4$	$>100$	$>0.99$

Table 2. Table 2 : New events with astrophysical probability > 50 % absent percent 50 >50\% in all of the BBH banks. The rate distributions used to compute p astro subscript 𝑝 astro p_{\rm astro} are shown in Fig.  4 , the maximum-likelihood rates in banks BBH 3 and BBH 4 are ℛ max = 8 / O2 subscript ℛ max 8 O2 \mathcal{R}_{\rm max}=8/{\rm O2} and 5 / O2 5 O2 5/{\rm O2} , respectively.

Name	Bank	${ℳ}^{\det }\left(M_{\odot }\right)$	${\chi }_{\mathrm{eff}}$	$z$	GPS time444The times given are the ‘linear-free’ times of the best fit templates in our bank in the Hanford detector; with this time as the origin, the phase of the template is orthogonal to shifts in time, given the fiducial PSD.	${\rho }_H^2$	${\rho }_L^2$	IFAR (O2)555The inverse false alarm rates, or IFARs, are computed within each bank, and do not include any additional trials factors; our BBH analysis has 5 chirp-mass banks. The IFAR is given in terms of “O2” instead of physical time, since the ranking statistic includes a time-dependent volumetric correction factor to account for the significant and systematic changes in the network’s sensitivity over the run. If the network’s sensitivity were constant during the observing run, the unit $``\mathrm{O2}"\approx 118$ days.	$\frac{W\left(\text{event}\right)}{ℛ\left(\text{event}\mid 𝒩\right)}$ (O2)	$p_{\mathrm{astro}}$
GW170121	BBH (3,0)	${29}_{-3}^{+4}$	$-{0.3}_{-0.3}^{+0.3}$	${0.24}_{-0.13}^{+0.14}$	$1169069154.565$	$29.4$	$89.7$	$2.8\times {10}^3$	$>30$	$>0.99$
GW170304	BBH (4,0)	${47}_{-7}^{+8}$	${0.2}_{-0.3}^{+0.3}$	${0.5}_{-0.2}^{+0.2}$	$1172680691.356$	$24.9$	$55.9$	377	$13.6$	$0.985$
GW170727	BBH (4,0)	${42}_{-6}^{+6}$	$-{0.1}_{-0.3}^{+0.3}$	${0.43}_{-0.17}^{+0.18}$	$1185152688.019$	$25.4$	$53.5$	370	$11.8$	$0.98$
GW170425	BBH (4,0)	${47}_{-10}^{+26}$	${0.0}_{-0.5}^{+0.4}$	${0.5}_{-0.3}^{+0.4}$	$1177134832.178$	$28.6$	$37.5$	15	$0.65$	$0.77$
GW170202	BBH (3,0)	${21.6}_{-1.4}^{+4.2}$	$-{0.2}_{-0.3}^{+0.4}$	${0.27}_{-0.12}^{+0.13}$	$1170079035.715$	$26.5$	$41.7$	6.3	$0.25$	$0.68$
GW170403	BBH (4,1)	${48}_{-7}^{+9}$	$-{0.7}_{-0.3}^{+0.5}$	${0.45}_{-0.19}^{+0.22}$	$1175295989.221$	$31.3$	$31.0$	4.7	$0.23$	$0.56$

Table 3. Table 3 : Sub-threshold candidates with astrophysical probability above 10% in all of the BBH banks. The rate distributions used to compute p astro subscript 𝑝 astro p_{\rm astro} are shown in Fig.  4 , the maximum-likelihood rates in banks BBH 3 and BBH 4 are ℛ max = 8 / O2 subscript ℛ max 8 O2 \mathcal{R}_{\rm max}=8/{\rm O2} and 5 / O2 5 O2 5/{\rm O2} , respectively.

Bank	GPS time666The times given are the ‘linear-free’ times of the best fit templates in our bank in the Hanford detector; with this time as the origin, the phase of the template is orthogonal to shifts in time, given the fiducial PSD.	${\rho }_H^2$	${\rho }_L^2$	IFAR (O2)777The inverse false alarm rates, or IFARs, are computed within each bank, and do not include any additional trials factors; our BBH analysis has 5 chirp-mass banks. The IFAR is given in terms of “O2” instead of physical time, since the ranking statistic includes a time-dependent volumetric correction factor to account for the significant and systematic changes in the network’s sensitivity over the run. If the network’s sensitivity were constant during the observing run, the unit $``\mathrm{O2}"\approx 118$ days.	$\frac{W\left(\text{event}\right)}{ℛ\left(\text{event}\mid 𝒩\right)}$ (O2)	$p_{\mathrm{astro}}$
BBH (4,1)	$1172487817.477$	$48.6$	$19.1$	$0.82$	$0.147$	$0.45$
BBH (3,0)	$1170914187.455$	$20.4$	$41.4$	$0.43$	$0.044$	$0.28$
BBH (3,1)	$1172449151.468$	$29.5$	$32.4$	$0.31$	$0.025$	$0.18$
BBH (4,0)	$1174138338.385$	$37.1$	$28.4$	$0.62$	$0.034$	$0.17$
BBH (3,0)	$1171863216.108$	$46.5$	$21.6$	$0.27$	$0.016$	$0.125$
BBH (3,1)	$1187176593.222$	$20.3$	$42.0$	$0.2$	$0.014$	$0.12$
BBH (3,0)	$1182674889.044$	$34.1$	$28.7$	$0.23$	$0.016$	$0.12$
BBH (3,1)	$1171410777.200$	$40.8$	$21.0$	$0.18$	$0.014$	$0.11$

Equations

$$p_{\mathrm{astro}}(\mathrm{event})$$

$$=\frac{\mathcal{R}(\mathrm{event}\mid \mathcal{S})}{\mathcal{R}(\mathrm{event}\mid \mathcal{N})+\mathcal{R}(\mathrm{event}\mid \mathcal{S})}.$$

$${\rho }_{\mathrm{th}\text{BBH 0}}^2$$

$${\rho }_{\mathrm{th}\text{BBH 1}}^2$$

$${\rho }_{\mathrm{th}\text{BBH 2}}^2$$

$${\rho }_{\mathrm{th}\text{BBH 3}}^2$$

$${\rho }_{\mathrm{th}\text{BBH 4}}^2$$

$$\mathcal{R}(\mathrm{event}\mid \mathcal{S})=W(\mathrm{event})\mathcal{R},$$

$$\mathcal{L}(\mathrm{data}\mid \mathcal{R})$$

$$\prod _{\mathrm{triggers}}[\mathcal{R}(\mathrm{event}\mid \mathcal{N})+W(\mathrm{event})\mathcal{R}],$$

$$p_{\mathrm{astro}}(\mathrm{event})={\int }_0^{\infty }{p}_{\mathrm{astro}}(\mathrm{event}\mid \mathcal{R})P(\mathcal{R})d\mathcal{R}.$$

$$\left[\begin{array}{cc}⟨f,f⟩& ⟨f,1⟩\\ ⟨f,1⟩& ⟨1,1⟩\end{array}\right]\left[\begin{array}{c}2\pi t_0\\ {\psi }_0\end{array}\right]$$

$$⟨a(f),b(f)⟩$$

$$\mathcal{L}(d\mid i,B,\mathcal{S})$$

$$\approx \frac{1}{(\Delta c_{\alpha ,i}{)}^{n_{c_{\alpha ,i}}}}\iint \mathrm{d}\mathbf{\Theta }\mathrm{d}Ap(A)p(\mathbf{\Theta })\mathcal{L}(d\mid A,\mathbf{\Theta }),$$

$$\mathcal{L}(d\mid A,\mathbf{\Theta })$$

$$\mathcal{L}(d\mid \mathbf{\Theta })\propto \exp \left[\frac{⟨d\mid t(\mathbf{\Theta }){⟩}^2}{2⟨t(\mathbf{\Theta })\mid t(\mathbf{\Theta })⟩}\right],$$

$$⟨A_{\ast }t(c_{\alpha })\mid A_{\ast }t(c_{\alpha ,0})⟩\approx {\rho }^2\left(1-\frac{\delta c_{\alpha }^2}{2}\right).$$

$$\mathcal{L}(d\mid i,B,\mathcal{S})$$

$$\propto \frac{(\rho /{\rho }_0{)}^{-4}}{N_{\mathrm{temp},B}(\Delta c_{\alpha ,i}{)}^{n_{c_{\alpha ,i}}}}\iint \mathrm{d}\delta \mathbf{\Theta }\exp \left[\frac{{\rho }^2}{2}(1-\delta {\Theta }^2)\right]$$

$$=\frac{(\rho /{\rho }_0{)}^{-4}\exp \left({\rho }^2/2\right)}{N_{\mathrm{temp},B}(\Delta c_{\alpha ,i}{)}^{n_{c_{\alpha ,i}}}}{\left(\frac{2\pi }{{\rho }^2}\right)}^{(n_{\mathbf{\Theta }}\equiv (n_{c_{\alpha ,i}}+2))/2}.$$

$$\begin{gathered} \mathcal{L}(d_{\mathrm{H}},d_{\mathrm{L}}\mid i,B,\mathcal{S}) \\ \propto \frac{p(\Delta t,\Delta \varphi ,{\rho }_{\mathrm{H}}^2,{\rho }_{\mathrm{L}}^2\mid {\rho }_{0,\mathrm{H}},{\rho }_{0,\mathrm{L}})e^{{\rho }^2/2}}{N_{\mathrm{temp},B}(\Delta c_{\alpha ,i}{)}^{n_{c_{\alpha ,i}}}}{\left(\frac{2\pi }{{\rho }^2}\right)}^{(n_{c_{\alpha ,i}}+4)/2}. \end{gathered}$$

$$P({\rho }_{\mathrm{H}}^2,{\rho }_{\mathrm{L}}^2,c_{\alpha },i)$$

$$=\frac{P({\rho }_{\mathrm{H}}^2\mid i)P({\rho }_{\mathrm{L}}^2\mid i)}{N_{\mathrm{temp},i}}P(i),$$

$$P(i)=\frac{N_{\mathrm{trig},i}}{N_{\mathrm{trig}}},$$

$$\frac{P({\rho }^2\mid i)}{P({\rho }_0^2\mid i)}\approx \frac{\mathrm{Rank}({\rho }^2\mid i)}{\mathrm{Rank}({\rho }_0^2\mid i)},$$

$$P({\rho }_0^2\mid i)$$

Full text

I Introduction

The LIGO–Virgo Collaboration (LVC) reported the detection of gravitational waves from ten binary black hole (BBH) and one binary neutron star mergers in their two latest observing runs, O1 and O2 Abbott et al. (2018, 2016a, 2016b, 2016c, 2017a, 2017b, 2017c, 2017d). The intrinsic properties of the mergers, namely the masses, the mass ratio and the spins of the black holes, are important observables that can inform us about how and where the binaries were assembled. All the BBH events are consistent with mergers of black holes with comparable masses (notably, GW170729 is mildly inconsistent with equal component masses Chatziioannou et al. (2019)); in two of the events, at least one of the components had a nonzero spin.

Currently, inference about the origin of the BBHs is limited by the small number of detected events. Future observations with improved sensitivity will enlarge the sample and map out the parameter space of the BBH population; meanwhile, it is important to adopt analysis techniques that maximize the yield of existing data.

A number of data analysis pipelines have been developed to search for transient events in LIGO–Virgo data. The two modeled searches used by the LVC are PyCBC Usman et al. (2016) and GstLAL Sachdev et al. (2019); these pipelines use matched filtering with a template bank of target compact binary coalescence signals. In addition, the LVC runs an unmodeled transient search with the coherent Wave Burst pipeline Klimenko et al. (2016). There are also groups external to the LVC running independent pipelines Nitz et al. (2018); Venumadhav et al. (2019b) on the public data released by the LVC, which now includes the O1 and O2 observing runs gwo (2 27); Vallisneri et al. (2015). In Ref. Venumadhav et al. (2019a), henceforth TV19, we presented a new and independent pipeline to analyze public LIGO data, which we applied to the public data from the O1 observing run. The cumulative improvements significantly increased the sensitive volume (at the same detection thresholds as those of Ref. Abbott et al. (2018)), and led to the detection of a new event, which is consistent with the merger of rapidly spinning and heavy black holes Zackay et al. (2019). In this paper, we present results from our search of coincident triggers in the Hanford (H1) and Livingston (L1) public data from the second observing run of advanced LIGO (O2) Vallisneri et al. (2015); gwo (2 27).

The paper is organized as follows: in Section II we begin with a preview of our final results to motivate what follows. We then discuss the details of the changes in the pipeline used in this paper compared to the one used to analyze the O1 data in Section III. Section IV summarizes the results of our analysis on the events that were previously detected by the LVC analysis pipelines in Ref. Abbott et al. (2018). We describe the new events we found in Section V, estimate our improvement in sensitivity in Section VI, and conclude with some remarks in Section VII. In Appendices A–C, we present the posterior distributions for the parameters of the new events, as well as some technical details.

II Preview of results

We improved upon the methods we previously developed to analyze the O1 data (see TV19 and Sec. III for details), and applied them to look for massive binary black hole mergers in the H1 and L1 data from the O2 data released by the LVC. We recover with high confidence all the BBH events that the LVC searches previously found in the bulk data release, with the notable exception of GW170818, which is one of our highest ranked candidates, but is effectively close to a single-detector trigger that would be hard to confirm using the methods adopted in this search (see Sec. IV for details).

This serves as an independent validation of the events presented in GWTC-1 Abbott et al. (2018), and a confirmation of the sensitivity of our pipeline. The fact that we can only place upper bounds on the false alarm rates (FARs) for these events, in particular, even for the least significant event in GWTC-1, GW170729, is consistent with an improvement in sensitivity over the range of parameters included in our search. Also consistent with this is the fact that we find six new events with (a) false alarm rates below the threshold of 1 in 30 days defined in GWTC-1 (even after applying appropriate trials factors), and (b) probabilities of astrophysical origin $p_{\rm astro}>0.5$, as effectively quantified by their consistency with the rates of occurrence of louder events.

Figure 1 places these new events in context of the population of events in the GWTC-1 catalog by showing their marginalized posterior probability distributions in the plane of the effective spin-parameter, $\chi_{\rm eff}=(m_{1}\chi_{1,z}+m_{2}\chi_{2,z})/(m_{1}+m_{2})$ (where $m_{i}$ and $\chi_{i,z}$, respectively, are the masses and projections of the dimensionless spin onto the orbital angular momentum), and the total source-frame mass. We chose to present the points for the new events colored according to their estimated values of $p_{\rm astro}$ to emphasize that these events come as a population, with these numbers being an integral part of their interpretation. In the rest of the paper, we will present the changes to our analysis methods, the procedure we use to assign values of FAR and $p_{\rm astro}$ to events, and a deeper view into the results.

III Changes to the O1 analysis pipeline

Our analysis pipeline is similar in overall structure to the one we used in the O1 analysis . The overall flow of the search is described in Section II of TV19: for completeness, we briefly repeat the description here. We divide our search space into a set of banks, and within each bank, construct a set of templates that approximate the waveforms generated using the IMRPhenomD approximant Khan et al. (2016) to a desired degree of fidelity in the presence of characteristic detector noise. We then generate a time-series of matched-filtering scores for each template against the data streams from the Hanford and Livingston detectors. When these scores cross a threshold value in a detector, we call this a trigger in that detector: we collect coincident triggers between H1 and L1 (with templates indexed by the same coefficients, and within $10\text{\,}\mathrm{ms}$ of each other). Subsequently we locally refine the triggers in each detector onto a finer template bank, apply vetoes based on signal quality, and pick the best coincident trigger from the subsets of refined H1 and L1 triggers. Finally we assign triggers (and hence the underlying signal/background event) to banks, compute a ranking statistic for the triggers in a bank, and assign a false alarm rate to triggers according to the empirically measured distribution of the ranking statistic.

The interested reader can find a full description of the details of the methods in TV19. In this section, we focus on the differences between the current and previous analysis. In detail, the pipeline for the O2 analysis differs in the following aspects:

Construction of the template bank: Our search covers the parameter space of compact binary mergers with component masses between 3 and 100 $M_{\odot}$, with mass-ratios q $\in[1/18,1]$ and the projections of the dimensionless spins of the constituents on the orbital angular momentum satisfying $|\chi_{i}|<0.99$. As in the O1 analysis, we partition this search space into five banks (BBH 0-4) based on the detector-frame chirp masses, and divide each bank into sub-banks based on the shape of the frequency-domain amplitude profile. Within each sub-bank, we use the procedure described in Ref. Roulet et al. (2019) to construct grids of templates. We now use the $k$-means algorithm to automatically divide each bank into sub-banks, each of which has a frequency-domain amplitude, $\overline{A}(f)$, that is the root-mean-square average of the amplitudes of the constituent astrophysical waveforms. In Fig. 2 we show $\overline{A}(f)$ for all the sub-banks in our banks covering detector chirp masses 20–40 ${\rm M}_{\odot}$ (bank BBH 3) and above $40\,{\rm M}_{\odot}$ (BBH 4). We show the profiles only for these banks to avoid overcrowding the figure, and since all the events we report in this paper fall in these banks. The curves for the other banks (BBH 2, BBH 1, and BBH 0) continue the trend that is visible in Fig. 2, in which the cutoff frequency shifts rightward as we go to lower masses.

The phases of the templates are linear combinations of basis functions, $\psi_{\alpha}(f)$, whose form depends on the noise power spectral density (PSD); we use the PSD estimated from a representative set of files from the run (instead of a model PSD as was done in the O1 analysis). We also restrict the templates to frequencies between $24$ and $512\text{\,}\mathrm{Hz}$. Cumulatively, these changes enable us to cover the same astrophysical parameter space using $\sim 30\%$ fewer templates, and achieve slightly better effectualness Roulet et al. (2019). 2. 2.

Preprocessing and flagging the data: Our analysis pipeline produces a stream of whitened data per $4096\text{\,}\mathrm{s}$ file. As in our O1 analysis, we perform several tests on the data to identify prolonged or transient disturbances in the detector that appear as segments with excess power, and discard these segments to avoid polluting our search. The overall nature and number of these tests is unchanged from the O1 analysis; a full description of these tests can be found in Section III C of TV19. What is different, however, is how we set thresholds to trigger these tests. If the thresholds are too low, they can be triggered by real astrophysical events. Conversely, if they are too high, the tests are not easily triggered by bad data segments. To mitigate the former effect, in the O1 analysis, we set the threshold for each test to the power achieved by signals with a fiducial signal-to-noise ratio $({\rm SNR})_{\max}$ in the absence of noise. This method of setting thresholds is problematic, since in the presence of additive Gaussian random noise, real signals with ${\rm SNR}<({\rm SNR})_{\max}$ can have excess power above these thresholds due to upward fluctuations. Hence, in our O1 analysis we used a relatively high value of $({\rm SNR})_{\max}=30$ to ensure completeness at lower values Venumadhav et al. (2019b).

In the analysis described in this paper, we instead set the thresholds for tests by demanding a given false-positive probability for signals with a desired signal-to-noise ratio. For whitened data, power in a given band, and on a given timescale, is distributed according to a non-central chi-squared distribution in the presence of a signal; hence, we can set thresholds such that signals with a given value of $({\rm SNR})_{\max}$ are flagged with a probability $<10^{-4}$. Given this promise, we are able to set thresholds with a lower target value $({\rm SNR})_{\max}=20$, without losing relevant signals in this search. By nature, the thresholds resulting from this criterion, and hence the amount of data flagged as bad, are dependent on the bank (i.e., at the same value of SNR, the thresholds on excess power are higher for a signal in BBH 4 compared to those in BBH 0). Regardless of bank, the total length of bad data never exceeds 2 % of the overall runtime of O2, and hence, the effective reduction in ‘runtime’ is not an important factor to take into account in spacetime volume estimates for the search. 3. 3.

Refining coincident triggers: The phases of our templates belong to a vector space, $\mathcal{V}$, spanned by the basis functions $\{\psi_{\alpha}(f):\alpha=1,2\dots n\}$; the template bank is a discrete subset of this space with basis coefficients that live on an $n$-dimensional grid. We determine the extent of the grid in all dimensions by projecting a large random sample of astrophysical waveforms into $\mathcal{V}$ and ensuring that every waveform has a nearby grid point. We allow comparatively large mismatches ($\lesssim 10\%$) between astrophysical waveforms and the best template in the bank, which enables us to work with coarser grids on $\mathcal{V}$ when generating triggers. We then reduce the mismatch for significant triggers by refining their coefficients on a finer local grid. In the O1 analysis, we chose this grid to be a uniform regular grid centered on a trigger of interest Venumadhav et al. (2019b).

For heavy BBH waveforms, the set $\mathcal{S}$ of projected astrophysical waveforms is typically thin and mildly curved in higher dimensions ($\alpha\gtrsim 2$) Roulet et al. (2019). The strategy used in the O1 analysis can cause us to step outside $\mathcal{S}$ in these dimensions, and introduce unphysical degrees of freedom that pick up noise but no signal. In this analysis, we change the spacing of the finer grid and excise unnecessary elements to ensure that we enumerate over templates within $\mathcal{S}$. Note that different choices of how this refinement is done can lead to different quoted SNRs for the same astrophysical signal depending on how closely the finer grid approaches it. Hence, we apply the same strategy to the background triggers (found via time slides) to avoid biasing the calculation of false-alarm rates. 4. 4.

Reducing cross-contamination between banks: High ${\rm SNR}$ triggers tend to appear in several of our chirp-mass banks, both in the time slides used to estimate our background, as well as in the set of coincident triggers. In the O1 analysis, we assigned triggers to the best sub-bank in a given chirp-mass bank (as determined by the incoherent network ${\rm SNR}^{2}=\rho_{\rm H}^{2}+\rho_{\rm L}^{2}$, where $\rho_{\rm H}^{2}$ and $\rho_{\rm L}^{2}$ are the incoherent squared ${\rm SNR}$s in Hanford and Livingston, respectively), but allowed them to appear in multiple banks Venumadhav et al. (2019a). This choice was conservative, in that it caused us to overestimate the FAR for real events. In this analysis, we improve upon this in two ways: we assign both background and coincident triggers to a unique bank (and sub-bank within), and instead of the incoherent network ${\rm SNR}^{2}$, we use a discriminator that better accounts for the different structures of the sub-banks. Appendix B contains the expression (see Eq. (13)), and outlines a derivation. 5. 5.

Computing the false-alarm rate:

After collecting triggers and assigning them to banks, we estimate the FAR for a given coincident trigger by comparing it to the background triggers (collected using time slides) within the same chirp-mass bank. Our statistic for comparing triggers is the coherent score Nitz et al. (2017), which for a given trigger $t$ is an approximation of the likelihood ratio under the signal ($\mathcal{L}(t\mid\mathcal{S})$) and noise ($\mathcal{L}(t\mid\mathcal{N})$) hypotheses. In order to estimate $\mathcal{L}(t\mid\mathcal{N})$, we approximate the distribution of background triggers in each detector (as a function of ${\rm SNR}^{2}$) by its survival function. The survival function (or more correctly, the complementary cumulative distribution function) has the virtue that we can estimate it empirically by just ranking triggers as a function of the value of ${\rm SNR}^{2}$ in that detector (see Section III J of TV19).

In the O1 analysis, when estimating $\mathcal{L}(t\mid\mathcal{N})$, we assumed that the background rate per template is flat over all sub-banks in a given bank Venumadhav et al. (2019a). We revisited the validity of this assumption while developing the analysis in this paper. Figure 3 shows the rank functions for the five sub-banks of the bank BBH 4 as a function of the ${\rm SNR}^{2}$ in L1, $\rho_{\rm L}^{2}$ (these curves are completely empirically determined and not fits; the quantity on the y-axis is twice the logarithm of the number of triggers above a particular value of $\rho_{\rm L}^{2}$, with an offset. See discussion in Section III J of TV19, and below.). The figure clearly shows that the sub-banks have substantially different background-rates per template. The distributions of $\rho_{\rm L}^{2}$ transition from exponential (chi-squared-like) to power-laws when glitches become important, and hence the rank functions flatten: this happens at progressively lower values of $\rho_{\rm L}^{2}$ for higher-mass sub-banks. If we do not account for this, a sub-bank with more glitches can disproportionately influence the background in the search.

In this search, when computing the false alarm rates of a trigger, we use the rank function of the sub-bank it is drawn from. This enables us to give a more local estimate of the false-alarm rate, which better accounts for the fact that templates in some regions of the search space are more prone to triggering on glitches (this effect has been noted before Dal Canton et al. (2014); Bose et al. (2016)). A technical point is that we compute FARs over all templates in a bank (which is the union of several sub-banks), and hence we need to properly fix the normalization of the rank functions when comparing triggers from different sub-banks. Appendix C describes the procedure; the rank functions shown in Fig. 3 were normalized in this way. Curiously, we see that the rank functions have different slopes in different sub-banks even at low values of $\rho_{\rm L}^{2}$ (in the Gaussian-noise-dominated regime); this is due to the maximization over templates when we collect triggers (see Appendix B for a derivation in a different context). 6. 6.

Computing $p_{\rm astro}$: We define $\mathcal{R}{({\rm event}\mid\mathcal{N})}$ and $\mathcal{R}{({\rm event}\mid\mathcal{S})}$ to be the rates of a given event under the noise ($\mathcal{N}$) and signal ($\mathcal{S}$) hypotheses. The probability that an event is astrophysical is

[TABLE]

We define a rate $\mathcal{R}$ for each bank to be the overall number of astrophysical events satisfying $\rho_{\rm H}^{2},\rho_{\rm L}^{2}>16$ and $\rho_{\rm H}^{2}+\rho_{\rm L}^{2}>\rho^{2}_{\rm th}$, where

[TABLE]

The rate $\mathcal{R}$ is for a hypothetical network consisting of two identical detectors, each having a sensitivity equal to the median Hanford sensitivity in the O2 run, which observe in coincidence for $118$ days. $\mathcal{R}$ is assumed to be uniform across templates within the bank, regardless of which sub-bank they might fall in. For a given event, we have

[TABLE]

where the factor $W$ depends on the instantaneous sensitivities of the detectors, as well as the extrinsic parameters of the event. Note that $W$ does not depend on the (unknown) astrophysical rate. We determine the terms in Eqs. (1) and (2) in a similar manner to that of our O1 analysis Venumadhav et al. (2019b), but with two changes. Firstly, we estimate the rate of producing triggers under the noise hypothesis, $\mathcal{R}{({\rm event}\mid\mathcal{N})}$, using only the background triggers in the respective sub-banks that the candidates belong to. Secondly, we determine the rate of astrophysical events, $\mathcal{R}$, from the data itself rather than assuming it from the loudest events detected in the search.

Given a particular value of the rate $\mathcal{R}$, the likelihood of the data is

[TABLE]

where the product is over all the triggers in the bank, including those detected originally by the LVC. Using this likelihood, we compute a posterior on $\mathcal{R}$, assuming a uniform prior $P(\mathcal{R})$ between 0 and 50. As earlier, this number is the overall number of triggers of astrophysical origin in a region in the $\rho_{H}^{2}$–$\rho_{L}^{2}$ plane for a hypothetical two–detector network with each detector having the median H1 sensitivity over the O2 run, i.e., a period of 118 days. We expect this number to be of the order of the number of detections, and determine it from the data itself; hence the most important requirement is that the prior be non-informative over the relevant range. For any value of the overall rate $\mathcal{R}$, we can calculate the probability that an event is of astrophysical origin; our final quoted values were obtained by marginalizing over $\mathcal{R}$:

[TABLE]

IV Results on the previously reported events

Table 1 summarizes our pipeline’s results for the O2 events published by the LVC Abbott et al. (2018). We detect all previously reported BBH events except for GW170608, for which the LVC did not provide the Hanford data in their bulk data release, and thus that time was not part of our coincidence search. We only report results from our BBH search in this paper, so we exclude the binary neutron star GW170817 from our results. Nearly all of the LVC events have only an upper limit for the FAR of $1/($20,000$\,{\rm O2})$. All of these events are certainly astrophysical sources with lower limit on $p_{\rm astro}\geqslant 0.99$.

An interesting special case is GW170818, which was not found by the PyCBC pipeline, and deemed potentially interesting but not confirmed by the GstLAL pipeline using Hanford and Livingston alone Abbott et al. (2018). It was subsequently detected with high confidence by GstLAL when Virgo data was included; we did not analyze Virgo for the search reported in this paper.

Due to the low score in H1, this is close to a single-detector trigger, and hence the detection heavily relies on our understanding of the background distribution in the vicinity of the L1 ${\rm SNR}^{2}$. Empirically, our background distributions do not have any L1 triggers that are louder than this event. In the regime relevant to triggers like GW170818, the background is dominated by only a few loud glitches, and hence the L1 ranking score, $\tilde{\rho}_{\rm L}^{2}$, that we use to define our test statistic saturates (see Fig. 3) and comes with significant error bars. Moreover, we cannot compute the local probability density of the background in the $\rho_{\rm H}^{2}$–$\rho_{\rm L}^{2}$ plane without extrapolating from lower values, and this is a vital ingredient in the calculation of $p_{\rm astro}$ (see Eq. (1)). Hence we need a more careful analysis to reliably assess the FAR and probabilities of astrophysical origin of such events. We defer this analysis to a follow-up paper, in which we will show that such events can be detected using the Livingston and Hanford detectors alone. In this subsequent work, we will also report the results of a search for similar events, i.e., events which are loud in one detector (and saturate the rank score), and which either have low scores in the other detectors, or have no coincident data.

The following technical details about the quantities we report in Table 1 (and subsequent tables) are worth noting for clarity:

• •

Our false alarm rates (FARs) are calculated per bank, and do not include any additional trials factors. Hence, the reader should interpret these numbers in the context of the set of banks used in the search. We adopt this convention since (a) FARs reported this way are invariant even if we subsequently expand the search space (i.e., add extra banks), and (b) it makes it convenient for readers to apply their own priors when computing the look-elsewhere correction, or the trials factor (i.e., the interpretation of a particular FAR is different in regions of parameter space known to have larger astrophysical rates).

• •

We report the ‘linear-free’ times of the best fit templates in the template bank in the Hanford detector. The ‘linear-free’ time is defined as follows: suppose that for a particular set of intrinsic parameters, we have a frequency domain template $h(f)=A(f)\exp{\left[i\psi(f)\right]}$, where the phase $\psi(f)$ is unwrapped as a function of frequency. Suppose through the inverse Fourier transform, we have a time-domain template $h(t)=\int{\rm d}f\,h(f)\exp(i2\pi ft)$, defined on the domain $[0,T]$. The ‘linear-free’ time of the template $h(t)$ occurs at $t=T-t_{0}$, where $t_{0}$ satisfies the equations

[TABLE]

In the above equation, $S(f)$ is the PSD of the detector noise. Intuitively, after applying the parameter-dependent shift $t_{0}$ to the template, the measurements of arrival time and intrinsic parameters are uncorrelated in a single detector (see related discussion in Ref. Roulet et al. (2019)). When using Eqs. (5)-(6), i.e., only to define $t_{0}$, we used the same fiducial PSD that we used to construct the template bank as described in Sec. III. For the short waveforms in the heavier binary black hole banks, the linear-free time is practically close to the merger time.

V Overview of the newly discovered BBH mergers

Table 2 summarizes the basic properties of the newly discovered events: their parameters, inverse FAR, or IFAR, and estimated probabilities of being of astrophysical origin, $p_{\rm astro}$. We need the rate of astrophysical events in the detector, $\mathcal{R}$, to calculate these values of $p_{\rm astro}$, and hence, we also report the values of the factor $W$ (from Eq. (2)), so that external sources of information about the rate (e.g., future runs) can be easily incorporated if needed.

We derive a distribution for $\mathcal{R}$ from all the events we detected (including those already detected by the LVC analyses), and at the same time estimate the $p_{\rm astro}$ of each event, using the procedure described in Section III. Figure 4 shows the posteriors on the rates for chirp-mass banks BBH 3 and BBH 4; the values of $p_{\rm astro}({\rm event})$ quoted in Tables 1 and 2 were marginalized over these posteriors. As is apparent from Table 2, all of our detections are in the high-mass region of parameter space, which is covered by the banks BBH 3 and BBH 4. Hence, curves analogous to those in Fig. 4 for the other banks only yield upper bounds on the respective rates (of the mergers detected in previous runs, GW151012 and GW151226 lie in banks BBH 2 and BBH 1, respectively).

As in Table 1, the IFARs in Table 2 are calculated per bank, and do not include any additional trials factors. The criterion in Table 2 is based on the value of $p_{\rm astro}$, but all the events have significantly higher values of the ‘per-bank’ IFAR than the threshold of 1 in 30 days, or $\sim 0.25$ O2, adopted in GWTC-1 Abbott et al. (2018) (i.e., their IFARs exceed 0.25 O2 even if divided by a trials factor of five).

Figure 5 shows the spectrograms for segments of data of length 1 s around the GPS times of the events in the Hanford and Livingston detectors. These spectrograms show no obvious evidence of glitches in the immediate vicinity of the events, and are visually consistent with how binary black hole merger signals with the reported values of SNR would look like in the presence of additive Gaussian noise. Note that visual inspection of spectrograms is not a particularly good, or even well-defined way to identify either glitches or signals in the data: our pipeline automatically flags bad segments of data using mathematically well-defined test statistics (see TV19 for details), and looks for signals using matched filtering.

Appendix A includes the posteriors for the parameters of all the new events. Figure 1 places these events in the context of the previous LVC events, as well as the one we reported in Ref. Zackay et al. (2019), by showing their distribution in the plane of source-frame total mass and the effective spin parameter $\chi_{\rm eff}$. In the remainder of this section, we briefly comment on the properties of each of the newly found events.

GW170121:

This event has the lowest FAR among those not reported by the LVC (${\rm IFAR}\approx 2.8\times 10^{3}\,\rm O2$), and $p_{\rm astro}>0.99$. The posterior distribution of $\chi_{\rm eff}$ has most of its support at negative values, $\chi_{\rm eff}>0$ being ruled out at the $95\%$ confidence level. The chirp mass, mass-ratio and redshift of this event are similar to those of the heavy BBHs reported by the LVC. Its redshift is $z\sim 0.24$.

GW170304 & GW170727:

These two events have ${\rm IFAR}\sim 370\,\rm O2$ and $p_{\rm astro}\approx 0.98$. Their masses and spins are similar to those of the heavy BBHs detected by the LVC. Both events are consistent with zero $\chi_{\rm eff}$, and are on the massive end of the population. They have relatively high redshifts $z\sim 0.5$ and $0.43$, respectively.

GW170425:

This candidate has $p_{\rm astro}\approx 0.77$ and ${\rm IFAR}\approx 29\,\rm O2$. Its inferred parameters are similar to those of the heavy BBHs reported by the LVC; the effective spin $\chi_{\rm eff}$ is consistent with zero. The posterior distribution has a tail extending to large values for the masses. Its inferred redshift is large, $z\sim 0.5$.

GW170202:

This candidate has $p_{\rm astro}\approx 0.7$ and ${\rm IFAR}\approx 6\,\rm O2$. The masses and the spins are similar to those of the heavy LVC BBHs. It is found in the bank with the largest number of secure detections (BBH 3). It has a bimodal posterior, in which the solution with lower masses has a more negative spin, and is located closer. The inferred redshift is $z\sim 0.27$.

GW170403:

This candidate has ${\rm IFAR}\approx 5\,\rm O2$ and $p_{\rm astro}\approx 0.55$; this is close to the threshold $p_{\rm astro}=0.5$ to make it into a list of detections (as defined in Ref. Abbott et al. (2018)). The inferred redshift is $z\sim 0.45$. Interestingly, the posterior for $\chi_{\rm eff}$ is inconsistent with positive values.

In addition to these events, we list in Table 3 the sub-threshold triggers of our search, defined as those with $0.1<p_{\rm astro}<0.5$. The sum of the $p_{\rm astro}$ of the events in this list exceeds unity; in fact, a candidate in bank BBH (4,1) has $p_{\rm astro}\approx 0.45$, which is close to the detection threshold (though it has a relatively high ${\rm IFAR}\approx 0.8\,$O2). It is possible that an improved analysis, or rate-estimate, can push some of these candidates above the detection threshold. Note that we have a threshold of $p_{\rm astro}>0.1$ for candidates to appear in Table 3. The set of loudest triggers in banks BBH 2, BBH 1, and BBH 0 did not pass this cut.

VI Sensitivity of our pipeline

In the previous section, we described several additional events we detected that are not in the catalog of events published by the LVC. All of these events pass the thresholds for detection in GWTC-1 Abbott et al. (2018) (their FARs are above the threshold of 1 in 30 days, even accounting for the five banks in our BBH search, or eleven banks in a hypothetical binary neutron star and neutron-star–black-hole search Roulet et al. (2019), in the extreme scenario in which searches for different kinds of astrophysical systems are considered together).

The values of the FAR and $p_{\rm astro}$ of the events come with error bars, and for the same event, their values can fluctuate even between two searches that use the same algorithms due to choices within. Given that searches operate with thresholds, the presence or absence of a near-threshold event in one search, by itself, is not proof of a difference in sensitivity. However, equally sensitive searches should agree on the overall population of events, and, on average, assign comparable values of FAR and $p_{\rm astro}$ to comparable events. Our results in Sec. IV show that we assign systematically lower values of FAR to the events discovered by the LVC pipelines in general, and in particular, for GW170729, the least secure event in GWTC-1. Moreoever, we have a substantial number of events that individually clear the thresholds for detection; while it is certainly possible that some of the events with lower values of $p_{\rm astro}$ could have benefited from upward fluctuations in their test statistics, it is hard to explain the entire population away this way. These two effects suggest that our analysis has a substantially larger sensitive volume within the search space that we defined in Sec. III.

In an idealized case, a simple way to compare sensitive volumes would be to inject a large number of signals in the data and measure the fraction recovered above the appropriate thresholds in the test statistics adopted. We will present the results of such a systematic injection campaign in a future paper; pending this, we would like to obtain a simple estimate of the additional sensitivity. In the rest of this section, we will consider the population of loud astrophysical events in this data akin to a common set of injected signals, and get a simple ‘back-of-the-envelope’ level estimate of the order of magnitude of the sensitivity change. Since there are relatively few events, this is necessarily highly uncertain, and hence the results should not be over-interpreted.

Figure 6(a) shows the background triggers we collected using $20\,000$ time slides in those BBH sub-banks in which all the events considered in this work, both from the LVC and our analysis, reside. This figure does not include the BBHs from the O1 run (GW150914, GW151012, GW151216, GW151226), nor GW170608, which was not included in the bulk data release we analyzed. This figure is not intended as a demonstration of how we compute the FAR or $p_{\rm astro}$ for particular events: firstly, it shows $\rho_{\rm H}^{2}$ and $\rho_{\rm L}^{2}$, i.e., the incoherent H1 and L1 ${\rm SNR}^{2}$, while we compute the FAR using a coherent score that takes into account the time-delays and the relative phases of the triggers, and the differing detector sensitivities; secondly, we estimate the FAR and the $p_{\rm astro}$ for a particular event using the background in its chirp-mass bank, and sub-bank, respectively. We include this figure only to easily visualize the sensitive volume.

It is clear from Fig. 6 that our pipeline has substantially lower background in the relevant region: for example, we see no background triggers within the sensitive region of the LVC pipelines. All but one of the LVC reported events have values of $\rho_{\rm L}^{2}$ that are so large that we do not have even single-detector background triggers at their level (this is a consequence of our data-cleaning procedure, as well as our signal-quality vetoes). The only exception to this is GW170729, with $(\rho_{\rm H}^{2},\rho_{\rm L}^{2})=(62,53)$. This event had FARs of $0.2\text{\,}\mathrm{y}\mathrm{r}^{-1}$ and $1.36\text{\,}\mathrm{y}\mathrm{r}^{-1}$ in the GstLAL and PyCBC pipelines, respectively, but we have no background in its vicinity even incoherently (i.e., allowing for arbitrary phases, time-delays, and sensitivity ratios of the two detectors).

In order to calculate sensitive volumes, we need to estimate the detection limits of the analysis. With the above caveat on the validity of thresholds based on incoherent ${\rm SNR}^{2}$, the solid and dashed lines show the approximate detection thresholds that we judged as appropriate for different analyses. The detection thresholds shown for the LVC catalog are approximate and conservative, they err on the side of reporting a better sensitivity. At the single-detector level, we set the threshold to $\rho_{\rm H}^{2},\rho_{\rm L}^{2}>30$ by the non-detection of GW170121. We set the minimum network ${\rm SNR}^{2}=\rho_{\rm H}^{2}+\rho_{\rm L}^{2}>90$ by scaling the reported FAR of GW170729 to 1/O2, and rounding down. We approximated our incoherent limit as $\rho_{\rm H}^{2}+\rho_{\rm L}^{2}>68$, and $\rho_{\rm H}^{2},\rho_{\rm L}^{2}>16$, based on our cuts, and the FAR we would assign to events at this level given our search background.

The difference in the detection limits, at the same threshold on FAR, maps to a difference in the sensitive volume between the searches. The size of this difference depends on the ratio of the sensitivities of the two detectors we analyze, i.e., H1 and L1. The sensitivities $n_{\rm H}$ and $n_{\rm L}$ are proportional to the SNR with which a gravitational wave signal with a given strain amplitude is measured by H1 and L1, respectively. Figure 6(b) shows the volume as a function of the sensitivity ratio (scaled such that a volume of $V_{\rm max}$ corresponds to detecting all events with $\rho_{\rm H}^{2}+\rho_{\rm L}^{2}>68$); the solid orange and the dashed blue lines are for the cuts shown in Fig. 6(a). The curves were obtained by randomizing the angular locations and inclinations of a large number of mergers on the sky and recording the incoherent scores (including the stochastic noise contribution) at the detectors with a given sensitivity-ratio. We see that the gain between the blue and range lines depends significantly on the relative sensitivities of the H1 and L1 detectors. The filled histogram in Figure 6(b) shows the distribution of the ratio of the Hanford and Livingston sensitivities in our analysis of the O2 run, weighted by the momentary space-time volume. We see that the orange curve gains approximately a factor of two in sensitive volume over the relevant range of sensitivities. As emphasized in the beginning of this section, this is only a rough, order-of-magnitude estimate of the relative sensitivity using a small number of events. This number comes with significant uncertainty, and as such, should not be over-interpreted.

Also of particular interest is the limit in which the ${\rm SNR}$ is much larger in one detector than in the other. For part of the O2 run, the Livingston detector was substantially more sensitive than the Hanford one, and hence there is a substantial phase-space volume for astrophysical signals to have disparate ${\rm SNR}$ in the detectors. Figure 6(a) shows that there is substantially less background in L1 in the high ${\rm SNR}$ regime (see the teal shaded region). The few background events in this region come from the same small number of loud events in L1 matching with Gaussian fluctuations in H1 at different time shifts, and hence the ranking function $\tilde{\rho}_{\rm L}^{2}$ saturates and is severely affected by Poisson noise (see Fig. 3). We need a different analysis to estimate a meaningful FAR in this regime; the green, dashed-dotted curve in Figure 6(b) shows that we can gain a non-trivial amount of sensitive volume if we open up this region. We will study this regime in more detail in a subsequent paper, in which we will introduce a formalism for searching for events and assessing the FAR in this region. Notably, the LVC event GW170818, which was detected using Virgo data, belongs to this category.

VII Conclusions

In this paper we presented the results of our search for BBHs in the data from the O2 observing run of advanced LIGO, using the methods introduced in our work in Ref. Venumadhav et al. (2019b). We report six new events above the detection thresholds defined by the LVC (in terms of FAR and $p_{\rm astro}$), three of whom have probability $p_{\rm astro}>0.98$ of being of astrophysical origin. Interestingly, all the new events are in banks BBH 3 and BBH 4 (our heavy chirp-mass banks), as are most of the ones reported by the LVC.

The most significant new event (GW170121) prefers negative $\chi_{\rm eff}$ and is inconsistent with positive values at the $95\%$ level. The most marginal candidate event (GW170403), with ${p_{\rm astro}}\sim 0.5$, is inconsistent with zero or positive $\chi_{\rm eff}$. The spin of the merging BBHs is an important discriminator between formation channels Farr et al. (2018). Hence, the new events presented in this work can throw light on the mechanisms by which BBHs are assembled.

More generally, with the increased number of events, the clear next step is to perform a population analysis that accounts for selection biases, which will map out the distribution of the intrinsic parameters of the mergers. In particular, including new events in population analyses can significantly inform us about the dependence of the merger-rate on mass and redshift.

The LVC recently started their third observing run (O3), and several new detections are expected in the near future. The new events we report in this paper show that there will be additional information in the LIGO and Virgo data in addition to what the pipelines used by the LVC currently extract. The development of our pipeline has been facilitated by access to the public O1 and O2 data, as well as the LIGO Algorithm Library LIGO Scientific Collaboration (2018). We thank the LVC for releasing the data and tools to the community. We hope that data from current and future runs can be made available quickly to incentivize external groups to develop new analysis methods and maximize the scientific yield of the LVC data.

Acknowledgements

This research made use of data, software and/or web tools obtained from the Gravitational Wave Open Science Center (https://www.gw-openscience.org), a service of LIGO Laboratory, the LIGO Scientific Collaboration and the Virgo Collaboration. LIGO is funded by the U.S. National Science Foundation. Virgo is funded by the French Centre National de Recherche Scientifique (CNRS), the Italian Istituto Nazionale della Fisica Nucleare (INFN) and the Dutch Nikhef, with contributions by Polish and Hungarian institutes.

TV acknowledges support by the Friends of the Institute for Advanced Study. BZ acknowledges the support of The Peter Svennilson Membership fund. LD acknowledges the support from the Raymond and Beverly Sackler Foundation Fund. MZ is supported by NSF grants AST-1409709, PHY-1521097 and PHY-1820775 the Canadian Institute for Advanced Research (CIFAR) program on Gravity and the Extreme Universe and the Simons Foundation Modern Inflationary Cosmology initiative.

Appendix A Properties of the new events

We performed full parameter estimation for all the new events presented in this paper. We coherently analyzed the data from all available detectors (including Virgo), in contrast to the search, which only used Hanford and Livingston data. We use a likelihood model that assumes the data is the sum of a gravitational wave signal and additive Gaussian random noise, and a prior that is uniform in the intrinsic source-frame parameters $m_{1}$ and $m_{2}$, the effective spin parameter $\chi_{\rm eff}$, and luminosity volume (as detailed in Ref. Zackay et al. (2019)). We use the IMRPhenomD approximant Khan et al. (2016) to generate waveforms, the relative binning method Zackay et al. (2018) to evaluate the likelihood, and PyMultiNest Buchner et al. (2014) to generate samples from the posteriors.

Figures 7 and 8 show the posterior distributions of detector-frame chirp mass, mass-ratio, effective spin, and redshift for the new events reported with $p_{\rm astro}>0.5$, marginalized over the extrinsic parameters of the mergers such as their sky location and inclination. Figure 9 shows the posteriors for the location of all the events on the sky, marginalized in turn over intrinsic parameters and inclination: we only show the 50% and 90% contours in order to present all the events on the same plot. Samples from the full posterior distribution (i.e., the joint distribution of the intrinsic and extrinsic parameters) are available at https://github.com/jroulet/O2_samples.

Appendix B Populating events in sub-banks and banks

In our analysis, we ultimately assign events to a single chirp-mass bank. We compute false alarm rates for an event within a particular bank by comparing its corresponding trigger to the triggers from background events (i.e., those generated by noise) in that bank. However, a single signal or noise transient produces triggers in several sub-banks, and even across chirp-mass banks. If we do not properly assign events (both real and background) to banks, we run the risk of losing sensitivity in our search, and hence we should carefully choose our criteria for bank assignment888This is not unique to our search. For an extreme example, GW150914 produced triggers across the entire binary black hole parameter space in the LVC search Abbott et al. (2016). If one of the fainter triggers had been picked as representing the event, it would have been assigned a different false alarm rate..

We would like to perform the analysis as blindly as possible, and hence would like to avoid tweaking criteria after touching the data. Since we are analyzing a new dataset in this work, this gives us an opportunity to update the bank assignment procedure from the one we used for the O1 data (step 6 of Section II of TV19). We begin by considering the distribution of the triggers that a particular signal in the data produces across the banks. This distribution depends on the number of templates placed in each bank (effectively their density in parameter space), as well as the effective dimensionality of the parameter space in each bank/sub-bank (the dimensionality has a concrete meaning in the $c_{\alpha}$ parameter space, see Ref. Roulet et al. (2019)). Given this understanding, we use the set of triggers that an event throws up in a given bank to compute a net likelihood, and use this as a discriminator for assigning events to banks. In this rest of this appendix, we describe the details of this procedure.

For simplicity, let us start with triggers in a single detector. Let $d$ denote the strain data, $A$ be the signal amplitude, and $\mathbf{\Theta}$ be the other parameters of the signal: these include the coefficients $c_{\alpha,i}$ in sub-bank $i$ of chirp-mass bank $B$, and the phase and the merger-time. The likelihood of the sub-bank $i$ under the signal hypothesis is

[TABLE]

where the sum runs over the grid points in $c_{\alpha,i}$, and the integral is over continuous parameters (time and phase). We can view the sum over the coefficients $c_{\alpha,i}$ as a Riemann sum for an integral over the underlying continuous space, $\mathcal{V}$, and thus approximate Eq. (7) as

[TABLE]

where $\Delta c_{\alpha,i}$ and $n_{c_{\alpha,i}}$ are the spacing and dimensionality of the template grid in sub-bank $i$. The factor in front of the integral is the volume per template in the discrete grid.

We adopt the following assumptions:

In the prior $p(\mathbf{\Theta})$, astrophysical signals are equally likely to occur in two different chirp-mass banks. 2. 2.

The astrophysical rate per template is uniform within a chirp-mass bank, i.e., given that a signal occurs within bank $B$, the probability that it occurs with a given template equals $1/N_{{\rm temp},B}=1/\sum_{i}N_{{\rm temp},i}$, where $N_{{\rm temp},B}$ and $N_{{\rm temp},i}$ are the number of templates in the bank $B$ and sub-bank $i$, respectively (we make the same assumption when estimating the detector rate of events in Section III). 3. 3.

The integral in Eq. (8) receives most its contribution from around the best fit parameters, $\mathbf{\Theta}_{*}$, and amplitude $A_{*}$, which is valid in the limit of high ${\rm SNR}$.

The merger-time and phase are uniformly distributed within their ranges, while the amplitude has a prior distribution $p(A)\sim 1/A^{4}$ in a Euclidean universe. The integrand $\mathcal{L}(d\mid A,\mathbf{\Theta})$ has the form

[TABLE]

where $t(\mathbf{\Theta})$ is the template, and the inner product is weighted by the inverse PSD. Next, we integrate over the amplitude in Eq. (8). Under assumption 3 above, this gives us a prefactor $A_{*}^{-4}$, and simplifies the integrand to

[TABLE]

where we have removed a term that does not depend on the templates (this step is identical to the standard derivation of the $\mathcal{F}$-statistic Jaranowski and Królak (2012)). We express the best fit amplitude as $A_{*}=(\rho/\rho_{0})A_{0}$, where $\rho_{0}$ is the ${\rm SNR}$ of a merger at a fiducial distance and orientation, which captures the instantaneous sensitivity of the detector.

By assumption 3 above, the integrand in Eq. (10) is sharply peaked around the bestfit parameters $\mathbf{\Theta}_{0}$, there it equals $\exp(\rho^{2}/2)$, where $\rho^{2}=\langle A_{*}t(\mathbf{\Theta})\mid A_{*}t(\mathbf{\Theta})\rangle$ is the ${\rm SNR}^{2}$. A nice feature of our template banks is that the deviation in coefficients, $\delta c_{\alpha}$, directly measures the degradation in the overlap between templates Roulet et al. (2019):

[TABLE]

A similar relation holds for the other continuous parameters (time and phase), i.e the degradation in overlaps is quadratic with displacement (the principal directions are some linear combinations of time and phase), and the width is inversely related to $\rho$. We simplify the integrand of Eq. (10), substitute it into Eq. (8), and use the prior on the templates from assumptions 1 and 2 to obtain

[TABLE]

When we refine coincident triggers between two detectors, the $\mathbf{\Theta}$ contains four extra parameters apart from the template coefficients: two times, and two phases. A complication is that astrophysical signals are not uniformly distributed in the space of the time delay and relative phase, and the distribution depends on the relative sensitivities of the two detectors, $\rho_{0,{\rm H}}/\rho_{0,{\rm L}}$. Thus there is a nontrivial prior in the space of parameters $\mathbf{\Theta}$.

We make progress by noting that the likelihood in the integrand in Eq. (8) is independent of the relative times and phases. In this case, we are greatly helped by our assumption 3 above, which tells us that all we need is to evaluate the prior at the bestfit parameters $\mathbf{\Theta_{*}}$. We combine the prior and amplitude prefactor together into a function $p(\Delta t,\Delta\phi,\rho_{\rm H}^{2},\rho_{\rm L}^{2}\mid\rho_{0,{\rm H}},\rho_{0,{\rm L}})$, which we evaluate by Montecarlo sampling methods.

The integral over the likelihood in Eq. (8) can be evaluated in a similar manner as above, with the difference that now $\rho^{2}=\rho_{\rm H}^{2}+\rho_{\rm L}^{2}$. Thus we finally obtain

[TABLE]

Note that the exponent of the last term is different due to the extra degrees of freedom in the two-detector case.

We use the ratio of the likelihoods given by Eq. (13) as a discriminator to compare triggers in the same location in different sub-banks (and possibly chirp-mass-banks).

Appendix C Definition of rank functions in different sub-banks

We compute FARs of events by comparing them to background triggers in their chirp-mass bank. An essential ingredient in this computation is the likelihood $\mathcal{L}(t\mid\mathcal{N})$ for a trigger $t$ under the noise hypothesis. Coincident triggers are produced with random time-delays and relative phases between the two detectors, and the likelihood $\mathcal{L}(t\mid\mathcal{N})$ depends only on the incoherent ${\rm SNR}^{2}$ in the two detectors, and the template in sub-bank $i$ (which we denote by the set of coefficients $c_{\alpha}$ of the basis phase functions), i.e., $\mathcal{L}(t\mid\mathcal{N})=P(\rho_{\rm H}^{2},\rho_{\rm L}^{2},c_{\alpha},i)$. We can write:

[TABLE]

where $P(i)$ is the probability that noise produces a coincident trigger in sub-bank $i$, regardless of the ${\rm SNR}^{2}$, and $N_{{\rm temp},i}$ is the number of templates in sub-bank $i$.

Equation (15) assumes that a) the background is flat over templates within each sub-bank, and b) the triggers in different detectors are independent of each other. Both these assumptions fail to some degree (the latter happens because we refine coincident triggers, and pick the best common template). One solution would be to directly estimate the probability in Eq. (14) from the background, but in practice, the many-dimensional distribution is hard to sample finely enough, and thus real coincident triggers can receive spurious penalties to their scores. We use the same assumptions to rank the time slides as well as the coincident triggers, and hence our FARs are not biased for the strategy we adopt; the price of the above assumptions is that our ranking is no longer strictly optimal.

The probability $P(i)$ is

[TABLE]

where $N_{{\rm trig},i}$ is the number of triggers in sub-bank $i$, and $N_{{\rm trig}}$ is the number of triggers summed over all sub-banks.

We approximate the probability $P(\rho^{2}\mid i)$ for each detector using the same ranking function that we adopted in our previous work Venumadhav et al. (2019b), but compute it separately for each sub-bank:

[TABLE]

where ${\rm Rank}(\rho^{2}\mid i)$ is the ranking of a given trigger in its sub-bank (with the lowest rank given to the loudest event) and $\rho^{2}_{0}$ is a normalization point that we set to $\rho_{0}^{2}\approx 30$. We estimate $P(\rho^{2}_{0}\mid i)$ by taking the ratio between the number of triggers in a bin around $\rho^{2}_{0}$ in sub-bank $i$ over the total number of triggers in the sub-bank, i.e.,

[TABLE]