key: cord-0298378-adpdpilo authors: Datta, Sayak; Phukon, Khun Sang; Bose, Sukanta title: Recognizing black holes in gravitational-wave observations: Challenges in telling apart impostors in mass-gap binaries date: 2020-04-13 journal: nan DOI: 10.1103/physrevd.104.084006 sha: a1ef0dd6cf2708c3da7560f51dca238bf02d8f9e doc_id: 298378 cord_uid: adpdpilo We study how by careful monitoring of the presence or absence of tidal deformability (TD) and tidal-heating (TH) in the inspiral signal of compact object binaries in ground-based gravitational wave (GW) detectors, one can test if its components are black holes or not. The former property (TD) is finite for neutron stars but vanishes for black holes (in general relativity), whereas the latter is finite for black holes and negligible for neutron stars, and affects the GW phase evolution of binaries in a distinctly different way. We introduce waveform parameters that characterize the strength of tidal-heating, and are zero when there is no horizon. We develop Bayesian methods that use TD and TH for distinguishing the presence or absence of horizons in a binary. This is timely owing to several claims that these stellar-mass objects, especially, with masses heavier than those of neutron stars, may not have a horizon but may be black hole mimickers or exotic compact objects (ECOs). It is also astrophysically important to have the tools to test the presence or absence of horizons in mass-gap binaries and, thereby, help detect the heaviest neutron star or the lightest black hole. A proper accounting of tidal-heating in binary waveform models will also be critical for an unbiased measurement of characteristics of the equation of state of neutron stars in GW observations of binaries containing them -- or even to probe the existence of ECOs. We show that purely based on GW waveforms it will not be possible to discern binary horizons in the mass gap in Advanced LIGO, Virgo and KAGRA detectors unless the binary is within a few tens of Mpc. However, third generation ground-based detectors will be able to do so for binaries a few hundred Mpc away. In recent times, the discovery by LIGO-Virgo detectors of several compact binary coalescences (CBCs) has ushered in the era of gravitational wave (GW) astronomy [1, 2] . The LIGO-Virgo Collaboration also observed the binary neutron (BNS) star merger GW170817 [3] . These observations provided a fillip to tests of general relativity (GR) in the strong-field regime [4, 5] ; e.g., stringent bounds on the mass of the graviton and violations of Lorentz invariance have been placed [6] [7] [8] . Significantly, it has also become possible to test the nature of the compact objects in binaries. The deduced compactness of the components has led to the conclusion that they are either black holes (BHs) or neutron stars (NSs). In the case of GW170817, radius measurements were made [9, 10] that strongly disfavor them as BHs. A similar claim may be posited for the other BNS contender GW190425 [11] . However, for the other LIGO-Virgo binaries (which are much heavier than GW170817 or GW190425) [1] , it remains to be conclusively proven that their components are indeed BHs of GR and not, say, some exotic compact objects (ECOs) [12] [13] [14] . If the binaries show up with measured components masses in the mass-gap ∼ 2 − 5M [15] , then it poses the immediate challenge of determining whether the component(s) with mass(es) in the gap are NSs or BHs. Either occurrence will be significant, for it will either raise the maximum known mass of a NS or lower the minimum known mass of a BH. These issues make it imperative that methods be devised to discern compact objects with horizon from those without. In this work we study if the presence of horizon can be detected in binaries in the mass-gap by LIGO-Virgo. We also include in our study the mass range 1 − 2M where neutron-star masses commonly occur. Apart from NS and BH, ECOs may also occur in the same mass range. Multiple models of ECOs have been proposed. These include Planck-scale modifications of BH horizons [16, 17] , gravastars, [18] , and boson stars [19] -to name a few. In light of such proposals, it becomes necessary to devise strategies to tell them apart from BHs. Several tests have been proposed to probe the black-holeness of the compact objects in a binary. Distinguishing binary merger remnants from BHs in the postmerger phase using echoes has initiated rigorous modelling and search for those features in GW data [13, [20] [21] [22] [23] [24] [25] [26] [27] . Measurement of tidal deformability (TD) [28] [29] [30] [31] and spin-induced multipole moments [32] [33] [34] [35] [36] from the late inspiral can also be used to test black-holeness. In this paper, we expand on past work to study how diffi-cult it is to perform a horizon test by using GWs emitted during the inspiral phase of binary coalescences. For this purpose we include terms in the binary waveform phase beyond the point-particle ones that arise due to the material characteristics of the objects or the presence of horizons. In particular, we introduce two new best measured horizon parameters for stellar-mass binaries in groundbased detectors. Their precise measurement in binary observations is useful in probing the existence of horizons in those systems. Owing to their causal structure, BHs in GR are perfect absorbers that behave as dissipative systems [37] [38] [39] [40] . The defining feature of a BH is the presence of its horizon, which is a null surface and a one-way membrane. It is due to the presence of the horizon that a BH in a binary absorbs energy and angular momentum from the orbit. This phenomenon is called tidal-heating (TH) [41] [42] [43] . Energy loss via TH backreacts on the binary's evolution, resulting in a shift in the phase of the GWs emitted by the system. Therefore, the absence of a horizon -or any kind of change in the near horizon structure that modifies this absorption -will leave its imprint in the phasing of GWs emitted. A careful observation thus has the potential to measure these differences in the GW phase. Indeed, TH has been proposed for probing the presence of horizons along with the existence of higher dimensions and quantum effects at horizon scale [44] [45] [46] [47] [48] . Its importance in identifying horizons of intermediate-mass and supermassive compact objects has been examined for the space-mission LISA [30, 33, 49] . In the current work, we study its usefulness for stellar mass binaries -of the type observable by ground-based GW detectors like LIGO, Virgo, and KAGRA [50] . The TH of a black hole or any other star can be expressed in similar mathematical forms if the viscosity coefficient (η) of a BH is identified with its mass [51] . For NS, one has η NS ∼ Since the correction in GW phase due to TH is proportional to η, for an NS that correction is 10 orders of magnitude smaller than BH [51] . While this distinction presents an interesting prospect for observational exploitation, as we show here the magnitude of TH for binary black holes (BBHs) remains small and is useful for discerning the presence of horizons for very large M (as for EMRI central objects in LISA) or for strong signals. This implies that for stellar-mass BHs, detection of TH in LIGO-Virgo will require the binary to be within tens of Mpc, as shown below. While the occurrence of such a golden binary is not impossible, GW observations to date rule it as improbable. Nevertheless, for completeness of TH analysis, we examine this case in this work. For more realistic BBH distances, a detection of TH and its utilization for discerning horizons will have to wait for third generation detectors. The formalism initiated here for accounting for TH will be relevant for those detectors as well. Another property of compact objects that leaves an imprint on GWs is tidal deformability (TD). A body immersed in an external tidal field, such as due to a binary companion, experiences an induced quadrupole moment. That moment is proportional to the tidal field, and the proportionality factor is the tidal deformability (λ). This tidal deformation in turn affects the binary's orbital motion and the emitted GWs. The GW phasing carries an imprint of the dimensionless tidal deformability (Λ i ≡ λ i /m 5 i ) of the two masses m i=1,2 [52] . Material bodies, such as NS, have substantial Λ i values [9, 11] , but black holes have a vanishing value [53, 54] . 1 Hence, using appropriate modeling it is possible to measure Λ i and probe the properties of the bodies. More than TH, it is TD that we find to have a dominating influence in recognizing the absence of horizons in a stellar-mass binary, particularly, when the components masses are around 1 − 2M . TD decreases with mass, and above this range is vanishingly small for realistic neutron star equations of state. Mass-gap objects can be as heavy as ∼ 5M . This is the reason we analyze binaries with component masses between 1 − 5M . Since TD has little influence above 2M it is left to TH to help recognize the presence of horizons. We find that it is highly improbable to do so for binaries with component masses in the range 2−5M in the current generation of detectors. It has been shown that in GR the Love number vanishes for BHs, but not for other compact objects like NSs [54, 57, 58] . However, recently it has been suggested that the Love number can be nonzero for nonaxisymmetrically perturbed rotating BHs [59] . In the current work we take that the tidal deformability of all BHs is zero. Thus, our results may need to be revisited depending on how this matter gets resolved. We begin by studying in Sec. II the TH terms that appear in the GW phase of a binary. There we identify two horizon parameters that are best measured for stellarmass binary signals in ground-based detectors. There we also show how the spin-induced quadrupole moment and tidal deformability of the binary components influence the waveform. In Sec. III we develop the method for weighing the evidence in data for the presence or absence of horizon, utilizing the aforementioned horizon parameters and phase terms in a Bayesian formalism. In Secs. IV and V we implement this formalism on a large population of simulated binary signals in noisy data simulated with Advanced LIGO and Advanced Virgo noise. We conclude with a discussion on future prospects in Sec. VI. Consider a compact binary with component masses m i (i = 1, 2), total mass m = m 1 + m 2 , and mass-ratio q = m 2 /m 1 , with m 2 ≤ m 1 . Let the dimensionless component spins be χ i . Under the adiabatic approximation the orbital evolution of the binary can be quantified in the post-Newtonian formalism with reasonable accuracy, especially, when it is far from merger [60] . The dynamics of the system is governed by energy and angular momentum loss from the orbit. Usually it has a contribution arising from taking the components as point particles (PP) and another one originating from their finite size. The latter can be decomposed into two parts, (i) tidal deformation of each component due to the gravitational field of the other and (ii) the amount of energy absorbed by individual components from the orbit, namely, tidal heating. The dynamics of the system and, therefore, the emitted GW depends on all of these contributions. Hence, the Fourier transformed GW waveform can be written as where f is the instantaneous GW frequency and A(f ) is the frequency-dependent amplitude. The phase terms -Ψ PP , Ψ TD , and Ψ TH -are the phase contributions arising from the point-particle approximation, TD, and TH, respectively. Since GW absorption is negligible for matter [51] , it is reasonable to exploit evidence of TH in binary waveforms to discern the existence of horizons [30, 33] . This expectation led us to introduce the horizon parameter H for extreme mass-ratio inspirals (EMRIs) that LISA may observe [33] . Till now horizon distinguishability employing TH has been addressed primarily for LISA sources, such as EMRIs and supermassive BH binaries. However, even in the case of supermassive BH binaries, it is the combined tidal heating of both binary components is what has been employed, which ignores the possibility that not both components may have horizons (or lack them) [30] . Such an approach is reasonable for initial forays in this subject but, in general, different values of H need to be considered. Here we apply the formalism to binaries with similarly massive components primarily to target the LIGO-Virgo population of stellar-mass binaries. As it is a broad subject, we keep the studies with third-generation for the future. For a near-equal-mass binary we define horizon parameters for each component, (H 1 , H 2 ), such that the value of H i is 1 (0) when the ith component has a horizon present (absent). In the case of circular orbits, the flux of energy at the horizon can be expressed as a PN expansion [39, 40, [61] [62] [63] [64] [65] . Since TH signifies presence of horizon, we multiply the energy flux absorbed by each component with the corresponding H i . In the case of partial absorption, one has 0 < H i < 1. Therefore, the absorbed flux is (2) where ν = m 1 m 2 /m 2 is the symmetric mass-ratio, v is the orbital velocity, andŜ i andL N are the unit vectors along the directions of the ith spin and the orbital angular momentum, respectively. New waveform parameters characterizing TH The horizon parameters H 1,2 appear in the GW phase in terms that also include mass and spin factors. This makes them degenerate with those parameters, in that it is more practical to measure the following effective observables instead of H 1,2 : These are analogous to the effective spin parameter χ eff that was introduced [66-68] to characterize spinning compact binary waveforms: While the spins of the individual binary components are themselves difficult to measure (like H 1,2 here), their combined impact on the waveform phase, captured by χ eff , lends itself to more precise measurements. Dependencies of H ef f 5 and H ef f 8 on component spins are shown in Fig. 1 and Fig. 2 , respectively. If the system is a binary black hole (BBH), as long as any one of the component spins is finite both H ef f 5 and H ef f 8 will be nonzero. By contrast, for the same spins a horizonless binary would have both H ef f 5 and H ef f 8 vanish. Therefore, it is easiest to discern between the presence and absence of horizons in BBHs that have at least one component with sufficiently large spin. On the other hand, when both component spins of a BBH tend to zero, one has H ef f 5 → 0 but H ef f 8 = 0; see the inset in Fig. 2 . Therefore, in the low-spin limit H ef f 8 emerges as a discriminator for the presence or absence of horizons. Here the measurement is helped for small mass-ratio (q), which ensures large H ef f 8 . It is important to note that our choice of waveforms, based on the stationary-phase approximation (SPA) [69] , is for illustrative purpose, essentially as a proof of principle that the method proposed here is promising for identifying binary components with horizons from those without. For making such classification in real data, it will likely be important to use more accurate templates, such as those based on the EOB-NR formalism [70] [71] [72] . We will present those results in future. Having said that, our choice of SPA-based inspiral waveforms is a reasonable one for illustrating the power of this method for the systems studied here. We deduce the GW phase involving TH by using Eq. (2.7) of Ref. [73] (see [74] for the details). We find the phase shift due to the associated horizon absorption to be Note that H ef f 5 and H ef f 8 arise at different PN orders in the phase. The two bodies in a coalescing compact binary can have spin. In case of a nonzero spin a body would develop a spin-induced quadrupole moment. The leading order contribution arises due to the mass quadrupole moment (M 2(i) ) of both bodies, (i = 1, 2), at 2 PN order. If the bodies are BHs, then M 2(i) = −χ 2 (i) m 3 (i) . If they are NSs or ECOs that moment may be modified as Measuring the quadrupole moment κ from observations can be used to probe the nature of the compact objects [32, 33, 75, 76] . Since in a binary the quadrupole moments κ i of both the bodies contribute at the similar order, they are degenerate. Usually, a combination of κ 1 and κ 2 are used for the measurement [32] . Here we define a new effective parameterQ as follows: . Then the phase can be expressed as [77] Note that once the phase has been expressed in terms of Q, it is not necessary forQ to be limited to Eq. (6): It is straightforward to incorporate other models of M 2 into it, such as for boson stars [25, 78] . We will use Eq. (7) in our modeling of phase due to nonzero quadrupole moment. In Refs. [79, 80] , observational constraint on spin induced quadrupole moment has been found. Individual measurement suffers from broad posterior distribution, pointing towards low measurability with current detectors. Our results below are consistent with these observations. In the presence of a GW signal strain h(t, θ θ θ), characterized by parameters θ θ θ, the detector strain time-series can be modeled as d(t) = n(t) + h(t, θ θ θ), where n(t) denotes the detector's noise. In the presence of a GW signal h(t, θ θ θ), described by a model H, the likelihood of the data is [69] : under the assumption of Gaussian and stationarydetector noise. The angular bracket in Eq. (8) defines a noise-weighted inner product between two real timeseries a(t), b(t), and is given as where S n (f ) is the one-sided power spectral density (PSD) of the detector noise, and f low and f high are the low-frequency cutoff and high-frequency cutoff, respectively [69, 81] . Using the inner product, one can also define the signal-to-noise ratio (SNR) ρ for the template h(t, θ θ θ) as where σ = h|h is the template normalization. We will assume that noncolocated detectors on the globe have uncorrelated noise; hence, the combined likelihood is given as [82] , where d ∈ {d 1 , d 2 , · · · , d N } represents combined data from all N detectors. Using the coherent network likelihood function, posterior probability density can be written as where P (θ θ θ|H) is the prior probability density function or prior of the parameters θ θ θ. In the denominator, P (d|H) is the marginalized posterior probability density over all parameters θ θ θ, and is also known as the evidence for the model H. The evidence P (d|H) serves as a normalization constant of the posterior probability for H. The evidence computed for two competing models or hypotheses can be used to determine which one is favored by the data. In this work, we compute Bayes factors for simulated signals to compare two hypotheses, namely, 1. The horizon hypothesis H BBH : Signal carries imprints of horizon absorption and spin-induced quadrupole moment, 2. The no-horizon hypothesis H BNS : Signal has no imprint of horizon absorption, but has TD and spininduced quadrupole moment. In Bayesian model selection, we compute the Bayes factor, If the Bayes factor is greater than some preset threshold, i.e., BF > BF Th then the hypothesis H BBH is preferred over the other hypothesis H BNS in the data. Moreover, we use the Dynesty sampler [83] , as implemented in the Bilby package [84, 85] , to compute the posterior probability densities for our simulated signals. We use a likelihood function marginalized over time t c and phase φ c at coalescences of binaries [86, 87] and distance d L [88, 89] , thus removing the need for sampling those parameters without affecting the posterior probability densities in the parameters of interest. The posterior probability densities for these parameters can be reconstructed analytically from the full set of posterior samples [87] . The posteriors of some of the parameters for the H BBH hypothesis are shown in Fig. 3 . To compute them, we considered the signal integration in a frequency range such that it ends at f ISCO , while the duration of the signal is 16s. f ISCO is the instantaneous GW frequency at the ISCO of the binary [90, 91] . In practice, it may be possible to begin the signal integration at a frequency as low as 10Hz, which is what aLIGO design targets. Similarly, when waveform modeling is available to accurately incorporate TH beyond the ISCO, the upper frequency cutoff will also be raised. Both these changes will improve parameter estimation as well as Bayes-factor based model discrimination. Before setting up signals simulation for BF computations, it is worthwhile to examine through computationally inexpensive, even if approximate, means how precisely the horizon parameters would be measurable in mass-gap binaries. Such a computation is afforded by the Fisher information matrix (FIM), as defined below. We estimate how large the noise-limited errors are of the horizon parameters ϑ, by modeling the measured values after the maximum likelihood estimators MLEs [92] . Owing to noise, the MLE will fluctuate about the respective true values, i.e.,θ = ϑ + δϑ, where δϑ is the random error. The extent of these fluctuations is estimated by the elements of the variance-covariance matrix, γ ab = δϑ a δϑ b [92] , which is bounded by the signal via the Cramer-Rao inequality, namely, where Γ is the FIM: Above, ∂ a is the partial derivative with respect to the parameter ϑ a . Therefore, ∆ϑ a ≡ ( δϑ a δϑ a ) −1/2 = Γ −1/2 aa gives the lower bound on the root-mean-square error in the estimate of ϑ a . The two are equal in the limit of large SNR [92] . The error estimates listed here are the ∆ϑ a obtained from the FIM. When one computes Γ ab for the binary parameters one typically finds that its offdiagonal terms are nonzero, which implies that there are covariances among the parameter errors. It is, however, possible to mitigate those covariances for a different set of parameters. In the twodimensional parameter subspace of (H 1 , H 2 ), we find that (H ef f 5 , H ef f 8 ) are such parameters. One can also use FIM to deduce errors (∆H ef f 5 , ∆H ef f 8 ) in the new horizon parameters for our binaries of interest. This is how we estimate that for a mass-gap binary at a distance of 10Mpc to a few tens of Mpc, it is possible to measure (H ef f 5 , H ef f 8 ) to a few tens of percent in a three detector LIGO-Virgo network with the aforementioned noise PSD. A similar FIM calculation for the third generation detector Einstein Telescope shows that the same measurement precision is achievable even when the same mass-gap BBH is pushed out to a few 100 Mpc. As mentioned above, in spite of the weak effect of TH in mass-gap binaries in current detectors, for the completeness of the waveforms used in our simulations we continue to retain the TH terms in their phases. The impact of those terms for third generation detectors and binaries not limited to the mass-gap will be studied elsewhere. The distributions and ranges of parameter priors of the simulated binary waveforms used in our Bayesian model selection studies are listed in Table I . The possible values of H ef f 5 and H ef f 8 are shown in Fig. 1 and Fig. 2 , respectively. In Fig. 3 , we show the posterior probability distributions of various parameters of a BBH injected signal obtained from a Bayesian analysis. The luminosity distance (d L ), chirp mass (M), mass ratio (q), and effective spin (χ ef f ) are well measured. The estimation recovers the injected values. Comparatively H ef f 5 and H ef f 8 are poorly measured. Although we recover the injected values, and the posterior is certainly different from the flat prior, the error is large. This is expected as TH is a higher-order effect. We now quantify how successfully one can discriminate between a BBH signal from a BNS one in noisy data. For this signal model selection test we simulated a population of 1250 binaries, which are distributed uniformly in comoving volume between 50Mpc to 250Mpc. Component masses were taken to be m 1,2 ∈ [1 − 5]M and spins chosen to be aligned or antialigned with the orbital angular momentum, and with dimensionless magnitude χ 1,2 ∈ [0, 0.9]. For model selection we constructed two families of waveforms -both for signals (for adding in simulated noisy data) and templates (for matchedfiltering that data) -namely: (a) TaylorF2 (TF2), modified with TD contribution (TidalTF2) for representing horizonless components with nonzero TD. Here, the GW phase is devoid of any contribution from H ef f 5 or H ef f 8 ; (b) HeatedTaylorF2 (HTF2), which is TF2 but with additional phase terms arising from TH, as described in Eq. (4). We have included the effect of the spin-induced quadrupole moment appropriately in both cases via the phase term in Eq. (7) . We used the Akmal, Pandharipande, and Ravenhall (APR) equation of state (EOS) [93] for this purpose [94] to model the new effective param-eterQ introduced in Eq. (7) . The injected values ofQ have the range [−1.60, 0]. From Ref. [94] we constructed the values of κ i of the ith body of mass m i and spin χ i . From these values we find the corresponding value ofQ using Eq. (6), which is used for injection. Using the aforementioned waveform models we performed simulated signal injection studies in simulated [95] . To keep computational costs manageable we limited all our signals (and the filtering and parameter estimation) to only 16sec, and till the innermost stable circular orbit (ISCO). In one study, the sources are taken to be CBCs, with BHs as components. Hence, the injected waveform used is HTF2. We then performed a Bayesian analysis to measure the parameters of these sources with both TidalTF2 and HTF2 templates and compared the natural-log of their Bayes factor ln BF, utilizing the definition in Eq. 13, for the same "horizon" injections to test if such an analysis has the power to identify the true signal model. In Fig. 4 we plot the ln BF with respect to χ ef f in the x-axis and M in color. Each point in this figure represents an HTF2 injection. The fact that for a large majority of them the values of ln BF are positive, suggests that for this injection set model selection strongly favors horizon injections. In an ideal case, all of the points should be above the ln BF = 0 line. Deviation from this expectation for a minority of the injections is due to their low SNRs or similarity of their signals with the TidalTF2 waveforms for the same parameters (as will be explored in more detail below). With longer duration waveforms and higher SNRs this result should get somewhat better, for lower masses. The origin of the high values is likely due to a combination of TD, quadrupole moment, which will be discussed below. B. Assessing the statistical significance of the horizon discriminator In the preceding section, we found that barring a small subset the model selection returns positive ln BF values for the injected sources, which is tantamount to saying that the observations favor the true signal model, namely, HTF2 here. However, for a small subset, with nega- tive ln BF values, the wrong signal model (TidalTF2) is preferred. This raises the possibility that the opposite can also happen, i.e., some TidalTF2 injections, searched with both types of templates, may return ln BF values favoring the HTF2 signal model. As with any statistical analysis, it also becomes important to interpret quantitatively the probability with which the nature of those sources will be identified correctly. BBH injection studies enable one to do precisely that. However, it is also important to assess the probability with which the nature of that source will be misidentified. For example, if the value of ln BF turns out to be 5 for a BBH signal, it is important to interpret that value in terms of how probable it is to be identified correctly as of BBH origin (the true hypothesis, H BBH ) and incorrectly from a BNS (the wrong hypothesis, H BNS ), in noisy data. As we discuss next, the former probability can be assessed from the above study of BBH injections, whose ln BF values form the foreground distribution p (ln BF | H BBH ), which is the probability distribution of the ln BF values given that the H BBH hypothesis is true, i.e., the (injected) signals belong to the HTF2 model. On the other hand, to assess how probable it is for the TidalTF2 signals to be misidentified as HTF2, we also study injections of horizonless signals generated using the TidalTF2 waveform model; these ln BF values form the background distribution p (ln BF | H BNS ) when the hypothesis being tested for a detected event is that it is from a BBH. To obtain the foreground distribution corresponding to HTF2 signals, we compute the BF values for HTF2 injections using Eq. (13) . We plot the distribution of the ln BF for these values in red in Fig. 5 . To construct the background distribution for the same signals, we compute the BF values -but now for TidalTF2 injectionsusing Eq. (13) . The distribution of the ln BF for these values is plotted in black in the same figure. The samples of ln BF from the foreground and background distributions are used to estimate the efficiency with which BBHs can be identified in GW events and assign a statistical significance to each identification. In Fig. 5 , we show the estimated foreground and background distributions of the ln BF for a subpopulation of the injected binaries discussed in Fig. 4 and the pre-ceding subsection. This subpopulation includes only those binaries that have component masses in the range m 1,2 ∈ [1−2]M . As noted earlier, the foreground distribution is constructed by injecting HTF2 and calculating the Bayes factor in favor of HTF2. The background distribution is estimated by calculating the Bayes' factor in support of HTF2 for the TidalTF2 injections. We use the background distribution p (ln BF | H BNS ) to compute the false-detection-probability (FDP) of a BBH claim, given that a BNS signal is actually present in the data. 2 The FDP is computed from p (ln BF | H BNS ) for a measured ln BF-value ln BF measured as follows: If the FDP is sufficiently low, then it is less likely that the event is consistent with the H BNS hypothesis. Often the FDP values are converted to equivalent significance levels, e.g., n σ deviation of a Gaussian random process. From the background distribution of ln BF, we can compute the threshold Bayes' factor ln BF threshold corresponding to a certain statistical significance. In our analysis, each GW signal is injected in a 16second-long simulated colored-Gaussian data of the two Advanced LIGO detectors. In Fig. 5 , the foreground distribution p (ln BF | H BBH ) and the background distribution p (ln BF | H BNS ) are shown in red and black colors, respectively. The blue vertical line denotes the threshold value ln BF threshold corresponding to 2σ significance. Above that 2σ threshold, the areas under the foreground and the background curves are ∼ 0.17 and ∼ 0.04, respectively. Therefore, around 17% of BBH signals will have ln BF greater than that threshold, and will be correctly identified as BBHs with a significance at the 2σ level. However, there is a 4% chance that a signal from a BNS will be mischaracterized as a BBH at the same level of significance. We notice that there is a significant overlap between the foreground and background distributions, which is expected to decrease somewhat with longer waveforms and more sensitive detectors. In Fig. 6 we perform the same exercise for heavier component masses, namely, m i > 2M . In this mass range interestingly we find that irrespective of the type of the source injected, model selection prefers HTF2. This is because for those masses neither the TD (andQ) phase terms in TidalTF2 nor the TH ones in HTF2 are large enough to induce phase difference between the TidalTF2 and the HTF2 waveforms that is significant enough to tell them apart. In fact, both waveforms are very similar to the point-particle waveform there. While with increasing total mass TH would eventually become large, nevertheless in the mass-gap it is weak enough to tell such binaries apart at realistic distances in the current generation of detectors. We further test the conclusion above with a study of any systematics that may be induced in the estimation of parameters when the wrong waveform is used to search for a signal. For this purpose, we focus on the most precisely measured binary parameter, namely, the chirp mass M. In Fig. 7 we investigate these systematics. When the injection is TidalTF2 but parameter estimation employs TidalTF2, on the one hand, and HTF2, on the other hand. We plot the values of M so measured in Fig. 7 . If there were no systematics, then the measured values should be highly correlated between TidalTF2 and HTF2 measurements. In an ideal case they should fall along the diagonal line, barring a spread owing to detector noise. In the figure we find this behavior as expected. The measured M values fall exactly on the diagonal line. This implies that the measured values are highly correlated, which implies the absence of systematics. The fact that even the presence of the TH terms in the phase of HTF2 waveforms are not able to effect any discriminatory power may not be suprising but is of special significance. While this may be a disappointing result for prospects of characterizing the nature of compact objects in the mass gap, it is important to note that this conclusion is reached with the TH terms in phase. In retrospect, this is not surprising since Fisher studies point to the same conclusion. These studies, with full waveforms, also indicate that such a distinction is possible in third-generation detectors for binaries within a few hundred Mpc. Bayesian studies with longer waveform simulations for those detectors are computationally expensive and will be pursued elsewhere. We have developed a method to search for and characterize TH in the inspiral phase of a binary. We have defined two new parameters that capture the effect of TH in the inspiral waveform. These parameters are robust enough that even partial absorption can be modeled with them -something we will pursue in detail in the future. To test for the presence of horizon we performed model selection using the Bayes factor. We constructed two sets of waveforms, one for BBHs, which incorporates TH but no TD, and the other for binaries of horizonless compact objects, which does not include TH but has TD. We also defined a new effective parameter for the quadrupole moment, namelyQ, which has been added in both waveform models appropriately. We showed that for m i > 2M it was not possible to distinguish between the two models. We did so by employing the Bayes factor in a full Bayesian analysis with simulated injections of both types of signals. We also checked our results with a Fisher analysis and found that in this mass range it will be hard to test the presence of horizon. It remains to be seen if it is possible to find better result with golden binaries. It is obvious that with in- creased SNR the error will decrease, resulting in better measurement of H ef f 5 and H ef f 8 . We can estimate the error reduction in such a case. Assuming the sources like in Fig. 3 at d L ∼ 50Mpc and taking account of the whole signal duration we can estimate the error reduction. Taking into account of these two, we find that the SNR will increase by a factor of ∼ 3.1. Hence, the estimation in such case would result in H ef f 5 ∼ 0.06 ± 0. 66 and H ef f 8 ∼ −1.88 ± 7.66. Hence, with golden binaries even with advanced LIGO detectors, it is possible to find some meaningful constraints on the sources. An immediate continuation of the current work will be to construct better and more complete waveform models than TidalTF2 and HTF2 that can be used for more precise parameter estimation and more accurate model selection for real signals in contemporaneous GW detector data. Since the waveforms used here were limited to a short part (16s) of the late inspiral phase, it is likely that utilizing more complete waveforms may improve the ability to distinguish BBH and BNSs signals. Another problem we plan to address is the challenge posed by mixed binaries (NSBH) in discerning the presence of horizons. Thirdly, future generation detectors may allow enough precision so that proper discrimination of NSBH binaries as well as the horizon parameter with intermediate values (0 < H i < 1) may be realizable, thereby, affording the possibility of probing the existence of ECOs, such as stellar-mass gravastars, boson stars, etc. [45, 49] . In the absence of systematics the measured values should be highly correlated between TidalTF2 and HTF2 measurements, and should fall along the diagonal line. In the presence of noise, however, we expect some scatter around that line -as evidenced here. As expected, the measured M values shown above are nicely clustered around the diagonal. This shows that any biases, if present, are much less significant than the statistical errors. We thank Samanwaya Mukherjee for providing useful inputs which helped us express our results better. It is a pleasure to thank N. V. Krishnendu, Andrea Maselli and Paolo Pani for useful discussions. We would also like to thank Richard Brito and Otto Hannuksela for carefully reading the manuscript and providing helpful inputs, and Bhaskar Biswas, Soumak Maitra and Niladri Paul for useful comments. We gratefully acknowledge the use of the IUCAA computing cluster, Sarathi, and the computational resources provided by the LIGO Laboratory (CIT) and supported by National Science Foundation Grants PHY-0757058 and PHY-0823459. SD would like to thank University Grants Commission (UGC), India, for financial support for a senior research fellowship. KSP acknowledges support of the Netherlands Organisation for Scientific Research (NWO). This work was done with partial support provided by the Tata Trusts. This paper has been assigned LIGO Document Number LIGO-P2000115. We would like to thank all of the essential workers who put their health at their risk during the COVID-19 pandemic, without whom we would not have been able to complete this work. VIRGO) Proc. Nat. Acad. Sci Black holes: the membrane paradigm Proceedings of the Second Marcel Grossmann Meeting ofGeneral Relativity Gravity: Newtonian, Post-Newtonian, Relativistic Progress of Theoretical and Experimental Physics Marginalisation of the Time Parameter in Gravitational Wave Parameter Estimation Advanced LIGO anticipated sensitivity curves On the criterion that a given system of deviations from the probable in the case of a correlated system of variables is such that it can be reasonably supposed to have arisen from random sampling On the problem of the most efficient tests of statistical hypotheses