Abstract

1 Introduction

Finding genetic variation associated with bacterial phenotypes such as antibiotic resistance, virulence and host specificity has great potential for a better understanding of the evolution of these traits, and may be able to inform new clinical interventions. Genome-wide association studies (GWAS) address this in a hypothesis-free manner. The recent availability of thousands of whole-genome sequences from bacterial populations has made this approach possible, though issues of strong clonal population structure and a variable pangenome (sequence variation of both core regions present in every isolate, and in variably present accessory regions) must be accounted for in any successful analysis (Power et al., 2017).

One method that addresses these issues in a scalable manner is SEER (Lees et al., 2016), which uses non-redundant variable length k-mers between 9-100 bases as a generalized variant to represent variation across the pangenome. Linear models with a control for population structure are then used to perform the association. Other methods include: bugwas, which uses a linear mixed model (LMM) and also checks for lineage effects (Earle et al., 2016); scoary, which tests clusters of orthologous genes (COGs) with and without accounting for population structure (Brynildsrud et al., 2016); and phylogenetic regression, which uses a known phylogeny to adjust for the covariance between samples (Garland and Ives, 2000).

Each of these methods has its own set of advantages and limitations, so recent GWAS analyses have used a combination of these techniques along with methods designed for human genetics, often in a somewhat ad-hoc manner that requires tedious file format conversion and familiarity with a wide range of methods and software packages (Lees et al., 2017).

As the number of large bacterial population datasets increases, we recognize the need to make these methods more accessible to the field. We have therefore re-implemented SEER in an easy-to-use and install python package, pyseer. pyseer also includes many additional features covering new association models, input sources and output processing. This brings together all the methods mentioned above into a single piece of software, and further enables new combinations of analysis. Our method is equally applicable to viral GWAS, and indeed any species with structured populations and/or genomes with significant variation in gene content.

2 Materials and methods

The foundation of pyseer is a direct python re-implementation of SEER (originally written in C++). K-mers of variable length counted from draft assemblies are used as the input, and their association with a phenotype of interest is assessed by fitting a generalized linear model to each k-mer. To control for population structure, multi-dimensional scaling of a pairwise distance matrix is performed and these components are included as fixed effects in each regression. In the case of a binary phenotype with a highly penetrant variant (those with high effect sizes), Firth regression is performed to maintain power. This adds a penalty to the logistic regression likelihood when there is a large error on the slope, which alleviates underestimation of p-values in cases with nearly separable data (Heinze and Schemper, 2002). Significant k-mers, after adjusting for multiple testing, can then be mapped to a reference annotation to find regions of the genome associated with the phenotype.

After re-implementation, results were the same as SEER (Supplementary Fig. S1) and pyseer is of comparable speed (Supplementary Table S1). We then expanded the features in pyseer to bring together the different methods mentioned above.

3 Discussion

Starting with a re-implementation of SEER, we have added the models and input types used by other microbial and human GWAS approaches into a single package. Analyses which were previously challenging to perform, such as association of all variable length k-mers with a LMM, can be performed in an efficient and user-friendly manner. We have also enabled new types of analysis, such as population structure corrected burden testing, and gene level summaries of k-mers. Our package includes comprehensive documentation and a tutorial which shows how to perform GWAS using the new input sources and both association models, as well as how to interpret significant k-mers. We have implemented unit tests in our code to ensure consistency of output as features are added.

With pyseer we have therefore reconciled many of the existing methods for regression-based microbial GWAS into a single package. Our focus on documentation and ease-of-use of pyseer will make GWAS more accessible to the microbial genomics community.

Supplementary Material

Acknowledgements

References