Gravitational Wave Statistics for Pulsar Timing Arrays: Examining Bias from Using a Finite Number of Pulsars

TL;DR

This study investigates how using a limited subset of pulsars in PTA data analysis can bias Bayesian upper limits on gravitational wave background amplitude, often leading to underestimations.

Contribution

It demonstrates that selecting specific pulsars can systematically bias upper limit estimates, highlighting potential pitfalls in PTA data interpretation.

Findings

Using fewer pulsars can produce lower upper limits than the true value.

Bias can occur in over 10% of simulated realizations when selecting pulsars to minimize upper limits.

Subset selection influences the distribution of recovered upper limits, often shifting them lower.

Abstract

Recently, many different pulsar timing array (PTA) collaborations have reported strong evidence for a common stochastic process in their data sets. The reported amplitudes are in tension with previously computed upper limits. In this paper, we investigate how using a subset of a set of pulsars biases Bayesian upper limit recovery. We generate 500 simulated PTA data sets based on the NANOGrav 11-year data set with an injected stochastic gravitational wave background (GWB). We then compute upper limits by sampling individual pulsar likelihoods, and combine them through a factorized version of the PTA likelihood to obtain upper limits on the GWB amplitude using different numbers of pulsars. We find that it is possible to recover an upper limit (95\% credible interval) \textit{below} the injected value, and that it is significantly more likely for this to occur when using a subset of pulsars to compute the upper limit. When picking pulsars to induce the maximum possible bias, we find that the 95\% Bayesian upper limit recovered is below the injected value in 10.6\% (53 of 500) of realizations. Further, we find that if we choose a subset of pulsars in order to obtain a lower upper limit than when using the full set of pulsars, the distribution of upper limits obtained from these 500 realizations is shifted to lower amplitude values.