Pre-DECIGO can get the smoking gun to decide the astrophysical or cosmological origin of GW150914-like binary black holes

Pre-DECIGO consists of three spacecraft arranged in an equilateral triangle with 100km arm lengths orbiting 2000km above the surface of the earth. It is hoped that the launch date will be in the late 2020s. Pre-DECIGO has one clear target: binary black holes (BBHs) like GW150914 and GW151226. Pre-DECIGO can detect $\sim 30M_\odot-30M_\odot$ BBH mergers up to redshift $z\sim 30$. The cumulative event rate is $\sim 1.8\times 10^{5}\,{\rm events~yr^{-1}}$ in the Pop III origin model of BBHs like GW150914, and it saturates at $z\sim 10$, while in the primordial BBH (PBBH) model, the cumulative event rate is $ \sim 3\times 10^{4}\,{\rm events~ yr^{-1}}$ at $z=30$ even if only $0.1\%$ of the dark matter consists of PBHs, and it is still increasing at $z=30$. In the Pop I/II model of BBHs, the cumulative event rate is $(3-10)\times10^{5}\,{\rm events~ yr^{-1}}$ and it saturates at $z \sim 6$. We present the requirements on orbit accuracy, drag free techniques, laser power, frequency stability, and interferometer test mass. For BBHs like GW150914 at 1Gpc, SNR$\sim 90$ is achieved with the definition of Pre-DECIGO in the $0.01-100$Hz band. Pre-DECIGO can measure the mass spectrum and the $z$-dependence of the merger rate to distinguish various models of BBHs like GW150914. Pre-DECIGO can also predict the direction of BBHs at $z=0.1$ with an accuracy of $\sim 0.3\,{\rm deg}^2$ and a merging time accuracy of $\sim 1$s at about a day before the merger so that ground-based GW detectors further developed at that time as well as electromagnetic follow-up observations can prepare for the detection of merger in advance. For intermediate mass BBHs at a large redshift $z>10$, the QNM frequency after the merger can be within the Pre-DECIGO band so that the ringing tail can also be detectable to confirm the Einstein theory of general relativity with SNR$\sim 35$. [abridged]


Introduction
The first direct detection of a gravitational wave (GW) has been done by O1 of aLIGO [1]. The event, called GW150914, was a binary black hole (BBH) of mass ∼ 30M -30M . Such high mass black hole (BH) candidates had not been confirmed, although several suggestions existed [2][3][4][5][6][7][8][9][10], so that its origin is not known at present, while for the merger of a neutron star (NS) binary there exist several systems with merger time less than the age of the universe so that the event rate has been estimated [11][12][13], and this was the most plausible GW source for the first direct detection of GWs before GW150914 [14]. For BBHs, no observation of electromagnetic counterparts exists so that analyzing theoretical models is the only way to provide methods to study their properties. One method is population synthesis, in which Monte Carlo simulations have been performed to evolve binaries starting from binary ZAMS (Zero Age Main Sequence) stars. The mass of the observed BH candidate in X-ray binaries is at most ∼ 15M , which is about half of the mass of GW150914. This suggests that the progenitors of GW150914 are low-metal stars with little or no mass loss such as Pop II or Pop III stars [15,16]. In particular, the predicted mass of Pop III BBHs by Kinugawa et al. [4] (see also Refs. [6,17,18]) agrees astonishingly well with GW150914 [19]. However, a single event is not enough to restrict the origin of BBHs like GW150914, around 6-8 of which will be found in O2. The cumulative chirp mass and the total mass as well as the distribution of the spin parameter a/M of the merged BH will help to distinguish the plausible model among the various population synthesis models. Here, in the population synthesis Monte Carlo simulations of Pop II and III stars, there are so many unknown functions and parameters such as initial mass function, initial eccentric distribution function, the distribution of initial separation, the distribution of mass ratio, the Roche lobe overflow parameters, the common envelope parameters, and so on, so that the observed cumulative chirp mass and the total mass as well as the distribution of the spin parameter will only give constraints among these undetermined functions and parameters.
There are other completely different formation scenarios of BBHs like GW150914, such as primordial BBH (PBBH) [20] 1 and the three-body dynamical formation model [7]. In particular, in the PBBH model the mass spectrum primarily reflects the primordial density perturbation spectrum. This means that it is difficult to distinguish models only by observations of small-z BBHs like GW150914 since the average detectable range of GW150914 with the design sensitivity of aLIGO is at most z ∼ 0.3, whose luminosity distance is ∼ 1.5 Gpc.
It is highly requested to observe GW150914-like events for larger z to distinguish the various models. For this purpose, DECIGO (DECihertz laser Interferometer Gravitational wave Observatory) is suitable although it was originally proposed by Seto, Kawamura, and Nakamura [22] to measure the acceleration of the universe through GWs from binary NS-NS at z ∼ 1 and observe the chirp signal 1-10 years before the final merger. The GW frequency from such NS-NS binaries is ∼ 0.1 Hz (= decihertz), which is the origin of the name of DECIGO. 2 DECIGO consists of at least a triangle shape of three spacecraft in a heliocentric orbit with a 1000 km Fabry-Pérot laser interferometer to particularly observe the GW background (GWB) from the inflation at frequency 0.1 Hz [23]. Pre-DECIGO 3 is a smaller DECIGO which consists of three spacecraft arranged in an equilateral triangle 1 They simply applied the method in Ref. [21] for the case with mass from 0.5M to 30M . 2 Another origin of the name is "DECIde and GO project." 3 Pre-DECIGO was initially the pathfinder for DECIGO. That is, Pre-DECIGO was a smaller version of DECIGO, and the design of Pre-DECIGO had not been defined including the selection of possible targets, while DECIGO has a definite design and clear targets. In particular, the sensitivity of Ultimate-DECIGO is limited only by the uncertainty principle so that it can detect the inflation origin background GWs even if GW = 10 −20 . with 100 km arm length orbiting 2000 km above the surface of the earth. The orbit of Pre-DECIGO is geocentric and is different from DECIGO, whose orbit is heliocentric. We are hoping that the launch date will be in the late 2020s. In this paper, we show that Pre-DECIGO can detect events like GW150914 up to z ∼ 30 where the cumulative event rate is ∼ 1.8 × 10 5 events yr −1 in the Pop III origin model of GW150914, while in the PBBH model, it is ∼ 3 × 10 4 events yr −1 even if only 0.1% of the dark matter consists of primordial BHs (PBHs). This paper is organized as follows: In Sect. 2, the requirements on the orbit accuracy, drag-free techniques, laser power, frequency stability, and interferometer test mass will be shown. In Sect. 3, we show that Pre-DECIGO can measure the mass spectrum and the z-dependence of the merger rate up to z ∼ 30 to distinguish various models such as Pop III BBH, Pop II BBH, PBBH, threebody dynamical formation models, and so on. For small z = 0.1, Pre-DECIGO can also predict the direction of BBHs with accuracy ∼ 0.3 deg 2 , and the merger time with accuracy ∼ 1 s at about a day before the merger so that the Einstein Telescope (ET) [24] and the enhanced version of aLIGO as well as electromagnetic follow-up observations can prepare for the detection of the merger in advance, like a solar eclipse. For large z > 10 the quasinormal mode (QNM or ringing tail) frequency after the merger can be within the Pre-DECIGO band so that the ringing tail can also be detectable to confirm or refute the Einstein theory of general relativity (GR) with signal-to-noise ratio (SNR) ∼ 35 for intermediate mass BBHs such as ∼ 640M -640M . Section 4 is devoted to discussions.

Design of Pre-DECIGO
Pre-DECIGO is a space-borne GW antenna with 100 km arm lengths. It consists of three spacecraft separated by 100 km (Fig. 1). Test-mass mirrors in each spacecraft have masses of 30 kg and diameters of 30 cm. With these mirrors, a triangle-shaped laser interferometer unit is formed by three 100 km Fabry-Pérot cavities. The finesse of the cavities is 100, resulting in a cavity cut-off frequency around 20 Hz. As a laser source, a frequency-doubled Yb:fiber DFB laser with a wavelength of 515 nm and an output power of 1 W is used. The frequency of the laser source is pre-stabilized in reference to the saturated absorption of iodine molecules, and also stabilized by the common-mode signal of the 100 km arm cavities. So as to avoid the external force fluctuation on the test mass caused by the spacecraft motion, the test masses are kept untouched inside the spacecraft. In addition, the spacecraft is drag-free controlled; the displacement and attitude of the spacecraft are controlled by using the test masses inside it as references.   [25]; sensitivity curve at http://gwcenter.icrr.u-tokyo.ac.jp/en/researcher/parameter), thirdgeneration antenna (ET [26]), and space antenna (eLISA [27]) are shown together for references. The dashed curve shows the signal amplitude from BBH merger with masses of 30M at a distance of z = 1. The target sensitivity of Pre-DECIGO is 2 × 10 −23 Hz −1/2 in strain in the current design (Fig. 2). It is fundamentally limited by the optical quantum noise of the interferometer: laser shot noise and radiation pressure noise in high and low frequency bands, respectively. The external force noise on the test-mass mirrors is also critical; the requirement is 1 × 10 −16 N/Hz 1/2 . With this sensitivity, mergers of BBHs at z = 10 will be within the observable range of Pre-DECIGO, assuming optimal direction and polarization of the source, and detection SNR of 8 (Fig. 3).
The candidate orbit of Pre-DECIGO is a record-disk (or cartwheel) orbit around the earth. The reference orbit, the orbit of the center of the mass of the three spacecraft, is a Sun-synchronized duskdawn circular orbit with an altitude of 2000 km. Each spacecraft has a slightly eccentric orbit around this reference orbit so as to minimize the natural fluctuation in the relative distance between the spacecraft during the orbital motion. Moreover, there will be no eclipse in these spacecraft, which is beneficial to avoid thermal shock and drift in the spacecraft. The formation flight of the three spacecraft is realized by continuous feedback control. The laser interferometers measure the cavity length changes, which are fed back to the positions of the test-mass mirrors. Since the spacecraft follows the test-mass position inside it by drag-free control, the 100 km triangular formation flight is realized. The orbital period of the formation flight interferometer unit around the earth is 124 min. With this orbital motion and the earth's annual orbital motion around the sun, the antenna pattern of Pre-DECIGO to observe GWs will change in time. Parameter estimation accuracy for the GW sources, such as sky localization, is improved because of this selection of the orbit.

Physical, astronomical, and cosmological targets of Pre-DECIGO
The guaranteed science target of Pre-DECIGO is the pre-merger phase of BBHs like GW150914. For this purpose, we rewrite Eqs. (3)-(5) in Ref. [22] expressing the time before the final merge t c , the number of cycle N c , and the characteristic amplitude h c as where z, M 1 , M 2 , M t , M c , ν, and d L (z) are the cosmological redshift, the mass of each BH, the total mass, the chirp mass defined by (M 1 M 2 ) 3/5 /(M 1 + M 2 ) 1/5 , the frequency of the emitted GW assuming a circular orbit, and the luminosity distance, respectively. One might think that it is difficult to determine z since, for fixed physical masses, we have h c ∝ z −1/6 for z 1, where d L (z) is almost proportional to z. However, from the phase evolution of GW we will get the redshifted chirp mass (1 + z)M c so that from h c we will get d L (z) accurately. Then, assuming the standard CDM cosmological model value of d L (z), we can determine the redshift z.
According to Dalal et al. [28], the maximum and the average SNR 4 are given by where We obtain for a BBH of mass 30M -30M at 1 Gpc (z ∼ 0.2) SNR ∼ 90 in the Pre-DECIGO band of 0.01-100 Hz. Since the characteristic strain amplitude weakly depends on the redshift z as (1 + z) −1/6 for z 1, near face-on ∼10% of BBHs have SNR ∼ 5 even at z ∼ 30.

Pop III model
Kinugawa et al. [4] showed that Pop III binaries generally become BBHs whose mass has a peak at 30M -30M . Figure 4 shows the Pop III BBH chirp mass distribution and the cumulative chirp mass distribution of the standard model defined in Kinugawa et al. [6] where the initial mass range PTEP 2016, 093E01 T. Nakamura et al. is 10M < M < 140M with a flat initial mass function (IMF). We also give a plot for the model with the Salpeter IMF (see Ref. [6] for details about the other model parameters). The main reasons for this peak of chirp mass distribution are as follows: First, Pop III stars tend to be born as massive stars whose typical mass is 10-140M , 5 i.e., more massive than Pop I and II stars because of the larger Jeans mass of the initial gas cloud due to the lack of coolant and the weaker radiative feedback due to the lack of absorption of photons by the metal. Furthermore, since the stellar wind mass loss of Pop III stars is not effective, the stellar evolution makes them more massive compared with Pop I and II stars. Second, Pop III star evolution is qualitatively different from that of Pop I and II. Pop I and II stars evolve as a red giant at the late phase. On the other hand, Pop III stars whose masses are less than 50M evolve as a blue giant. If a star in a binary becomes a red giant, the mass transfer from such a giant star to the companion star often becomes unstable and they form a common envelope phase, during which the companion star plunges into the envelope of the giant star and spirals in. In this phase, the friction between the envelope and the companion star yields a loss of orbital angular momentum to decrease the binary separation, while the envelope will be evaporated by the energy liberated through the friction. Consequently, the binary becomes either a close binary that consists of the core of the giant or a single star absorbing the companion star during the common envelope phase. Therefore, all Pop I and II binaries initially close enough to become a close compact binary that merges within the age of the universe, evolve via a common envelope phase. However, the mass transfer of a blue giant Pop III star is so stable that Pop III stars whose mass is less than 50M do not experience a violent mass loss process such as the common envelope phase. They evolve via stable mass transfer phases and their mass loss is smaller than that in the evolution via a common envelope phase. They tend to lose only a minor fraction (1/10-1/3) of their initial mass so that they tend to end up with 20-30M BHs, while Pop III stars whose mass is larger than 50M are also likely to form a common envelope phase, and they lose 1/2-2/3 of their mass so that they tend to end up with 25-30M BHs, too. Therefore, the peak of the Pop III BBH chirp mass distribution at 20-30M can be understood as a natural consequence of the evolution of Pop III stars (see Fig. 4).
Kinugawa et al. [4,6] showed the merger rate of Pop III BBHs using their Pop III binary population synthesis result and the Pop III star formation rate (SFR) of de Souza et al. [30]. The 5 A Pop III star with mass larger than 140M explodes as a pair instability supernova or fragments.   Fig. 5 shows the cumulative merger rate of Pop III BBHs for the standard model and the Salpeter IMF model. The peak of Pop III BBH merger rate density is at z ∼ 8 because the Pop III star formation peak is at z ∼ 9 [30]. Thus, the cumulative merger rate is saturated at z ∼ 10.

PopII/I model
It is usually believed that Pop I and II stars are lighter than Pop III stars because the minimum Jeans mass of Pop I and II stars is lighter than that of Pop III, which is caused by the difference in cooling mechanism; that is, metal cooling and dust cooling are at work for Pop I and II, while only hydrogen atom cooling and molecular cooling apply for Pop III [31]. Although the Pop II IMF is not known well, the Pop I IMF is known rather well by observations so that the typical mass of Pop I is ≤ 1M [32,33]. Furthermore, Pop I and II stars lose mass via stellar wind due to the absorption of photons by the metal. Therefore, Pop I and II compact binaries tend to be lighter than those of Pop III.
It is considered that the stellar wind mass loss rate is an increasing function of metallicity [34], i.e., Pop II binaries tend to be more massive compact binaries than Pop I. Thus, the maximum mass  of Pop II BBHs is greater than Pop I BBHs. 6 In the case of Pop II BBHs, however, the peak of the chirp mass distribution is almost the same as that of Pop I, assuming that the Pop II IMF is the same as the Pop I IMF. The reason is that in the late evolution stage, Pop II stars evolve not as a blue giant like a Pop III star with mass less than 50M but as a red giant like a Pop I star. Figure 6 shows the chirp mass distribution and the cumulative chirp mass distribution of Dominik's standard model of submodel B for Z = Z and Z = 0.1Z [3] (see also Fig. S6 in "Supplementary online text" of Belczynski's latest paper [10]).
The top right panel of Fig. 5 shows the cumulative merger rates of Belczynski's latest model [10] and previous model [36]. In these models they use a grid of 11 metallicities such as Z = 0.03, 0.02 (= Z ), 0.015, 0.01, 0.005, 0.002, 0.0015, 0.001, 0.0005, 0.0002, 0.0001. The cumulative merger rates of Pop I and II saturate at z ≤ 6. Note that the cumulative merger rate depends on the SFR and the metallicity evolution models.

PBBH model
Nakamura et al. [21] 7 first pointed out that the PBH can be the source of GWs since many binary BHs can be formed by the existence of the third BH. In that paper, the prime interest was the possibility of PBH as MACHO (MAssive Compact Halo Object) of mass ∼ 0.5M [37,38] as dark matter and the source of GW for TAMA300 [39,40]. Sasaki et al. [20] applied this to GW150914-like BBHs (see also Refs. [41,42]).
If there exists random density perturbation in the radiation-dominated era, in some part of the universe the amplitude of the density perturbation can be high enough to collapse to a BH. Since the sound velocity is comparable to the light velocity in the radiation-dominated era, the Jeans length is comparable to the particle horizon size so that the mass of PBH is comparable to the horizon mass at that time. Then PBH with mass ∼ 30M was formed at ∼ 0.005 s and the temperature of the universe was ∼ 14 MeV. Once PBH is formed, it behaves like dust. Note here that the fraction of radiation which collapses to PBH is only ∼ f × 10 −7.5 at the formation time, where f is the fraction 6 This result depends on the stellar wind mass loss model [35]. 7 There are so many typos in this ApJ paper. of PBH to the total mass of dark matter at present. Ricotti, Ostriker, and Mack [43] argued that if PBH exists there will be an accretion disk around PBH and the CMB will be distorted. For ∼ 30M PBH, f is at most 10 −3 to be compatible with the data of COBE FIRAS [44]. After the equal time between the local BH energy density and the radiation energy density, some pairs of PBHs with separation x smaller than f 1/3x , wherex is the mean separation of PBHs, will merge in a free fall time during the globally radiation-dominated era. Two PBHs acquire angular momentum by the tidal force of the third PBH and they become a binary BH. Following the argument of Nakamura et al. [21], Sasaki et al. [20] obtained the distribution function of the semi-major axis a and the eccentricity e given by dP(a, e) = 3 2 The merger time t due to the emission of the GW is approximately given by Refs. [45,46]: We can integrate Eq. (8) fixing the merger time t to get the probability of merger from t to t + dt as 8 while there is an upper bound of e as Now let us calculate the cumulative event rate (R(z)) as a function of z as where we fixed the starting time t 0 = 13.7 Gyr from the big bang. The bottom panel of Fig. 5 shows the cumulative event rate per year of PBBH mergers as a function of z, corresponding to the fraction f = 0.001; that is, PBH contributes only 0.1% of the dark matter. Note that the rate is still increasing logarithmically even at z = 30. 8 When the BH fraction is even smaller, the upper bound of the integral over e is determined differently. As a result, we obtain a slightly different expression for the merger probability [20].

QNM, or ringing tail
For the LISA detector, a detailed study of ringdown GWs from massive BHs was done by Berti, Cardoso, and Will [47]. Here, we model the ringdown waveform as where f R and f I denote the real and imaginary parts of the QNM frequency, respectively, and t 0 and φ 0 are the initial ringdown time and phase, respectively. For non-spinning, equal mass BBH mergers, f R and f I are given by where we used fitting formulas given in Ref. [48] (see also Refs. [49][50][51]) to obtain the remnant BH parameters, i.e., the mass and spin parameter of the Kerr BH, and the BH parameters were converted to f R and f I by using Ref. [47].
As for the initial amplitude of the ringdown GWs, we consider a simple inspiral-merger-ringdown waveform for BBHs in the frequency domain, based on Ref. [52]. In practice, thanks to the breakthrough of numerical relativity for BBHs [53][54][55], we can have precise inspiral-merger-ringdown waveforms for BBHs (see, e.g., Refs [56][57][58]). When we focus on the amplitudes for the three phases, these are described by for f < f 1 , where we set f 1 = c 3 /(6 3/2 πG(1 + z)M t ), which is twice the frequency of the innermost stable circular orbit (ISCO) of a test particle in the Schwarzschild spacetime, and f 2 = f R . The constant C is chosen to be consistent with the amplitude of the inspiral phase, and A mer and A ring are calculated by imposing the continuity of A(f ) at f = f 1 and f = f 2 , respectively.
Using the amplitude A ring , we evaluate the SNR for the ringdown phase. In Fig. 7, we show the total mass M t dependence of SNR in the Pre-DECIGO noise curve. It is noted that for non-spinning, equal mass binaries, the ringdown event with SNR = 35 is enough to confirm Einstein's GR for the strong gravity space-time near the event horizon, i.e., whether the compact object emitting the ringdown GWs is a BH predicted by GR or not (see Refs. [59,60]). As shown in Fig. 7, we have the possibility of testing general relativity to that accuracy for M t > 640M and 1280M at z = 0.2 and 10, respectively.

Accuracy of direction and time of the merger of BBH
In this section, we investigate how accurately we can determine the waveform parameters of GWs from GW150914-like BBHs using Fisher analysis. Parameter resolutions are estimated by the inverse of the Fisher matrix ij for large SNR as [61,62]  where ij is defined as The symbol (·|·) denotes the noise-weighted inner product. Assuming stationary Gaussian detector noise, the inner product is expressed by We define the accuracy of sky localization as the measurement error in the solid angle, where α and δ denote the right ascension and declination of the GW source, respectively. We consider a GW from a coalescing binary system consisting of two point non-spinning masses with M 1 and M 2 . We use the restricted post-Newtonian (PN) waveform up to the 1.5PN order [63]: where the PN parameter x ≡ πG (1 + z) (M 1 + M 2 ) f /c 3 2/3 was introduced. F + and F × denote the "detector beam-pattern" coefficients, ι and φ c are the inclination angle of the binary and the phase at the coalescence time t c , respectively, and η = M 1 M 2 /(M 1 + M 2 ) 2 . The vectors n and r appearing in the Doppler phase φ D are the unit vector pointing from the geometrical center of Pre-DECIGO to the source and the vector pointing from the SSB (Solar System barycenter) to the geometrical center of Pre-DECIGO, respectively. We investigate the angular resolution defined by Eq. (25) and the accuracy of the time of coalescence using Fisher analysis. The GW signal is characterized by the nine waveform parameters θ = {α, δ, t c , φ c , d L , ln M c , ln η, ψ, cos ι} (see Ref. [63] for details; the polarization angle ψ is in F + and F × ). We set the cutoff frequencies in Eq. (24) in the following way: where is the GW frequency at the ISCO, and f 1yr is the frequency at one year before the coalescence. Substituting Eqs. (26)-(32) into Eq. (22), we evaluate the accuracy of the parameter estimation. When we perform the Fisher analysis, we assume GW signals from 30M equal mass BBHs at distances of z = 0.1, 1, and 10. The polarization phase ψ is set to be 0.5 radians. The sky location (α, δ) and the inclination cos ι are assumed to obey uniform distributions and are averaged over by using 10 3 Monte Carlo simulations. As a result, we obtain Fig. 8.    (10), about 10% of them can be observed by Pre-DECIGO with SNR > 5 (7), thanks to the fact that the overall amplitude only depends on the redshift as z −1/6 for z 1. Figure 8 for z = 1, where the detector trajectory is assumed to be a circular orbit around the Sun. Figure 9 shows the time evolution of the SNR and the angular resolution of GWs from z = 0.1. The angle parameters are set to be α = δ = 1.0 rad, ψ = 0.5 rad, and cos ι = 0.5. This figure indicates that we can identify the sky location of GW150914-like BBHs with an accuracy of about 0.3 deg 2 at about a day before the coalescence. Note that since SNRs are not so large in the case of z = 10 and 30, the estimation accuracies of the parameters would be overestimated [65].

Discussion
One of the observed binary NS PSR2127+11C (see Ref. [66] and references therein) has a similar parameter to the Hulse-Taylor binary pulsar, although PSR2127+11C is in the globular cluster (GC) M15. This suggests the possibility of the formation of BBHs in the GC. A BH of mass ∼ 30M is much larger than the typical mass of the constituent stars, ∼ 1M , so that it will sink down to the center of the GC or star cluster due to dynamical friction. Then, BBHs can be formed in the central high density region of GCs. Since the escape velocity from GCs is 10 km s −1 or so, the kick velocity in the formation process of BHs or the kick when BBHs are formed by three-body interaction is high enough for BBHs to escape from GCs. Rodriguez, Chatterjee, and Rasio [67] performed such a simulation to show that the event rate is at most ∼ 1/7 of Pop I and II origin BBHs. If we take their result as it is, the dynamical formation of binaries in GCs gives only a minor contribution of Pop II origin of BBHs. From only the chirp mass, total mass, and spin angular momentum, it will be difficult to distinguish the origin of GW150914-like BBHs. This is because the number of parameters that can be determined by the distribution function of the GW data is much smaller than that of the unknown model 13  parameters and the distribution functions assumed in each model. However, the redshift distribution of GW events varies robustly among the models. Namely, the maximum possible redshift is ∼ 6, 10, and > 30 for Pop I/II, Pop III, and PBBH models, respectively (see Fig. 10). In Fig. 10, we show the event rates for each model. These event rates are derived by differentiating the cumulative event rate in Fig. 5 with respect to ln z. To observe the maximum redshift as a smoking gun to identify the origin of GW150914-like events, the construction of Pre-DECIGO seems to be the unique possibility.
Pre-DECIGO can observe NS-NS and NS-BH mergers. However, no detection of GWs from the merger of these systems has been done, though many simulations exist. For the same distance of the source, the SNR for NS-NS and NS-BH (30M ) are 0.08 and 0.25 times smaller than for 30M -30M BBHs. We will here postpone discussing what we can do using Pre-DECIGO about these sources until the first observations of GWs from these systems, since the event rates are still uncertain and might be very small compared with BBH mergers.