Mapping the Universe Expansion: Enabling percent-level measurements of the Hubble Constant with a single binary neutron-star merger detection.

Current parameter estimation of BNS signals makes use of waveform templates including only the quadrupole mode. This causes a degeneracy between the inclination and distance parameters that fundamentally limits our ability to measure each. Several works have shown that inclusion of HMs in templates can break this degeneracy for sources with sufficiently-loud HMs in the detector sensitive band, as the observed combination of modes will depend on the orientation of the binary via the $Y_{\ell,m}$ factors (Graff et al., 2015; London et al., 2018; Calderón Bustillo et al., 2018; Pang et al., 2018; Calderón Bustillo et al., 2019). For Advanced LIGO and Virgo observations, unfortunately, this is only possible for large mass and asymmetric BBHs (Graff et al., 2015; London et al., 2018; Chatziioannou et al., 2019; LIG, 2020) for which the merger and ringdown emission, rich in HMs, is strong in the detector sensitive band. In contrast, the merger and post-merger of BNSs is un-observable due to its large frequency. This emission will, however, be observable with future high-frequency (Martynov et al., 2019; Ackley et al., 2020) and third-generation (Punturo et al., 2010a; Reitze et al., 2019a) detectors and, in addition, the pre-merger GW emission will last for several minutes in the sensitive detector band, allowing to accumulate the effect of weak HMs.

We test our ability to measure the source distance and $H_{0}$ using the inspiral and post-merger emission of BNSs. To this, we perform parameter inference on two kinds of simulated GW signals. First, we use 80 ms-long numerical-relativity simulations for the post-merger emission of BNS (Dietrich & Hinderer, 2017). These have mass ratios $q=1$ and $q=1.5$ and implement two different equations of state (EOSs): a soft one (SLy) and a stiff one (MS1b ¹¹1While a stiff EOS like MS1b is disfavoured by the observation of GW170817 and its EM counterparts, e.g., (Abbott et al., 2019, 2018; Dietrich et al., 2020), it provides a good test case for our study to show the effect of two different EOSs.). The left panel of Fig. 1 shows the time domain modes for the equal-mass SLy case located at a distance of $40$Mpc. The right panel shows in blue and green the spectra of full waveforms observed face-on (blue) and edge-on (green), together with the contribution of the HMs and the $(\ell,m)=(2,1)$ mode alone for the edge-on case. It can be noted not only how the face-on signal is stronger, but how the presence of HMs in the edge-on signal leads to noticeable morphological differences. Parameter inference on these short waveforms provides an idea of how the post-merger emission breaks the distance-inclination degeneracy. Restricting to this and ignoring the long inspiral signal, however, would greatly underestimate the signal-to-noise ratio (SNR) accumulated throughout the full minutes-long signal observable by the detector, reducing the accuracy of our measurements.

To obtain signals that can cover the full inspiral and post-merger emission from a BNS, we combine (Bustillo et al., 2015) our short numerical-relativity waveforms (Dietrich & Hinderer, 2017) with long analytical waveforms covering the early-inspiral minutes before the merger (e.g., computed using tidal effective-one-body model of Ref. (Nagar et al., 2018)). Unfortunately, parameter inference using these waveforms is computationally prohibitive. As a solution, we implement a two-step approach. First, we consider 128s-long phenomenological (Phenom) waveforms (Khan et al., 2016; Santamaria et al., 2010), constructed with the IMRPhenomHM model (London et al., 2018), covering the inspiral-merger and ringdown stages of non-precessing binary black hole (BBH) mergers. While computationally inexpensive, these waveforms omit the two main characteristic aspects of BNSs: tidal-deformability effects and the post-merger emission shown in Fig. 1. Since, as argued, the latter can improve our distance measurements, results obtained using this BBH signal model are rather conservative. Finally, in order to obtain improved and more realistic results for BNSs, we combine these BBH parameter estimates with those obtained by solely analysing the short numerical-relativity waveforms covering the post-merger stage of BNSs, following the procedure described in (Zimmerman et al., 2019).

We inject signals $h(\theta_{true})$ with source parameters $\theta_{true}$, that include HMs, in zero noise and estimate the source parameters using waveform templates $h(\theta)$ that omit and include HMs. We compute the posterior Bayesian probability of the parameters $\theta$ as

with $\pi(\theta)$ denoting the prior probability of the parameters $\theta$ and ${\cal{L}}(\theta|h(\theta_{true}))$ denoting their likelihood. As usual, the latter is defined as the standard frequency-domain likelihood for GW transients (Finn, 1992; Romano & Cornish, 2017)

where N runs over the different detectors of our network. As we discuss later, we work with two and three-detector networks. As usual, $(a|b)$ represents the inner product (Cutler & Flanagan, 1994)

where $\tilde{a}(f)$ denotes the Fourier transform of $a(t)$ and ^∗ the complex conjugate. The factor $S_{n}(f)$ denotes the one-sided power spectral density of the detector. In this work, we consider a network of detectors, all with noise sensitivity equivalent to that of the proposed 2.5-generation Neutron star Extreme Matter Observatory NEMO (Ackley et al., 2020). We choose a lower frequency cutoff of $f_{min}=20$Hz and a sampling frequency of $16$kHz so that $f_{max}=8$kHz. The NEMO detector has a proposed sensitivity similar to Cosmic Explorer and Einstein Telescope in the kHz regime; we could equally use those third-generation detectors and achieve similar results for the late inspiral and post-merger, albeit with larger signal-to-noise ratios for the full signal given the better low-frequency ($\lesssim 500$ Hz) sensitivity.

In all cases, we assume standard prior probabilities for the sky-location, source orientation and polarisation angle, together with a prior uniform in co-moving volume and an uniform prior on the time-of-arrival, with a width of $0.2\rm\ s$, centered on the true value. For our analyses making use of 128s-long Phenom waveforms, we impose uniform priors on the individual masses $m_{1,2}\in[1,2]M_{\odot}$ and on the components of the individual spins along the orbital angular momentum $\chi^{z}_{1,2}\in[-0.15,0.15]$. Since numerical-relativity waveforms, however, are only produced for a discrete set of intrinsic parameters, we assume the masses and spins to be known in this case. We find this is a reasonable assumption as the individual masses and the effective spin parameter (Santamaria et al., 2010) $\chi_{eff}=(\chi^{z}_{1}m_{1}+\chi^{z}_{2}m_{2})/(m_{1}+m_{2})$ are very well measured from the long inspiral. As an example, with a triple-detector network and using HMs, we determine the total mass, chirp mass and effective-spin parameters of our face-on unequal-mass source with respective uncertainties of $<1\%$, $<0.01\%$ and $<0.015$ at the $68\%$ level. We perform our parameter inference runs with the software Bilby (Ashton et al., 2019; Smith et al., 2020), sampling the parameter space with the algorithm CPNest (Veitch et al., 2020).

We consider three network configurations. The first (denoted HV) assumes two detectors with the location and orientation of Advanced LIGO Hanford and Virgo. Such a network has each detector sensitive to one of the two independent GW polarisations. The second network (HL) assumes LIGO Hanford and Livingston location and orientations, almost anti-aligned, so that both detectors are sensitive to roughly the same GW polarisation, missing the other one. Finally, we consider an HLV network that is sensitive to both GW polarizations and can pinpoint the sky-location of the source.

The accuracy of our distance measurement is of course limited by the loudness of the injected signals in our detector network, quantified by the optimal network signal-to-noise ratio (SNR), given by

which is inversely proportional to the source luminosity distance $d_{L}$. As we show in Appendix III, for face-on cases we obtain optimal SNRs of $\sim 190$ for HL and HLV networks and $\sim 140$ for HV. In contrast, for the weaker edge-on cases we obtain respective values $\rho_{opt}\in[30,65]$ depending on the mass ratio and network considered. For comparison, the SNR of GW170817, whose source was rather face-on, was only $\simeq 32$. We note that while in this study we restrict to sources placed at distances $d_{L,true}=40$Mpc, our results can be extended to different reference distances. In particular, in our SNR-regime, given fixed source and detector network, the uncertainty of our distance estimates roughly depends on the optimal SNR (and therefore on the reference distance) as $\Delta(d_{L})\propto 1/\rho_{opt}^{2}\propto d_{L,true}^{2}$ (Cramer, 1999; Rao, 1992; Li, 2013). In addition, uncertainties in the estimated value of the Hubble constant $H_{0}=cz/d_{L,true}$ grow linearly with distance.

We choose four target sources with total total mass $M=2.75M_{\odot}$, mass ratios $q=1$ and $q=1.5$, and oriented both face-on ($\iota=0$) and edge-on ($\iota=\pi/2$). Signals emitted in these two angles differ in three aspects: morphology, polarisation, and loudness. Face-on signals contain solely the $(\ell,m)=(2,2)$ mode, are circularly polarised (i.e., both $h_{+,\times}$ have the same amplitude), and are louder than edge-on signals. On the contrary, edge-on signals have contributions from higher-order emission modes that confer a richer structure, but are weaker in amplitude than face-on ones (see Fig. 1). In addition, edge-on signals contain only one of the two polarisations. Consequently, it has been shown that measuring the ratio of the two polarisations is key to infer correctly the inclination of the source, provided that the detector network can observe both polarisations (Usman et al., 2019). For each of these sources, we consider two EOSs, namely SLy and MS1b, which we assume to be known. Different EOSs trigger post-merger HMs in different ways, varying the accuracy of the distance estimate. Finally, we would ideally consider a wide range of distances and sky locations for all of our target sources. However, given the extreme computational cost of our parameter inference runs, and to allow for a direct comparison of our results, we place all of our sources at the same distance and sky location. For the former, we considered that the most reasonable choice was a value of $d_{L}=40$Mpc, consistent with that of GW170817, the only conclusive BNS observed to date through GWs. Finally, we placed all of our target sources at the same, arbitrary sky location.

Figure 2 shows the two-dimensional posterior distributions for the luminosity distance and inclination of face-on oriented binaries with mass ratios of $q=1$ (left) and $q=1.5$ (right), using an HLV detector configuration. The contours denote the $90\%$ credible regions and the legend provides median estimates with symmetric $68\%$ credible intervals. For unequal masses, and using 128s-long Phenom waveforms, the omission of HMs in the templates $h(\theta)$ leads to a biased estimate of $d_{L}^{\text{Phenom,22}}=36.2^{+2.2}_{-2.5}$ Mpc. Their inclusion corrects this bias and reduces the uncertainty by $\sim 70\%$, yielding a distance measurement of $d_{L}^{\text{Phenom,HM}}=39.4^{+0.4}_{-0.9}$. For equal masses, the impact of HMs is significantly milder, as the (usually strongest) odd-$m$ emission modes are suppressed (Pan et al., 2011, 2014; Blanchet, 2014). This leads to a biased estimate of $d_{L}^{\text{Phenom,HM}}=37.9^{+1.4}_{-2.3}$ even if HMs are included. We also note that this remains true even if one would assume the intrinsic source parameters (masses and spins) to be known. This situation changes, however, when using $80ms$-long numerical-relativity waveforms that can account for the rich post-merger signal morphology. For both mass ratios, and considering a SLy equation of state, results are better than those making use of 128s-long Phenom waveforms omitting HMs. Moreover, for equal-mass, the estimate even improves on that including HMs, yielding a non-biased estimate $d_{L}^{\texttt{SLy}}=38.6^{+1.8}_{-2.7}$ with smaller uncertainties.

We combine the distance estimates obtained using numerical-relativity waveforms with those obtained analysing Phenom waveforms restricted to frequencies not covered the former. To this, we multiply the respective posterior distributions for the the distance, dividing by one instance of the prior (Zimmerman et al., 2019) ²²2Note that this ignores that stronger constraints can be obtained for the inclination and the sky-location by combining both measurements, making our results rather conservative.. We obtain joint estimates of $d_{L}^{\text{joint}}=38.6^{+0.9}_{-1.3}$ Mpc for the $q=1$ case and $d_{L}^{\text{joint}}=39.4^{+0.4}_{-0.7}$ for the $q=1.5$ case. The reason behind this improvement is that matter effects arising during the post-merger of BNSs trigger HMs, helping to break the degeneracy between distance and inclination and even allowing to measure the azimuthal angle (see Appendix II). For the unequal-mass case, these only add a small contribution with respect to the integrated effect of the HMs during the 128s of the signal. For $q=1$, however, odd-m modes are suppressed during the inspiral and are only triggered during the post-merger of our numerical-relativity simulations due to an effect known as one-armed spiral instability or 21-mode instability (East et al., 2016; Radice et al., 2016; Lehner et al., 2016) (see Suppl. Material). As a consequence, the inclusion of the post-merger emission can have a large impact. In our study, the distance estimates are much better for SLy compared to MS1b, due to the stronger HM emission.

Finally, while we have discussed results assuming an HLV detector network, in the Suppl. Material we show results for all studied network configurations. In all cases we obtain similar results, both qualitatively and quantitatively.

Previous work has shown that the degeneracy between distance and inclination can be broken by computing the ratio $h_{\times}/h_{+}$, as this evolves from $1$ to $0$ as the inclination varies from face-on to edge-on (Usman et al., 2019). We find that for face-on binaries, all HL, HV, and HLV configurations yield almost equivalent distance measurements despite the differences in signal loudness across the network, so that HMs have a much larger impact on the measurement than the polarisation ratio. In contrast, in Fig. 3 we show that for edge-on cases it is key not only to access both signal polarisations but also having a third detector that can pin-point the sky-location of the source. For an HLV configuration, unbiased estimates with uncertainties lower than $4\%$ are obtained regardless of the usage of HMs, while biased estimates are obtained using both two-detector configurations. The reason is that such configurations cannot pin-point the sky location of the source using timing information. This way, Bayesian inference places the source at those patches of the sky, consistent with the two-detector timing, where the detector network is most sensitive to, biasing the distance toward large values. We obtain identical qualitative results for the $q=1$ case, and when analysing our 80 ms-long numerical-relativity simulations.

Combining our distance estimates with simulated redshift estimates, we can infer $H_{0}$ via $H_{0}=cz/d_{L}$, with $c$ the speed of light. We assume redshift estimates consisting of Gaussian posterior distributions centered at a value of $z_{0}=0.00897$, corresponding to $d_{L}=40$ Mpc in a $\Lambda$CDM cosmology with Hubble parameter $H_{0}=67.9$ km s^-1 Mpc^-1. We assume that the redshift has an uncertainty with standard deviation $\Delta z$ ranging from a realistic value of $10^{-3}$ consistent with that for GW170817 (Abbott et al., 2017c) to an improved value of $10^{-4}$ which would require higher resolution spectrograms and also the possibility to measure the internal motion of sources within individual galaxies, e.g., (Davis et al., 2019).

Figure 4 shows the $H_{0}$ estimates derived from the distance measurements of our face-on BNSs with mass ratios $q=1$ and $q=1.5$ as a function of $\Delta z$ using an HLV network. Once again, we quote results in terms of symmetric $68\%$ credible intervals centred at the median value. For unequal mass, and for $\Delta z=10^{-3}$, the omission of HMs never leads to biased estimates. Their inclusion, however, leads an important improvement from $H_{0}=73.8^{+9.0}_{-9.8}$ to $H_{0}=68.3^{+7.5}_{-7.5}$. More spectacularly, for the most optimistic $\Delta z=10^{-4}$ HMs improve the measurement from $H_{0}=73.8^{+5.3}_{-4.1}$ to $H_{0}=68.0^{+1.6}_{-1.0}$, enabling a $2\%$-level measurement. This shows that with 2.5G detectors, together with improved detector calibration and waveform models, $H_{0}$ estimates will be limited by EM redshift measurements and not by GW distance ones. Conversely, if HMs are omitted, a significant reduction of $\Delta z$ will not translate into an improved $H_{0}$ estimate. Moreover, we find that $H_{0}$ estimates would be biased when $\Delta z\lesssim 0.5$. Consistently with the previous section, the inclusion of the post-merger emission in the analysis does not lead to any relevant improvement.

As expected, the situation is different for equal-mass cases. For these, the inclusion of post-merger effects is crucial to obtain visible improvements in the $H_{0}$ measurement. For $\Delta z=10^{-3}$ we obtain a mild improvement from $H_{0}=71.2^{+8.6}_{-8.3}$ to $H_{0}=69.6^{+7.8}_{-7.7}$. When $\Delta z$ is reduced to $\Delta z=10^{-4}$, the HMs present in the post-merger emission allow for an estimate $H_{0}=69.3^{+2.5}_{-1.5}$, with uncertainties at the $\simeq 4\%$ percent level, while their omission doubles the uncertainty and biases the measurement.

Almost identical results hold for the HL and HV networks, as shown in the Suppl. Material. For the (weaker) edge-on binaries, we find percent-level measurements are possible using the HLV network regardless of the usage of HMs.