Probing Primordial Black Hole Mergers in Clusters with Pulsar Timing Data

The gist: The analysis systematically addresses whether the stochastic gravitational wave background suggested by Pulsar Timing Array (PTA) datasets can be substantial enough to compete with that from scalar induced GWs, finding that scalar-induced GWs dominate in a broad PBH mass distribution, although specific PBH models with enhanced clustering could yield a background dominated by PBH mergers.

PBH Mass Function Derivation

The paper derives the PBH mass function from a power law primordial spectrum of curvature perturbations using the Press-Schechter formalism. The mass is approximated by the horizon mass at re-entry, given by MH(k) = 4πM2p/H(1) ≃ 20 M⊙ k106 Mpc−1 g4∗,s(Tcr)g−3∗(Tcr) (2). The PBH mass function fPBH(m) is then expressed in terms of the primordial power spectrum Pζ (k) and the critical threshold δc. For a broad mass distribution obtained for values of ns ∼ 1, the behavior of the equation of state w(T) with temperature, especially during the QCD epoch, can play an important role.

Scalar Induced Gravitational Waves (SIGWs)

A stochastic spectrum of second-order GWs is induced by scalar perturbations. For Gaussian perturbations, this induced GW spectrum can be sizeable for enhanced primordial curvature perturbations at small scales, typically around Pζ ≃ 10−2. The observable induced GW spectrum today is defined as omegagw,0 ≡ (1/ρc) dρ0gw/d ln f, which can be expressed in terms of the induced GW spectrum in the radiation era through [24, 86]. For the spectrum defined in Eq. (3), one can derive an approximate analytical expression for the induced GW spectrum today given by omegaSIGW gw,0(f) h2 = 7.5 × 10−5 g∗(T) (g∗,s(T))−4/3 Q(ns) (f f⋆2(ns−1)) Θ(f − fmin)Θ(fmax − f).

GW Background from PBH Binaries

The GW background from PBH mergers is considered through two main channels: i) early Universe PBH binaries that decoupled from the Hubble flow before matter-radiation equality, and ii) late-Universe PBH binaries that can dynamically form in clusters via capture. The merger rate from both these channels has many uncertainties as discussed in section IV. For the SGWB from unresolved PBH binaries, the spectrum follows the standard relationship to the comoving number density of compact binary sources. The differential GW spectrum emitted in Eq. (25) can be expressed by plugging Eq. (24) and Eq. (14) in Eq. (25), giving omegaPBHB gw,0(f) = f 2/3(π G)2/3ρc Z dz d ln m1d ln m2 M5/3c with V = Vearly+Vlate.

Results and Interpretation of PTA Data

The total GWB is given by summing Eq. (13) and Eq. (29), i.e., h2omegagw(f) = h2omegaSIGW gw (f) + h2omegaPBHB gw (f). In the fixed clustering analysis, the PTA signal is always dominated by the scalar induced GWB, omegaSIGW gw (f) with the PBH binaries giving a much sub-leading contribution to the total GW spectrum. To explain the PTA common-spectrum process, a spectral index ns ∼ 1.3 would yield the required increasing frequency dependence of the GW spectrum from SIGWs, omegaSIGW gw ∝ f2(ns−1).

Constraints on Primordial Power Spectrum

The analysis imposes the CMB µ-distortion constraints on the curvature power spectrum at scales k ≲ 105 Mpc−1. Since gravitational waves act as radiation, they are thus constrained by current bounds on the effective number of neutrino species limit, ∆Neff < 0.28. The LVK constraint imposes that, at fLVK = 25 Hz, omegagw ≤ 1.7 × 10−8. The posterior values for the fixed clustering analysis in Fig. 8 do not violate this constraint.

Conclusion and Discussion

The main results from the Bayesian inference show that the contribution to the GWB from PBH binaries (both early and late channels) is largely subdominant with respect to the GWB from scalar induced perturbations for power law PBH mass distributions (see Sec. V A). It is not possible to explain the PTA signal with only PBH binaries in late time PBH clusters, under the assumption of standard PBH formation scenario with PBH clustering from Poissonian fluctuations. The analysis reveals that the values of Rcl required to explain the PTA signal with late PBH binaries is ∼ O(104) Gpc−3 yr−1 larger than the Rcl values arising from Poissonian fluctuations.

How it works

The analysis performs Bayesian searches using Monte Carlo Markov Chains (MCMC) in PTArcade [112, 113] to implement the GW signal h2omegagw(f) in Eq. (30). The search is divided into two parts: Section V A computes the clustering factor Rcl appearing in the late-time PBH merger rate in Eq. (24) as a function of redshift z from primordial power spectrum parameters, such that at any z, Rcl ≡ Rcl(Aζ, ns, kmin, kmax) is fixed. Section V B considers the clustering factor Rcl appearing in Eq. (24) as an extra free parameter in addition to the primordial spectrum parameters.

Summary of Key Findings

**- The analysis shows that scalar-induced GWs dominate the GWB from PBH mergers for a broad PBH mass distribution derived using Gaussian curvature perturbations. The IPTA DR2 strongly favors an astrophysical origin for the reported common-spectrum process over the PBH models considered in this analysis. **

**- The PTA signal is always dominated by the scalar induced GWB, omegaSIGW gw (f) with the PBH binaries giving a much sub-leading contribution to the total GW spectrum. **

**- The required clustering factor Rcl to explain the PTA signal with late PBH binaries is ∼ O(104) Gpc−3 yr−1 larger than the Rcl values arising from Poissonian fluctuations. **

**- The analysis reveals that PBH mass distributions with ns ≲ 1 are statistically preferred, because they naturally suppress low-mass PBHs (from overclosing the Universe) and do not require a small-scale cut-off kmax. **

**- The resulting maximum likelihood GW spectrum for this analysis is shown Fig. 5 (right) where the total GW background is a combination from scalar induced (shown in blue solid) till the cut-off kmax(fmax) ∼ 107 Mpc−1 (10−8 Hz). **

**- The best fit values used for the GW spectrum shown in Fig. 5 are in violation of this LVK bound. **

How it works

The analysis uses a broad PBH mass function including features from the QCD epoch (as detailed in sec. II). The broad PBH mass function used in our analysis has been derived assuming the direct collapse of enhanced density perturbations obtained from a power-law power spectrum of curvature perturbations at small scales, with a large-mass cut-off (or small k, kmin) introduced to evade the constraints from CMB µ-distortions and a low-mass cut-off (or large k, kmax) to prevent overclosing the Universe with light PBHs.

PBH Constraints

In Fig. 6, we show the posterior predictive distribution for the PBH abundance fPBH(m) (as defined in) compared to existing (most stringent) constraints on PBHs. The two other analyses including SIGWs, however, evade all these limits and are viable. The posterior values for the fixed and free Rcl analysis, including the contribution from both SIGWs and PBH binaries to the GWB peak at fPBH ∼ 10−5 − 10−4, at PBH mass of around 4−10 M⊙.

Conclusion

In this work, we have searched for a GW background arising from inspiralling PBH binaries in the International PTA dataset (IPTA DR2).

1. Bold header: Enhanced Bayesian Inference for GW Source Discrimination

This system can perform a Bayesian search in the International PTA Data Release 2 (IPTA DR2) for a combined GW background arising from scalar perturbations and unresolved PBH mergers, allowing it to systematically assess whether the observed signal is sourced by scalar-induced GWs or PBH mergers.

2. Bold header: Parameter Space Mapping with Clustering Nuisance

The improved system can utilize the free clustering factor Rcl as an extra free parameter in addition to the primordial spectrum parameters to determine if a GWB contribution from PBH binaries is substantial enough, allowing it to investigate models where a much more dense PBH cluster or enhanced PBH clustering could possibly lead to much enhanced merger rates.

3. Bold header: Constraint Evaluation Against Multi-Wavelength Probes

The AI can rigorously compare its derived posterior distributions against a wide array of constraints, including CMB µ-distortions, the LVK-Constraint, and microlensing surveys (such as EROS, ICARUS, OGLE), to provide a comprehensive assessment of viability.

4. Bold header: Mass Function Sensitivity Analysis

The system can analyze how different primordial power spectrum parameters affect the resulting PBH mass distribution, noting that the final PBH abundance fPBH is exponentially sensitive to the primordial power amplitude Aζ (through the variance σ2k defined in Eq. (6)).

5. Bold header: Model Comparison via Bayes Factors

The AI can calculate and compare the relative evidence between competing physical models, such as SIGWs + PBHBs [Fixed Rcl] log10(BSMBHB,PBH) = 1.88 versus other scenarios to determine which model is preferred by the PTA data.

We consider the possibility that the stochastic gravitational wave (GW) background suggested by Pulsar Timing Array (PTA) datasets is sourced by Primordial Black Holes (PBHs). Specifically, we perform a Bayesian search in the International PTA Data Release 2 (IPTA DR2) for a combined GW background arising from scalar perturbations and unresolved PBH mergers, assuming a broad PBH mass distribution. In our analysis, we incorporate constraints on the curvature power spectrum from CMB μ-distortions and the overproduction of PBHs, which significantly suppress the contribution of PBH mergers to the total GW background. We find that scalar-induced GWs dominate the nHz frequency range, while PBH mergers alone cannot account for the observed signal under the standard PBH formation scenario involving Gaussian perturbations, and including only Poissonian PBH clustering. However, specific PBH models, such as those with enhanced clustering, could yield a GW background dominated by PBH mergers. Overall, we find that the IPTA DR2 strongly favors an astrophysical origin for the reported common-spectrum process over the PBH models considered in this analysis.