Correcting for Bias in Estimates of θ w and Tajima's D From Missing Data in Next-Generation Sequencing

Population genetic analyses use information from the site frequency spectrum to infer evolutionary processes. Two summary statistics, Watterson's estimator ( ${\theta }_w$ ) of genetic diversity, and Tajima's $D$ , used for detecting non-neutral evolution, are among the most frequently computed statistics utilising this information. However, missing information in genomic data, particularly as encoded in the Variant Call Format (VCF), can bias these estimates, leading to incorrect evolutionary inferences. We assessed the impact of missing data on the estimation of these statistics using various population genetic software packages (VCFtools, PopGenome, pegas and scikit-allel). By simulating neutral genomic data with varying levels of missing genotypes and sites, we found consistent underestimation of ${\theta }_w$ across programs. We found a consequent bias in estimates of Tajima's $D$ , though the direction varied by software. We developed and implemented correction methods as functions in an update of the popular pixy software, significantly reducing these biases. Our findings highlight the need for accurate data handling in population genomics to avoid misinterpretations of evolutionary phenomena.