Selection in coastal Synechococcus (cyanobacteria) populations evaluated from environmental metagenomes

Abstract

Environmental metagenomics provides snippets of genomic sequences from all organisms in an environmental sample and are an unprecedented resource of information for investigating microbial population genetics. Current analytical methods, however, are poorly equipped to handle metagenomic data, particularly of short, unlinked sequences. A custom analytical pipeline was developed to calculate dN/dS ratios, a common metric to evaluate the role of selection in the evolution of a gene, from environmental metagenomes sequenced using 454 technology of flow-sorted populations of marine Synechococcus, the dominant cyanobacteria in coastal environments. The large majority of genes (98%) have evolved under purifying selection (dN/dS<1). The metagenome sequence coverage of the reference genomes was not uniform and genes that were highly represented in the environment (i.e. high read coverage) tended to be more evolutionarily conserved. Of the genes that may have evolved under positive selection (dN/dS>1), 77 out of 83 (93%) were hypothetical. Notable among annotated genes, ribosomal protein L35 appears to be under positive selection in one Synechococcus population. Other annotated genes, in particular a possible porin, a large-conductance mechanosensitive channel, an ATP binding component of an ABC transporter, and a homologue of a pilus retraction protein had regions of the gene with elevated dN/dS. With the increasing use of next-generation sequencing in metagenomic investigations of microbial diversity and ecology, analytical methods need to accommodate the peculiarities of these data streams. By developing a means to analyze population diversity data from these environmental metagenomes, we have provided the first insight into the role of selection in the evolution of Synechococcus, a globally significant primary producer.

Figure 1. Schematic of environmental metagenome tiling and dN/dS calculation.

Sequence reads were tiled to a reference genome. A majority-rule consensus sequence was generated from the aligned reads. The number of non-synonymous (NS) and synonymous (S) sites were calculated from the consensus sequence. The observed NS and S polymorphisms (polys) in the aligned reads relative to the consensus sequence are shown in bold. An example of aligned reads and dN/dS calculation is shown for a small section of a gene, but the same methodology would be used for an entire gene.

Figure 2. Number of metagenome sequence reads tiled to the CC9311 and CC9902 genomes.

Total (A) and additional (B) metagenome sequence reads tiled to the CC9311 and CC9902 genomes as the % identity (id) and % coverage (cov) was lowered.

Figure 3. Histograms of the average depth of read coverage for core and accessory genes.

The average depth of read coverage was calculated from tiling metagenome sequences to the CC9311 (A) and CC9902 (B) reference genomes. The tilings used criteria of 80% identity and 80% coverage.

Figure 4. Histograms of dN/dS ratios.

dN/dS ratios were calculated from gene alignments using metagenomic sequences from a Synechococcus clade I population (A), metagenomic sequences from a Synechococcus clade IV population (B), simulated CC9311 sequences (C) and simulated CC9902 sequences (D). The dashed line indicates where dN/dS = 1 and selection is neutral.

Figure 5. dN - dS versus the average depth of read coverage for each gene.

For genes from the core genome of a clade I population (A), the accessory genome of a clade I population (B), randomized alignments of reads tiled to CC9311 (C), simulated CC9311 sequences (D), the core genome of a clade IV population (E), the accessory genome of a clade IV population (F), randomized alignments of reads tiled to CC9902 (G), and simulated CC9902 sequences (H). The gray line is the linear regression and m is the slope of the line. Note that dN - dS is plotted instead of dN/dS.

Figure 6. dN/dS ratios from sliding windows of 201 bp from genes with transport, secretory, or motility functions.

These genes have regions of elevated dN/dS that were well-represented by metagenome reads.

Figure 7. Scatter plot of dN/dS ratios of homologous genes from a clade IV population versus a clade I population.

The dN/dS ratios were square root transformed prior to plotting to obtain a more normal distribution of values. The solid gray line is the linear regression (r² = 0.38, y = 0.58x+0.14). The dotted line has a slope of 1. The boxed outlier corresponds to an ABC transporter for sugars, solute-binding protein (sync_1402 and syncc9902_1078). Sync_1402 was very poorly covered by metagenome reads (Table 1). Not plotted are 96 genes with no observed polymorphisms and 8 genes with no observed synonymous polymorphisms based on metagenome tilings to either the CC9311 or CC9902 genomes.