Complex exchanges among plasmids and clonal expansion of lineages shape the population structure and virulence of Borrelia burgdorferi

Abstract

Background: In the United States, Borrelia burgdorferi (Bb) is the principal etiologic agent of Lyme disease. The complex structure of Bb genomes has posed challenges for genomic studies because homology among the bacterium's many plasmids, which account for ~40% of the genome by length, has made them difficult to sequence and assemble.

Results: We used long-read sequencing to generate near-complete assemblies of 62 isolates of human-derived Bb and collected public genomes with plasmid sequences. We characterized genetic diversity and population structure in the resulting set of 82 plasmid-complete Borrelia burgdorferi sensu stricto genomes. The Bb core genome is encoded by a chromosome and the conserved plasmids cp26, lp54, and lp17; the accessory genome is encoded by all other plasmids and the distal arm of the chromosome. Near-complete genomes reveal that the most granular Bb genotypes are clonal expansions of complex rearrangements among accessory genome elements. Ribosomal spacer types (RST) represent multiple collections of such genotypes, whereas OspC types are usually clonal. Structural rearrangements are non-randomly distributed throughout the genome, with cp32 plasmids undergoing dense exchanges and most linear plasmids, except lp54, sharing blocks among themselves and with the distal arm of the chromosome. OspC type A strains, known to possess greater virulence in humans, are distinguished by the presence of lp28-1 and lp56. Rearrangements among plasmids tended to preserve gene content, suggesting functional constraints among gene networks. Using k-partite graph decompositions, we identified gene sets with correlation patterns suggestive of conserved functional modules.

Conclusions: Long-read assemblies reveal that Bb population genetic structure results from clonal expansion of lineages that have undergone complex rearrangements among plasmid-encoded accessory genome elements. Genetic structure is preserved among genes even when plasmid rearrangements occur, suggesting that selection among epistatic loci maintains functional genetic networks. The analysis of near-complete genomes assembled using long-read sequencing methods advances our understanding of Bb biology and Lyme disease pathogenesis by providing the first detailed view of population variation in previously inaccessible areas of the Bb genome.

Keywords: Borrelia burgdorferi; Lyme disease; Pathogenicity; Virulence; Whole Genome Sequencing; genomics; lipoproteins; recombination.

Figure 1:. The phylogenetic structure of genetic markers and plasmid presence.

A phylogenetic tree assembled from core genome sequences is shown. Tips are colored by RST and tip labels are colored by OspC type. To the right of the tree, the results of genetic marker and cluster annotation from several typing schemes (RST, OspC, MLST, BAPS, PopPUNK core clusters, PopPUNK accessory clusters, DBScan). Plasmids are annotated by color and ordered by prevalence.

Figure 2:. Phylogenetic structure and plasmid occupancy patterns of Bb orthologs.

A phylogenetic tree assembled from core genome sequences is shown at left. Tips are colored by RST and tip labels are colored by OspC type. To the right of the phylogenetic tree, orthologs present in each isolate are shown in a matrix. For visual clarity, only orthologs present in 26 or more isolates are shown, with this threshold being selected due to the partitioning of isolates into 28 RST1, 28 RST2, and 26 RST3. Columns of the matrix are grouped by best hit replicon, with replicons ordered by presence rates across our sample. The plasmid encoding the ortholog is annotated by color.

Figure 3:. Matrix of pairwise average nucleotide identity, clustered by phylogeny.

(A) A phylogenetic tree constructed from core genome sequences is shown at left. Genetic markers and cluster membership are shown next to the phylogenetic tree. At the right, the lower triangle portion of a matrix whose elements correspond to the pairwise calculation of average nucleotide identity (ANI) is shown. Each element in the matrix is colored by the ANI value. (B) The boxplots depict genomic similarity (full average nucleotide identity) within each OspC type.

Figure 4:. Genetic structure and patterns of divergence in the core and accessory Bb pangenome.

(A) Principal component analysis shows the projection of Sorensen distances for the collection of isolates. (B) Core vs. accessory divergences. (C) Core (top panel) and accessory (lower panel) divergence grouped by RST. (D) COG2020 categories and their occupancy according to RST. (E and F) Venn diagrams of pangenome orthologs by RST (left panel) and pangenome orthologs with 100% prevalence in at least one RST (right panel). *** p < 0.001

Figure 5:. Gene group networks provide a structural view of variation by sequence type.

(A) Ternary plot of gene groups by RST types. Points are sized by gene frequency across the sample, and coloring is determined with a 3-dimensional RGB color cube mapping based on gene distribution amongst RST types. Surface lipoproteins are circled in black. For visibility, a small jitter has been applied to the points. (B) Positive correlation networks for (i) all gene groups and (ii) gene groups associated with lipoproteins. Nodes are sized by gene frequency across the sample, and coloring is determined with a 3-dimensional RGB colorcube mapping based on gene distribution amongst RST types. Edge weights indicate correlation strength. Gene groups encoding for surface lipoproteins are circled in black.

Figure 6:. K-partite analysis reveals 82 distinct modules within the Bb genome.

K-partite sets were determined using the negative correlation network, then induced subgraphs were created for each group of nodes using positive correlations.