Binary black hole merger rate constraints using GWTC-3 and full-O3 stochastic background constraints

README This dataset contains posterior measurements of the redshift-dependent merger rate, mass distribution, and spin distribution of binary black holes following the O3b observing run of the LIGO-Virgo-KAGRA network, including both direct compact binary detections and constraints on the astrophysical gravitational-wave background. In particular, the goal of this work is to constrain a more complex model for the black hole merger rate, with the comoving rate density evolving as $$R(z) \propto \frac{(1+z)^\alpha}{1 + \left(\frac{1+z}{1+z_p}\right)^{\alpha + \beta}}.$$ At redshifts \(z < z_p\), the merger rate grows approximately as \(R(z) \propto (1+z)^\alpha\), whereas at \(z>z_p\) it falls as \(R(z) \propto (1+z)^{-\beta}\). The analysis was performed as described in Callister et al (2020) and Abbott et al (2021), now using binary black hole detections from the GWTC-3 catalog (Abbott et al 2023a, 2023b). The parameter estimation samples used are those provided by the LIGO-Virgo-KAGRA collaboration at https://zenodo.org/records/8177023 Selection effects are calculated and mitigated using the suite of pipeline injections available at https://zenodo.org/records/7890398 Cross-correlation measurements of the stochastic gravitational-wave background are available at https://dcc.ligo.org/LIGO-G2001287, and correspond to results presented in Abbott et al (2021). As discussed in Abbott et al (2023b), the results of this combined BBH + stochastic analysis are categorically unchanged relative to results previously obtained using GWTC-2 events (Abbott et al 2021); sensitivities are not yet sufficient to resolve the redshift at which the black hole merger rate peaks and turns over. Contents processed_emcee_samples_together_r00r01.npy: File containing posterior samples run_emcee_plPeak.py: Script performing hierarchical inference; used to generate posterior samples Accessing posterior samples Posterior samples are contained in the file processed_emcee_samples_together_r00r01.npy. This is loaded via python as, e.g. >>> import numpy as np >>> dataset = np.load('processed_emcee_samples_together_r00r01.npy') Contained in this file is a single numpy array of size (# of posterior samples, # of hyperparameters): >>> dataset.shape (1152, 13) The 13 hyperparameters are defined as follows: Column Name Definition dataset[:, 0] xeff_mu Mean effective inspiral spin dataset[:, 1] xeff_sig Standard deviation of effective inspiral spin dataset[:, 2] R0 Total BBH merger rate at redshift \(z=0\) dataset[:, 3] mMin Minimum black hole mass dataset[:, 4] mMax Maximum black hole mass dataset[:, 5] lmbda Power-law index on primary mass distribution dataset[:, 6] mu_peak Mean of Gaussian peak in primary mass distribution dataset[:, 7] sig_peak Standard deviation of Gaussian peak dataset[:, 8] frac_peak Fraction of events occupying Gaussian peak dataset[:, 9] bq Power-law index on mass ratio distribution \(p(q\|m_1)\) dataset[:, 10] alpha Slope of \(R(z) \propto (1+z)^\alpha \) at low redshifts dataset[:, 11] beta Slope of \(R(z) \propto (1+z)^{-\beta}\) at high redshifts dataset[:, 12] zpeak Redshift at which \(R(z)\) peaks and turns over The exact usage of the above parameters can be seen in the included script run_emcee_plPeak.py, with which the inference was performed.