with \(\alpha =1.1, \beta =0.75, R_0 = 20 \text{Gpc}^{-3} \text{yr}^{-1}, \kappa =4\) , between redshift 0 and 5 and for \(m_1 \in [2, 500] M_{\odot }\) , \(m_2<m_1\) . We sample the angles drawing the source position in the sky \((\theta , \phi )\)  as well as the time of the day of the event from uniform distributions and the inclination angle \(\iota \)  from a distribution flat in \(\cos {\iota }\) . Then, we compute the network SNR of each sampled event for a three detector network composed by the Virgo and the two LIGO detectors, taking into account their position on the Earth at the given time and the correct orientation of the interferometers' axes. We use PSDs indicative of the O1-O2, O3a, and O3b for each observation period, as given in Figures 2 of [1]} for O2, [1]} for O3a, and [3]} for O3b.The corresponding data are available at https://dcc.ligo.org/P1800374/public/ and https://dcc.ligo.org/LIGO-P2000251/public for O2 and O3a, while we compute the PSD with the \(\tt {pycbc}\)  software around the GPS time indicated in Fig. 2 of [3]} for O3b.
Finally, we take into account the duty cycles of each detector during the O1-O2, O3a and O3b observing run accordingly to the sensitivity status for each run.Available at https://www.gw-openscience.org/summary_pages/detector_status/O2/, https://www.gw-openscience.org/detector_status/O3a/, https://www.gw-openscience.org/detector_status/O3b/ We use the \(\tt {ImrPhenomPv2}\)  waveform for O1-O2 and O3a, and the \(\tt {ImrPhenomXPHM}\)  waveform for O3b to match the ones used in the data releases when possible (this is the case for O1-O2 and O3b, while for O3a a mixture of waveforms was used; we checked that the effect of using different waveforms is small).
This procedure still contains some approximations. One is that we use only the average duty cycle of the detectors during the observing runs. This results in an average network duty factor that gives a total time with no interferometer active which is lower than the one reported in the sensitivity status pages (e.g. \(\sim 1.5 \%\)  in our simulation vs \(\sim 3 \%\)  in the sensitivity status report for O3b). This is due to the fact that, not knowing the exact detector status at each time, the best we can do is to “shut down” each detector randomly at a fraction of the detection time corresponding to that given in the sensitivity status reports.
The second approximation is the use of a fixed PSD. Incorporating the actual noise of the detectors throughout all the observing runs requires significantly larger computational efforts. However, using a relatively high SNR threshold such as 12, used in our main analysis, prevents us from being too sensitive to fluctuations in the detectors' PSD. A related point is that our threshold is imposed only on the network SNR and does not include conditions on the False Alarm Rate. However, we note that, once the SNR threshold of 12 is used, all the events considered have a very small \(\text{FAR}<2 \times 10^{-4} \text{yr}^{-1}\) . We also checked that, even when using a lower SNR threshold of 8, the population analysis of GWTC-2 with a broken power-law mass function perfectly agrees with the result of [5]} despite all the differences discussed so far.We compared with the posterior samples reported in the data release associated to the LVC population analysis paper [5]}, available at https://dcc.ligo.org/LIGO-P2000434/public. In particular, all the mass scales, relevant for this analysis, match perfectly, while the only parameter with a slight difference is the overall rate \(R_0\) , for which we get \(R_0=12.1^{+5.5}_{-3.9}\)  at \(68\%\)  C.L., while in [5]} \(R_0=19.1^{+16.2}_{-9.0}\)  at \(90\%\)  C.L. is found. This slightly lower rate is consistent with the fact that we find a smaller fraction of the observational time where no interferometer is active. The two results remain in any case statistically indistinguishable. Finally, we find very good agreement also with the result for \(H_0\)  obtained in [8]}, see App. .
