Genome-wide identification of bacterial plant colonization genes

<div><p>Diverse soil-resident bacteria can contribute to plant growth and health, but the molecular mechanisms enabling them to effectively colonize their plant hosts remain poorly understood. We used randomly barcoded transposon mutagenesis sequencing (RB-TnSeq) in <i>Pseudomonas simiae</i>, a model root-colonizing bacterium, to establish a genome-wide map of bacterial genes required for colonization of the <i>Arabidopsis thaliana</i> root system. We identified 115 genes (2% of all <i>P</i>. <i>simiae</i> genes) with functions that are required for maximal competitive colonization of the root system. Among the genes we identified were some with obvious colonization-related roles in motility and carbon metabolism, as well as 44 other genes that had no or vague functional predictions. Independent validation assays of individual genes confirmed colonization functions for 20 of 22 (91%) cases tested. To further characterize genes identified by our screen, we compared the functional contributions of <i>P</i>. <i>simiae</i> genes to growth in 90 distinct in vitro conditions by RB-TnSeq, highlighting specific metabolic functions associated with root colonization genes. Our analysis of bacterial genes by sequence-driven saturation mutagenesis revealed a genome-wide map of the genetic determinants of plant root colonization and offers a starting point for targeted improvement of the colonization capabilities of plant-beneficial microbes.</p></div

Putative operons containing 3 or more colonization genes.

<p>All genes selected for validation in these operons were successfully validated.</p

Selected candidates are validated using secondary luciferase-based screen.

<p>Twenty-two mutant strains were retrieved from an arrayed clone library for further validation. These strains were competed against an engineered wild-type (WT) <i>P</i>. <i>simiae</i> strain that produces luciferase. (A) False-color image of competition between the nonengineered WT WCS417r strain versus the luciferase-producing (Lux+) strain. Each group of 5 roots represents a different ratio of Lux+ to WT (0.0–1.0). Luminescence intensity is false colored in green. White rectangle indicates approximate region of the root tip used to measure luciferase activity. (B) All 22 insertion mutant strains (described in <a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.s012" target="_blank">S1 Data</a>) were competed with the Lux+ on roots and empty phytagel/mesh plates after inoculation as a 1:1 mixture. A ratio of each mutant strain to the Lux+ strain on the roots was estimated by interpolating the luminescence intensity of the root tip onto a standard curve (<a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.s005" target="_blank">S5 Fig</a>, <a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.s012" target="_blank">S1 Data</a>). The estimated fitness score (Lux+ colonization index [CI] in <a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.s012" target="_blank">S1 Data</a>) was then derived from the log-transformed root ratio minus the log-transformed mesh (no root final [NRF]) ratio (see <a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.s012" target="_blank">S1 Data</a>). Error bars are ± standard error (<i>n</i> = 3 biological replicates for y-axis, <i>n</i> = 15 for x-axis). Abbreviation: RB-TnSeq, randomly barcoded transposon mutagenesis sequencing.</p

Genome-wide map of root colonization.

<p>(A) Inner to outer tracks: transposon insertion density (per 1 kb); fitness score for genes with enhanced colonization ability when mutated; dominant cluster of orthologous group (COG) category for operons with 3 or more colonization-enriched genes (<a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.t001" target="_blank">Table 1</a>), gene density (for each strand); dominant COG category for operons with 3 or more colonization-reduced genes (<a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.t001" target="_blank">Table 1</a>); fitness score for genes with reduced colonization ability when mutated; chromosomal position. (B) Color legend of dominant COG categories and highlights are shown. (B) Distribution of genes significantly depleted or enriched among COG categories (see <a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.s012" target="_blank">S1 Data</a>).</p

Selected groups of genes with reduced or enhanced fitness scores with putative functions highlighted by in vitro fitness data.

<p>Colonization-depleted (A and B) and colonization-enriched (C) genes were selected for their functional characteristics determined by in vitro growth assays. (A) Genes and in vitro conditions with at least 1 strong phenotype (|fitness score| > 2), excluding any gene with significantly reduced fitness in motility assays. (B) Operon (ID = 1,092) has 5 colonization-depleted genes, and is also required for resistance to antibiotics, including polymyxin B. (C) Colonization-enriched genes that also have significantly reduced fitness in many in vitro assays (see <a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.s008" target="_blank">S8 Fig</a>) with conditions in which the amino acids noted are the only carbon source. These genes have profiles consistent with amino acid auxotrophy. For (A), (B), and (C), gene names are presented on the right, data can be found in <a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#pbio.2002860.s012" target="_blank">S1 Data</a>. Conditions shown are labeled on the bottom. The color scale corresponding to the fitness score is shown at the top right. “Root colonization” (in green) refers to the root fitness score as described earlier (<a href="http://www.plosbiology.org/article/info:doi/10.1371/journal.pbio.2002860#sec002" target="_blank">Results</a>). Abbreviation: COG, cluster of orthologous group; RB-TnSeq, randomly barcoded transposon mutagenesis sequencing.</p

Overview of root colonization randomly barcoded transposon mutagenesis sequencing (RB-TnSeq) screen.

<p>Wild-type <i>Pseudomonas simiae</i> WCS417r (A) was mutagenized with a mariner transposon system, generating 110,142 insertion mutants (B). Most insertion mutations do not significantly alter the growth phenotype on the root (black), while some insertion mutations make these mutant strains more (blue) or less (red) likely to colonize plant roots, while not significantly affecting their ability to grow in liquid culture or on the nylon substrate in the absence of plants (C). This mutant strain library was exposed to vertically-oriented phytagel plates with (D, left) or without (D, right) <i>Arabidopsis thaliana</i> seedlings. After colonization, surviving mutant strains from the root and from the plant-free mesh were collected and the abundance of insertion mutant strains within each population was quantified by RB-TnSeq. (F) Genes with under-represented insertion counts in the root population compared with the control population (shown here in red) were given low fitness scores (efficient colonization), while genes with over-represented counts were given high fitness scores.</p

NQS filtering improves fit of probability model to data.

<p>(<b>A</b>) Quantile-quantile (q-q) plots under NQS filtering show good fit of the probability model to the observed distribution of errors. Since the probability model is discrete, p values are projected onto a uniform distribution, and the distribution of projected p values is compared with the expected null distribution. See <a href="http://www.ploscompbiol.org/article/info:doi/10.1371/journal.pcbi.1002417#s4" target="_blank">Materials and Methods</a> section for details. (<b>B</b>) In contrast, q-q plots under no filtering show that no filtering skews the calibration of the probability model used by <i>V-Phaser</i>. Q-q plots of models based on subsets of the reads demonstrate that this effect becomes more pronounced with increasing coverage (see <a href="http://www.ploscompbiol.org/article/info:doi/10.1371/journal.pcbi.1002417#pcbi.1002417.s001" target="_blank">Figure S1</a>). Q-q plots are scaled to fit curve, so y = x line is not at a 45 degree angle.</p

Phase information increased sensitivity, and base quality scores increased specificity.

<p>We compared <i>V-Phaser</i> to alternate versions of <i>V-Phaser</i> with specific components disabled. In the No Phase version, <i>V-Phaser</i> called variants without phase information. In the Uniform Errors version, <i>V-Phaser</i> estimated uniform error rates within homopolymer and nonhomopolymer regions without regard to assigned base qualities. In the No Filtering version, <i>V-Phaser</i> did not filter out low quality bases. (<b>A</b>) Phase information increased sensitivity. The version without phase information attained a sensitivity of 90%, but all other versions of <i>V-Phaser</i> used phase information and attained a sensitivity of 97% or more. We calculated sensitivity as the percentage of known variants correctly identified. Data are from WNV mixed population control dataset. (<b>B</b>) Individual base quality scores increased specificity. Among loci with mismatches, the Uniform Errors version had only 91% specificity, but all other versions incorporated base quality scores in their probability model and attained 97% specificity or more. We calculated specificity as the percentage of loci in the control sample correctly identified as having no variants among loci that had at least one candidate variant. Data are from infectious clone (HIV NL4-3) control dataset.</p

Error rates were not uniformly distributed.

<p>Error rates varied by (<b>A</b>) read position, (<b>B</b>) base transition, and (<b>C</b>) base quality score. We counted as errors any mismatches to the consensus assembly for each of the two runs in the control read set under the assumption that the NL-43 infectious clone had no diversity. We defined the read position relative to the beginning or end of the read, whichever was closer. We defined a base transition as a dinucleotide representing the transition from the preceding base to the current base, and we scored a transition as an error if the current base was a mismatch. Base quality scores came from the sequencing process.</p

Phase increased sensitivity to detect variants.

<p>Phase increased sensitivity to detect variants, as seen over a range of error rates at coverages of 100-fold, 250-fold, and 500-fold. The <i>phased variant detection threshold frequency (VDTF)</i> is the lowest frequency of reads with variants at two specific loci that <i>V-Phaser</i> can distinguish from error among reads that span both loci. The <i>unphased VDTF</i> is the lowest frequency of one variant that <i>V-Phaser</i> can distinguish from error among reads that cover that locus. 100-fold <i>phased</i> sequence coverage achieves comparable detection thresholds as 500-fold <i>unphased</i>. We use Equation 7 to calculate the <i>phased</i> and <i>unphased VDTFs</i>. (See the <a href="http://www.ploscompbiol.org/article/info:doi/10.1371/journal.pcbi.1002417#s4" target="_blank">Materials and Methods</a> section for Equation 7 and its derivation.)</p