dTOURS: Dense-region tagging for outbreak detection using ratio statistics

Abstract

Surveillance for food safety in the United States of America is a collaborative effort among public health agencies with additional partners worldwide contributing sequence data. Assemblies in GenBank and sequence reads in the Sequence Read Archive for surveilled species are received, rapidly analyzed, and results published publicly by an automated pathogen detection pipeline at the National Center for Biotechnology Information. The pipeline detects close isolates with a recent common ancestor by finding single nucleotide polymorphisms (SNPs) in genomes for pairs of isolates. Very few vertically transmitted SNPs are expected between a pair of close isolates; any genomic region with many SNPs compared to the number of SNPs in the rest of the genome is indicative of a horizontally transferred region that needs to be excluded for counting vertically transmitted SNPs. We developed dTOURS that adapted the ratio statistic for finding outliers to the problem of finding regions of high SNP density in a pair of genomes where isolates typically have fragmented genome assemblies. Simulations for deciding the dTOURS parameter are presented. We illustrate correctness of dTOURS using five published outbreaks, one each for five bacterial species that cause many foodborne outbreaks or lead to a high mortality rate. Comparison to Gubbins shows that while both Gubbins and dTOURS use the ratio statistic, the implementation in dTOURS is more robust for finding close isolates in outbreak analysis. Comparison with the method used by the Food and Drug Administration shows that their method is simple and fast but not sensitive.

Fig 1. Flowchart for dTOURS algorithm.

Input is a set of SNPs on a pair of genomes and output is HDRs.

Fig 2. SNPs between GCA_011671145.1 and GCA_010507495.1 before filtering.

SNPs are color coded by contig. Coordinates on the X-axis are by concatenating aligned blocks of contigs in the order in which they appear in the subject genome GCA_010507495.1.

Fig 3. Distribution of RS score for some simulated combinations.

With N={10, 25, 50, 75, 100} and G={1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 12.7} million, combinations shown are (A) 10 positions from all counts in G (B) 100 positions from all counts in G (C) all positions in N from 1.5 million (D) all positions in N from 12.7 million.

Fig 4. Clades in Listeria outbreak.

Number of members in each clade and subclade are shown in round brackets. Number of pairs in each subclade where the unfiltered SNP count is at least 100 more than the filtered SNP count in our analysis is shown in square brackets.

Fig 5. Plasmid alignment and SNPs on contig 67 in SRR3476797 assembly with respect to SRR3476805 assembly.

Length of contig 67 in SRR3476797 assembly by SKESA and plasmid LR890577.1 is 41,027 bp and 92,280 bp, respectively. Each vertical pink bar shows a single SNP.