Search tips
Search criteria 


Logo of procbThe Royal Society PublishingProceedings BAboutBrowse by SubjectAlertsFree Trial
Proc Biol Sci. 2016 March 16; 283(1826): 20152563.
PMCID: PMC4810848

Ecological genomics of mutualism decline in nitrogen-fixing bacteria


Anthropogenic changes can influence mutualism evolution; however, the genomic regions underpinning mutualism that are most affected by environmental change are generally unknown, even in well-studied model mutualisms like the interaction between legumes and their nitrogen (N)-fixing rhizobia. Such genomic information can shed light on the agents and targets of selection maintaining cooperation in nature. We recently demonstrated that N-fertilization has caused an evolutionary decline in mutualistic partner quality in the rhizobia that form symbiosis with clover. Here, population genomic analyses of N-fertilized versus control rhizobium populations indicate that evolutionary differentiation at a key symbiosis gene region on the symbiotic plasmid (pSym) contributes to partner quality decline. Moreover, patterns of genetic variation at selected loci were consistent with recent positive selection within N-fertilized environments, suggesting that N-rich environments might select for less beneficial rhizobia. By studying the molecular population genomics of a natural bacterial population within a long-term ecological field experiment, we find that: (i) the N environment is indeed a potent selective force mediating mutualism evolution in this symbiosis, (ii) natural variation in rhizobium partner quality is mediated in part by key symbiosis genes on the symbiotic plasmid, and (iii) differentiation at selected genes occurred in the context of otherwise recombining genomes, resembling eukaryotic models of adaptation.

Keywords: symbiosis, mutualism, microbe, deposition, adaptation, bacteria

1. Introduction

Mutualisms play critical roles in natural and managed ecosystems. Recent work on context-dependency in mutualisms has shown that changes in the economy of benefits, such as a shift in the availability of an important traded resource in the environment, can alter the symbiont community composition and/or the forms of natural selection acting on mutualists, even causing the evolution of decreased partner quality or abandonment of the mutualism [15]. Despite these recent advances, however, we generally lack empirical evidence for the selective agents that act on mutualist partner quality and thus how partner quality (co)evolves in natural mutualist populations [6,7]. Moreover, few mutualist systems to date have afforded both the ecological relevance and genomic resolution necessary to study the mechanistic underpinnings of mutualism evolution.

The legume–rhizobium symbiosis, in which leguminous plants house nitrogen (N)-fixing bacteria in root nodules, is responsible for the vast majority of non-anthropogenically fixed N in terrestrial systems and thus plays a major role in the N cycle [8,9]. This symbiosis has become a key model for mutualism evolution, given the relative ease with which rhizobium genetic variants can be isolated and subsequently manipulated in experiments to test key predictions of evolutionary theory [1014]. This recent boom in studies of phenotypic evolution in rhizobia, a long history of molecular genetic investigation resulting in well-annotated reference genomes (e.g. [15,16]), and most recently the ability to collect high-quality population genomic data with relative ease (e.g. [1720]), all make N-fixing rhizobia ideal systems for addressing both mechanistic questions of mutualism evolution and more general questions about the mode of bacterial adaptation.

In a recent quantitative genetic study of Rhizobium leguminosarum, rhizobium symbionts of clover, we demonstrated that rhizobia from field soil populations fertilized with N for 22 years were, on average, less beneficial for their host plants, compared to those from control populations [5]. In particular, rhizobium strains isolated from the N-fertilized plots resulted in host plants with 17–30% less biomass, 10–28% fewer leaves, and 6–17% lower leaf chlorophyll (depending on the clover host species), supporting a major prediction of mutualism theory—that changes in the resource environment can lead to the evolutionary decline of mutualism benefits [2123]. The observed phenotypic differentiation in response to long-term experimental N-fertilization, however, does not address the rhizobium loci targeted by selection or the evolutionary mechanisms driving reduced rhizobium quality. For example, although numerous genes required for symbiosis with legumes have been characterized [24], whether those same loci underlie naturally occurring partner quality variation remains largely unknown [25]. Similarly, partner quality could decline in high N environments if (i) positive selection has favoured lower quality rhizobium strains (i.e. rhizobia have adapted to sustained high N), or (ii) selection for high-quality strains has been relaxed (i.e. strains have simply become defective mutualists in response to sustained high N). These alternative hypotheses point out a critical distinction between defective mutualists, which are simply less beneficial for their partners, versus true ‘cheaters’, which gain a fitness benefit from being less beneficial [26,27]. The evolutionary importance of low-quality partners to the dynamics of contemporary mutualism dynamics remains controversial [14,26,2831]. Thus, identifying whether selection favours cheating can shed light on this important issue in mutualism evolution.

Using natural populations of mutualist bacteria that have become phenotypically differentiated in response to a manipulated variable (N fertilization) provides the rare opportunity to address the genetic underpinnings of mutualism evolution in response to long-term changes in the resource environment in a well-replicated manner. Here, we compare the genome sequences of N-fertilized and control strains of R. leguminosarum to investigate patterns of genomic variation between N-fertilization treatments to identify the genomic basis of mutualism decline in rhizobia. Next, we ask whether patterns of nucleotide variation at putatively selected genes suggest adaptation of rhizobium mutualists, versus relaxed selection for mutualism, in high-resource environments. Our study merges long-term ecological field experiments and bacterial population genomics to answer key questions about evolutionary mechanisms in nature.

2. Material and methods

(a) Study system

Here, we studied genomic variation among 63 strains of R. leguminosarum bv. trifolii, N-fixing rhizobial mutualists of clover (Trifolium spp.), that were isolated and characterized phenotypically by Weese et al. [5]; complete methods detailing the long-term ecological experiment, strain isolation, and phenotypic experiments are contained there. Briefly, rhizobium strains were isolated from both N-fertilized (22 years of fertilization with 12.3 g N m−2 per year granular ammonium nitrate) and control (unfertilized) plots located at the Kellogg Biological Station Long Term Ecological Research Site (KBS LTER; Weese et al. [5] found that rhizobium populations from the N-fertilization treatment had evolved to be lower quality partners for their hosts (i.e. resulted in lower plant above-ground biomass and chlorophyll content), compared with rhizobium populations from control plots (see electronic supplementary material, Dataset S1).

(b) Whole genome sequencing and annotation

Strain isolates were grown on solid tryptone yeast (TY) media plates, incubated at 30°C until bacterial colonies became visible (3–4 days), and isolated strains were grown in 5 ml TY liquid media [32] in a shaking incubator for 1–2 days at 30°C. Genomic DNA was extracted using the MasterPure™ DNA Purification Kit (Epicentre). We used a NanoDrop™ 2000c (Thermo Scientific, Waltham, MA, USA) to ensure that all samples had A260/280 readings between 1.8 and 2.0. We generated multiplexed libraries (no bead normalization step) using the Nextera XT DNA Sample Prep Kit (Illumina, San Diego, CA, USA). Libraries were quantified with a Qubit fluorometer 2.0 (Life Technologies, Carlsbad, CA, USA) and then diluted with molecular grade water for normalization. Normalized libraries were pooled for paired end sequencing in a single lane on a HiSeq2000 sequencer (Illumina) at the W.M. Keck Center for Comparative and Functional Genomics at the University of Illinois at Urbana-Champaign. The lane produced 467 million reads, resulting in an average of 720 megabase (MB) and 97× coverage per genome.

We used the A5 pipeline [33] with default parameters to automate sequence data cleaning, error correction, de novo assembly, and quality control. Genomes were assembled de novo rather than using reference-based assembly to avoid potential biases in synteny and to better assess highly variable (i.e. flexible) genomic content often found in bacteria [20,34,35]. A single control strain (155_C) could not be assembled properly due to poor sequence quality and was excluded from subsequent analyses. A portion of gaps were filled in silico with GapFiller v. 1.11, which re-uses paired reads to fill gaps left in draft assemblies [36]. We annotated the genomes by uploading A5 assemblies to the ‘Rapid Annotation using Subsystem Technology’ server [37]. Genes from annotated assemblies were grouped into recognized homologues using the ITEP toolkit [38], which uses the MCL program (Markov Chain Clustering) to cluster BLASTp results (0.4 cut-off, 2.0 inflation parameter [39]). Assembled genomes had a median length of 7.50 Mbp, a median scaffold number of 45, and a median N50 of 465 604 bp (electronic supplementary material, Dataset S2).

We defined the ‘core genome’ as orthologous, single-copy ITEP gene clusters found within all strains and the reference R. leguminosarum bv. trifolii strain WSM1325 [15]. The remaining (non-core) genes comprise the ‘flexible genome’ (genes with variable copy number and/or shared among subsets of strains). Analysis of the flexible genome revealed high variation in gene content among strains, but no candidates strongly associated with partner quality decline (i.e. no genes were completely absent in control strains and present in all N-fertilized strains or vice versa). See electronic supplementary material, Dataset S3 and text S1, for full discussion of presence-absence variation (PAV) results. Using an ITEP pipeline, core gene sequences were aligned and reverse-translated with the auto flag in mafft v. 7.123b [40] and PAL2NAL v. 14 [41], respectively. The ITEP GBlocks wrapper script was used for automated alignment curation (-c 24 -f 24 -n 10 -m 5 -g a) for each core gene [42]. There were a total of 2 433 core (i.e. single copy, present in all strains) genes among all 62 strains. For analyses using the 56 strains harbouring symbiosis genes on the symbiotic plasmid (pSym) (see below), there were a total of 3 173 core genes.

(c) Genome-wide patterns of variation

To determine evolutionary relationships among the experimental strains, we concatenated aligned core gene sequences on the chromosome and five plasmids, including the pSym. We used the GTR-Γ model with 100 bootstrap replicates in RAxML v. 8.0.19 [43] to generate a core gene phylogeny and visualized the distribution of partner quality on the tree by mapping mean trait values using iTOL v. 2 [44]. To test whether nucleotide variation in the rhizobium genome was structured geographically and/or between N fertilization treatments, we used Arlequin v. [45] to carry out separate hierarchical analyses of molecular variance (AMOVA) for core gene single nucleotide polymorphisms (SNPs) on the chromosome (2 521 single-copy genes present in all strains) and the core plasmid genomes (238, 195, 8, 113, and 98 genes in plasmids pR132501–pR132505, respectively). Strain, plot within treatment (six plots per N treatment), and N treatment were included as factors. The first 100 000 SNP sites were used for AMOVA on the chromosome due to data input limitations of Arlequin. AMOVA on remaining chromosomal SNP sites (data not shown) did not alter the results.

(d) Patterns of gene loss at a key symbiosis gene region of the pSym

Many genes known to be critical for the formation and maintenance of symbiosis with legumes [46] exist within a symbiosis island roughly centred at 0.303 Mega base pair (Mbp) on the pSym of the reference strain, R. leguminosarum bv. trifolii WSM1325 [15], and we were interested in testing whether genes at this region were differentiated between N-fertilized and control rhizobia. We found that six strains (four N-fertilized strains, 40, 538, 643, and 773 and two control strains, 160 and 3) lacked many of the genes in this region and that the average partner quality for these six strains was below the average for all other strains (significantly so for plant biomass; t = −2.048, p = 0.045; electronic supplementary material, Dataset S2). This left a large gap in the core gene alignment at this critical region of interest (‘Results’, and gap on the pSym in electronic supplementary material, figure S1 and figure 1 of text S1); therefore, we used the 56 strains that possessed the region of the pSym in population genetic analyses of the core genome to test for patterns of differentiation at the symbiosis gene region of the pSym.

Figure 1.
Patterns of genetic structure (FST) and gene presence and absence across the genome for N-fertilized and control R. leguminosarum. (a) Gene-by-gene FST between N-fertilized and control rhizobia for all (3 173) core genes using the subset of 56 strains ...

(e) Genetic differentiation between N-fertilized and control strains

We expected that genes underlying phenotypic differentiation between N-fertilized and control rhizobia should be genetically differentiated (have elevated FST) between N-fertilized and control groups, and also be associated with symbiotic partner quality (i.e. differentiate low- and high-quality strains). Thus, we calculated FST for core genes in two ways. First, we identified genes showing elevated FST between N-fertilized and control strains (i.e. between the groups of rhizobia originally collected from N-fertilized versus control plots). Second, we ordered strains according to rank partner quality (rather than N treatment) by summing the ranks of all plant growth traits (see ranks in electronic supplementary material, Dataset S2) and split the strains into two groups: ‘high-quality’ (eight N-fertilized and 20 control strains) versus ‘low-quality’ (16 N-fertilized and 12 control strains). We then identified genes showing elevated FST between these two quality groups. In conjunction with the ITEP toolkit, we used a collection of ‘pop_genome’ scripts [47] ( to carry out population genetic analyses with core gene alignments. Population pairwise FST was calculated with default parameters for each core gene using Arlecore v. [45]. The 99th percentile of FST values was calculated with the ‘quantile’ function in R v. 3.0.3 (R Core Team, 2014; We performed a sliding window analysis in R to calculate the mean of pairwise FST values using a window size of 13 genes (to roughly match the size of the region spanning the between-treatment differentiated symbiosis genes; see ‘Results’) and a step-size of seven genes (452 windows total).

(f) Patterns of nucleotide variation at candidate selected genes

Forms of selection acting differently in N-fertilized (versus control) plots should generate distinct patterns of nucleotide variation at candidate genes. If positive directional selection has favoured low-quality strains in the N plots, for example, N-fertilized strains should have reduced nucleotide diversity (π) at candidate genes compared to control strains, relative to the rest of the genome [4851]. Alternatively, relaxed selection in the N plots should lead to increased π at candidate genes compared to control strains. We calculated π within N-fertilized and control populations for each core gene on the chromosome or pSym using ProSeq v. 3.5 [52] and then performed a sliding window analysis of π as described above (394 windows total). The other four plasmids were excluded due to a relative deficiency of core genes with significantly elevated FST. R v. 3.0.3 was used to plot the regression of sliding window results for N-fertilized and control populations, and to generate the 95% confidence and prediction intervals to identify regions of the genome (sliding windows) with anomalous π.

3. Results and discussion

We found that the vast majority of genetic variation occurred among individual strains (AMOVA; electronic supplementary material, table S1), consistent with other bacterial populations that are quite variable even at fine spatial scales [53,54]. We also detected spatial genetic structure (significant among-plot variance; electronic supplementary material, table S1), suggesting some limits to gene flow at a rather fine geographical scale in these rhizobia [55,56]. We found little evidence of genome-wide differentiation between N-fertilized and control rhizobia, as evidenced first by an intermixed core genome phylogeny (electronic supplementary material, figure S2), and second by the lack of between-treatment variance in the AMOVAs (electronic supplementary material, table S1A). This finding is consistent with a previous analysis of the internal transcribed spacer (ITS) locus [5] and facilitates the identification of genomic outliers that are particularly structured due to differential selection in N-fertilized versus control plots [57,58]. We posited that, compared to the rest of the genome, rhizobium genes underlying mutualism decline should: (i) be genetically differentiated between N-fertilized and control rhizobia, (ii) be associated with partner quality phenotypes (rhizobium effects on plant growth and chlorophyll content), and (iii) show anomalous patterns of nucleotide diversity consistent with recent selection.

(a) Genetic differentiation between N-fertilized versus control and high- versus low-quality rhizobia

The top 1% of FST values in our comparison of N-fertilized and control rhizobium populations (FST ≥ 0.074) included 20 chromosomal genes, nine pSym genes, and a single gene on both plasmid pR132502 and pR132505 (electronic supplementary material, Dataset S4). Twenty-five of these genes were significantly differentiated (17 chromosomal genes: mean FST = 0.090, s.d. = 0.016; eight pSym genes: mean FST = 0.125, s.d. = 0.037). These outliers stand out against a backdrop of low overall genomic differentiation (mean FST = 0.022, s.d. = 0.019, 99th percentile = 0.074), and the sliding window analysis revealed a particularly strong FST peak at the symbiosis gene region of the pSym (blue vertical line in figure 1a). Genes among the top 1% of FST values were significantly enriched for pSym genes (32% of genes in the top 1% were on the pSym, versus the null expectation of 8.8% from the core genome, two-tailed Fisher's exact test, p = 0.001). Evidence from flexible genome analyses (i.e. genes with variable copy number and/or presence among strains) also revealed that six relatively low-quality strains were missing several symbiosis genes on the pSym (see pSym gap in electronic supplementary material, figure S1), further implicating the symbiosis island as an area of interest.

Consistent with the above findings, genes in the top 1% of FST values comparing low versus high partner quality rhizobia were also enriched for pSym genes (48% on the pSym, p < 0.0001), and sliding window analysis revealed an island of high differentiation at the symbiosis region on the pSym (figure 1b), compared to a low genomic background (mean FST = 0.027, s.d. = 0.022, 99th percentile = 0.086). Specifically, the top 1% of FST values (FST ≥ 0.086; 34 genes) in this high- versus low-quality comparison included 19 on the chromosome (14 significant; mean FST = 0.100, s.d. = 0.017) and 14 on the pSym (13 significant; mean FST = 0.189, s.d. = 0.053), and a single gene on plasmid pR132502 (p > 0.05; electronic supplementary material, Dataset S4).

Putting these two analyses together, only six genes were among the most differentiated (top 1%) in both FST comparisons: nifH, nifA, nifD, fixC, nodB, and an electron transfer flavoprotein subunit (Rleg_4928)—all on the pSym (blue vertical line in figure 1). Together, these results implicate genes in this region as players in the evolution of reduced cooperation in response to N fertilization. A strong history of molecular genetic work using rhizobial mutants and changes in gene expression during various stages of symbiosis has revealed the necessity of many rhizobium genes (particularly those in this region) for the establishment and maintenance of symbiosis with legumes (e.g. [15,5964]). Ecological genetic studies such as ours, which focus on natural variation and the loci underlying ecologically relevant traits, can complement mutant analyses by helping to resolve the loci targeted by contemporary natural selection [6,7,25,65].

(b) Patterns of nucleotide diversity in N-fertilized versus control strains

The genome-wide average π was higher among N-fertilized strains compared with control strains, resulting from a few N-fertilized strains (209, 231, and 717) that dominated a clade more distantly related to the majority of other strains (electronic supplementary material, figure S2). Nevertheless, π covaried between the two groups for both the chromosomal and pSym windows (chromosome R = 0.81, p < 0.0001; pSym R = 0.56, p < 0.001; figure 2b,c), suggesting that the relative selective constraints among genes are similar across most of the genome for N-fertilized and control groups. Two chromosomal sliding window regions fell below the 95% prediction interval in the π correlation (figure 2b), suggesting positive selection in N-fertilized environments; however, these regions did not possess genes with FST in the top 1% of either FST analysis. Two pSym sliding window regions also fell below the 95% prediction interval (blue triangles in figure 2c), and these regions contained the six pSym loci from the top 1% of the FST outlier analyses above. Reduced nucleotide diversity at the six pSym genes revealed in the FST outlier analyses suggests positive selection on these genes in N-fertilized environments. In other words, these results are consistent with natural selection favouring less beneficial rhizobia in high N conditions. Manipulative experiments are required to test definitively whether low-quality N strains have higher fitness than control strains in high N environments; these experiments are underway.

Figure 2.
Patterns of within-treatment nucleotide diversity (π) across the chromosome and symbiotic plasmid (pSym) for 56 R. leguminosarum strains. (a) Points (N = 394) represent sliding window means (window size = 13 genes), and dashed lines represent ...

One selective agent favouring lower partner quality in high N plots might be a shift towards the saprophytic lifestyle at the expense of symbiosis. Among other ecological changes in the Long Term Ecological Research (LTER) study since N fertilization began, the availability of clover hosts has declined in N-fertilized plots [5,66], likely leading to more frequent/prolonged saprobic phases for R. leguminosarum. Genetic trade-offs between symbiotic and free-living lifestyles could drive positive selection for less beneficial rhizobia in high N environments. For example, selection could favour mutualists that invest more in storage compounds such as poly-3-hydroxybutyurate (PHB), which can increase rhizobium fitness in the soil but likely divert carbon away from N fixation and thus from host benefits [67]. Both molecular and experimental evolution studies suggest the potential for such a trade-off between saprotrophy and mutualism. Rhizobium etli (PHB synthase) and R. tropici (glycogen synthase) mutants have higher rates of nitrogen fixation in symbiosis and hoard less carbon for life outside the nodule [68,69], while Bradyrhizobium japonicum strains that evolved the most significant increase in host-free fitness after 450 generations in culture showed the largest decline in partner quality [70]. We are currently investigating whether the N-fertilized strains in our study are better saprophytes.

(c) Potential functions of key candidate genes

Only six symbiosis-related genes on the pSym (nifH, nifA, nifD, fixC, nodB, and an electron transfer flavoprotein subunit) matched all three criteria for identifying the genes underlying the response to N fertilization. Because these genes are in close proximity, it is unclear which particular genes have been targeted by selection; nevertheless, given the large body of genetic work on symbiosis genes in this region [46], several interesting candidates stand out. For example, nifH, nifD, nifA, nodB, and fixC are among the top 1% of FST values in both FST analyses, and additionally nifK, nifE, nifN, and nodF are in the top 1% of FST values in the second analysis of high- and low-quality groups (a full list of candidate loci is in the electronic supplementary material, Dataset S4). The N-fixing (nif) genes nifH, nifD, and nifK are responsible for the formation of the nitrogenase enzymatic complex (Fe and FeMo proteins), which reduces atmospheric di-nitrogen (N2) to ammonia for plant use [71]. The nifB, nifE, and nifN genes are involved in the synthesis of the FeMo cofactor that is located at the active site of the FeMo protein, where N2 binds and is reduced [72,73]. Additionally, nifA helps regulate the expression of various nif genes [74], while fixC is part of the fixABCX complex involved in gene transcription regulation in low oxygen conditions (e.g. inside the nodule) and is posited to have a central role in electron transfer to nitrogenase [75,76]. In general, nod genes code for enzymes that generate Nod factors, rhizobium signalling molecules between rhizobia and plant required for nodulation initiation; nodB in particular is part of the nodABC operon, which generates the main Nod factor structure [24].

Aside from the six symbiosis-related pSym genes discussed above, 12 window regions on the chromosome and one window on the pSym possessed decreased π in control, compared with N-fertilized, strains (and thus fall above the upper bounds of the prediction interval in figure 2b,c). There were no outliers from both FST analyses within these 13 windows; however, they did contain five chromosomal genes that were significantly differentiated between N-fertilized and control rhizobia (electronic supplementary material, Dataset S4), including adenylate/guanylate cyclase and FecR. It is unclear how these genes might govern partner quality traits per se, since they did not differentiate low- and high-quality rhizobia in the second FST comparison; nevertheless, they might have an important role in uncharacterized phenotypic differences between N- and control strains. Adenylate cyclases have recently been implicated in infection thread formation [77], as well as in the negative regulation of secondary infections after initial infection thread formation [78], potentially controlling the number of successful nodulation events. N-fixing rhizobia inside the nodule have a high demand for iron (e.g. nitrogenase and cytochrome synthesis), and iron uptake could occur by means of ferric citrate [79], which requires the cytoplasmic membrane protein FecR in Escherichia coli [80]. A full understanding of how chromosomal and pSym genes together determine mutualist partner quality and other rhizobium phenotypes will require a much larger effort, ideally implementing emerging techniques for association mapping in bacteria [81,82].

(d) Bacterial adaptation in nature

Bacterial symbionts are key players in many plant and animal mutualisms [8385], yet the historical dearth of ecological genetic studies in bacteria has meant that the rules governing adaptation in natural populations of bacteria are still being elucidated [86]. The balance between natural selection and recombination within natural bacterial populations is one important unresolved aspect. With little recombination, selection acting at one or more loci should result in genome-wide sweeps that purge previously existing genetic variation within the population [19,87]. Many recent microbial studies, however, have discovered adaptive regions of reduced genetic diversity amidst high background polymorphism [18,8891]—implying abundant recombination. Using a population genomic approach, we found a genomic island of differentiation between rhizobium populations that have experienced different fertilization histories and that this region stands out against a genome-wide background of recombination/admixture between the two groups. These findings are congruent with a recent study of a natural R. leguminosarum population from the UK, in which high rates of recombination sufficient to prevent divergence by drift were found [20]. Thus R. leguminosarum in nature appear to form dynamic, diverse populations that are unified by gene flow despite selection acting at one or more loci. In our case, where genomic islands of differentiation (this study) and large phenotypic differences [5] were detected, selection due to N addition must be exceptionally strong.

Because the rhizobium populations were not assayed when the LTER experiment began, we cannot know with certainty whether particular nucleotide-level differences underlying our phenotypic differentiation have arisen through mutation since that time (22 years ago), or instead pre-dated the establishment of the N plots. These alternatives represent two distinct mechanisms of evolution [92]. On the one hand, new mutations could have arisen and increased in frequency atop a diverse phylogeny of rhizobium strains (selection on new mutations); on the other hand, pre-existing alleles could have increased in frequency (selection on standing genetic variation). Given the levels of nucleotide diversity found throughout the genome, in addition to our observations of substantial standing variation in partner quality variation within the control plots [5], we suspect that the alleles at the key outlier region on the pSym pre-dated the establishment of the N plots—i.e. that selection has acted on standing genetic variation. Evolution from standing variation should proceed more quickly than evolution from new mutations [93,94], potentially allowing populations to adapt more quickly. Indeed, how quickly N-fertilized populations of rhizobia might return to control levels of partner quality after a return to low N conditions is an interesting open question.

Supplementary Material

Table S1:

Supplementary Material

Figure S1:

Supplementary Material

Figure S2:

Supplementary Material

Text S1: Assessing Gene Presence Absence Patterns:

Supplementary Material

Datasets S1-S4:


We thank B. Gordon for laboratory assistance, N. Youngblut and M. Friesen for sequencing and bioinformatics advice, and R. Whitaker, Z. Cheviron, P. Tiffin, and the Heath laboratory for comments on the manuscript. This is KBS publication no. 1903.

Data accessibility

Illumina reads and de novo genome assemblies are available via DRYAD (

Competing interests

We declare we have no competing interests.


This work was funded by the National Science Foundation (NSF DEB-1257938), the NSF LTER Program at the Kellogg Biological Station (DEB-1027253), and Michigan State University AgBioResearch.


1. Sachs JL, Simms EL 2006. Pathways to mutualism breakdown. Trends Ecol. Evol. 21, 585–592. (doi:10.1016/j.tree.2006.06.018) [PubMed]
2. Heath KD, Stock AJ, Stinchcombe JR 2010. Mutualism variation in the nodulation response to nitrate. J. Evol. Biol. 23, 2494–2500. (doi:10.1111/j.1420-9101.2010.02092.x) [PubMed]
3. Johnson NC. 2010. Resource stoichiometry elucidates the structure and function of arbuscular mycorrhizas across scales. New Phytol. 185, 631–647. (doi:10.1111/j.1469-8137.2009.03110.x) [PubMed]
4. Kiers ET, Palmer TM, Ives AR, Bruno JF, Bronstein JL 2010. Mutualisms in a changing world: an evolutionary perspective. Ecol. Lett. 13, 1459–1474. (doi:10.1111/j.1461-0248.2010.01538.x) [PubMed]
5. Weese DJ, Heath KD, Dentinger BTM, Lau JA 2015. Long-term nitrogen addition causes the evolution of less-cooperative mutualists. Evolution 69, 631–642. (doi:10.1111/evo.12594) [PubMed]
6. Hoeksema JD. 2010. Ongoing coevolution in mycorrhizal interactions. New Phytol. 187, 286–300. (doi:10.1111/j.1469-8137.2010.03305.x) [PubMed]
7. Heath KD, Stinchcombe JR 2014. Explaining mutualism variation: a new evolutionary paradox? Evolution 68, 309–317. (doi:10.1111/evo.12292) [PubMed]
8. Sprent JI, Sutherland JM, De Faria SM, Dilworth MJ, Corby HDL, Becking JH, Materon LA, Drozd JW 1987. Some aspects of the biology of nitrogen-fixing organisms. Phil. Trans. R. Soc. Lond. B 317, 111–129. (doi:10.1098/rstb.1987.0051)
9. Cleveland CC, et al. 1999. Global patterns of terrestrial biological nitrogen (N2) fixation in natural ecosystems. Glob. Biochem. Cycles 13, 623–645. (doi:10.1029/1999GB900014)
10. Simms EL, Taylor DL, Povich J, Shefferson RP, Sachs JL, Urbina M, Tausczik Y 2006. An empirical test of partner choice mechanisms in a wild legume-rhizobium interaction. Proc. R. Soc. B 273, 77–81. (doi:10.1098/rspb.2005.3292) [PMC free article] [PubMed]
11. Heath KD, Tiffin P 2009. Stabilizing mechanisms in a legume–rhizobium mutualism. Evolution 63, 652–662. (doi:10.1111/j.1558-5646.2008.00582.x) [PubMed]
12. Barrett LG, Broadhurst LM, Thrall PH 2012. Geographic adaptation in plant-soil mutualisms: tests using Acacia spp. and rhizobial bacteria. Funct. Ecol. 26, 457–468. (doi:10.1111/j.1365-2435.2011.01940.x)
13. Regus JU, Gano KA, Hollowell AC, Sachs JL 2014. Efficiency of partner choice and sanctions in Lotus is not altered by nitrogen fertilization. Proc. R. Soc. B 281, 20132587 (doi:10.1098/rspb.2013.2587) [PMC free article] [PubMed]
14. Porter SS, Simms EL 2014. Selection for cheating across disparate environments in the legume-rhizobium mutualism. Ecol. Lett. 17, 1121–1129. (doi:10.1111/ele.12318) [PubMed]
15. Reeve W, et al. 2010. Complete genome sequence of Rhizobium leguminosarum bv. trifolii strain WSM1325, an effective microsymbiont of annual Mediterranean clovers. Stand. Genomic Sci. 2, 347–356. (doi:10.4056/sigs.852027) [PMC free article] [PubMed]
16. Schneiker-Bekel S, et al. 2011. The complete genome sequence of the dominant Sinorhizobium meliloti field isolate SM11 extends the S. meliloti pan-genome. J. Biotechnol. 155, 20–33. (doi:10.1016/j.jbiotec.2010.12.018) [PubMed]
17. Cadillo-Quiroz H, Didelot X, Held NL, Herrera A, Darling A, Reno ML, Krause DJ, Whitaker RJ 2012. Patterns of gene flow define species of thermophilic Archaea. PLoS Biol. 10, e1001265 (doi:10.1371/journal.pbio.1001265) [PMC free article] [PubMed]
18. Shapiro BJ, Friedman J, Cordero OX, Preheim SP, Timberlake SC, Szabó G, Polz MF, Alm EJ 2012. Population genomics of early events in the ecological differentiation of bacteria. Science 336, 48–51. (doi:10.1126/science.1218198) [PMC free article] [PubMed]
19. Kopac S, Wang Z, Wiedenbeck J, Sherry J, Wu M, Cohan FM 2014. Genomic heterogeneity and ecological speciation within one subspecies of Bacillus subtilis. Appl. Environ. Microbiol. 80, 4842–4853. (doi:10.1128/AEM.00576-14) [PMC free article] [PubMed]
20. Kumar N, Lad G, Giuntini E, Kaye ME, Udomwong P, Shamsani NJ, Young JPW, Bailly X 2015. Bacterial genospecies that are not ecologically coherent: population genomics of Rhizobium leguminosarum. Open Biol. 5, 140 133 (doi:10.1098/rsob.140133) [PMC free article] [PubMed]
21. West SA, Kiers ET, Simms EL, Denison RF 2002. Sanctions and mutualism stability: why do rhizobia fix nitrogen? Proc. R. Soc. Lond. B 269, 685–694. (doi:10.1098/rspb.2001.1878) [PMC free article] [PubMed]
22. Denison RF, Kiers ET 2004. Lifestyle alternatives for rhizobia: mutualism, parasitism, and forgoing symbiosis. FEMS Microbiol. Lett. 237, 187–193. (doi:10.1016/j.femsle.2004.07.013) [PubMed]
23. Akçay E, Simms EL 2011. Negotiation, sanctions, and context dependency in the legume-Rhizobium mutualism. Am. Nat. 178, 1–14. (doi:10.1086/659997) [PubMed]
24. Jones KM, Kobayashi H, Davies BW, Taga ME, Walker GC 2007. How rhizobial symbionts invade plants: the Sinorhizobium-Medicago model. Nat. Rev. Microbiol. 5, 619–633. (doi:10.1038/nrmicro1705) [PMC free article] [PubMed]
25. Heath KD, Burke PV, Stinchcombe JR 2012. Coevolutionary genetic variation in the legume–rhizobium transcriptome. Mol. Ecol. 21, 4735–4747. (doi:10.1111/j.1365-294X.2012.05629.x) [PubMed]
26. Ghoul M, Griffin AS, West SA 2014. Toward an evolutionary definition of cheating. Evolution 68, 318–331. (doi:10.1111/evo.12266) [PubMed]
27. Jones EI, et al. 2015. Cheaters must prosper: reconciling theoretical and empirical perspectives on cheating in mutualism. Ecol. Lett. 18, 1270–1284. (doi:10.1111/ele.12507) [PubMed]
28. Kiers ET, Denison RF, Kawakita A, Herre EA 2011. The biological reality of host sanctions and partner fidelity. Proc. Natl Acad. Sci. USA 108, E7 (doi:10.1073/pnas.1014546108) [PubMed]
29. Friesen ML. 2012. Widespread fitness alignment in the legume-rhizobium symbiosis. New Phytol. 194, 1096–1111. (doi:10.1111/j.1469-8137.2012.04099.x) [PubMed]
30. Frederickson ME. 2013. Rethinking mutualism stability: cheaters and the evolution of sanctions. Q. Rev. Biol. 88, 269–295. (doi:10.1086/673757) [PubMed]
31. Friesen ML, Heath KD 2013. One hundred years of solitude: integrating single-strain inoculations with community perspectives in the legume–rhizobium symbiosis. New Phytol. 198, 7–9. (doi:10.1111/nph.12173) [PubMed]
32. Somasegaran P, Hoben HJ 1994. Handbook for rhizobia: methods in legume-Rhizobium technology, p. 510. New York: Springer-Verlag.
33. Tritt A, Eisen JA, Facciotti MT, Darling AE 2012. An integrated pipeline for de novo assembly of microbial genomes. PLoS ONE 7, e42304 (doi:10.1371/journal.pone.0042304) [PMC free article] [PubMed]
34. Welch RA, et al. 2002. Extensive mosaic structure revealed by the complete genome sequence of uropathogenic Escherichia coli. Proc. Natl Acad. Sci. USA 99, 17 020–17 024. (doi:10.1073/pnas.252529799) [PubMed]
35. Boucher Y, Cordero OX, Takemura A 2011. Endemicity within global Vibrio cholerae populations. MBio 2, 1–8. (doi:10.1128/mBio.00335-10.Editor) [PMC free article] [PubMed]
36. Boetzer M, Pirovano W 2012. Toward almost closed genomes with GapFiller. Genome Biol. 13, R56 (doi:10.1186/gb-2012-13-6-r56) [PMC free article] [PubMed]
37. Aziz RK, et al. 2008. The RAST server: rapid annotations using subsystems technology. BMC Genomics 9, 75 (doi:10.1186/1471-2164-9-75) [PMC free article] [PubMed]
38. Benedict MN, Henriksen JR, Metcalf WW, Whitaker RJ, Price ND 2014. ITEP: an integrated toolkit for exploration of microbial pan-genomes. BMC Genomics 15 8 (doi:10.1186/1471-2164-15-8) [PMC free article] [PubMed]
39. Enright AJ, Van Dongen S, Ouzounis CA 2002. An efficient algorithm for large-scale detection of protein families. Nucleic Acids Res. 30, 1575–1584. (doi:10.1093/nar/30.7.1575) [PMC free article] [PubMed]
40. Katoh K, Standley DM 2013. MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol. Biol. Evol. 30, 772–780. (doi:10.1093/molbev/mst010) [PMC free article] [PubMed]
41. Suyama M, Torrents D, Bork P 2006. PAL2NAL: robust conversion of protein sequence alignments into the corresponding codon alignments. Nucleic Acids Res. 34, W609–W612. (doi:10.1093/nar/gkl315) [PMC free article] [PubMed]
42. Talavera G, Castresana J 2007. Improvement of phylogenies after removing divergent and ambiguously aligned blocks from protein sequence alignments. Syst. Biol. 56, 564–577. (doi:10.1080/10635150701472164) [PubMed]
43. Stamatakis A. 2014. RAxML version 8: a tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics 30, 1312–1313. (doi:10.1093/bioinformatics/btu033) [PMC free article] [PubMed]
44. Letunic I, Bork P 2011. Interactive tree of Life v2: online annotation and display of phylogenetic trees made easy. Nucleic Acids Res. 39, W475–W478. (doi:10.1093/nar/gkr201) [PMC free article] [PubMed]
45. Excoffier L, Lischer HEL 2010. Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Mol. Ecol. Resour. 10, 564–567. (doi:10.1111/j.1755-0998.2010.02847.x) [PubMed]
46. Black M, Moolhuijzen P, Chapman B, Barrero R, Howieson J, Hungria M, Bellgard M 2012. The genetics of symbiotic nitrogen fixation: comparative genomics of 14 rhizobia strains by resolution of protein clusters. Genes 3, 138–166. (doi:10.3390/genes3010138) [PMC free article] [PubMed]
47. Youngblut ND, Wirth JS, Henriksen JR, Smith M, Simon H, Metcalf WW, Whitaker RJ 2015. Genomic and phenotypic differentiation among Methanosarcina mazei populations from Columbia River sediment. ISME J. 9, 2191–2205. (doi:10.1038/ismej.2015.31) [PMC free article] [PubMed]
48. Lewontin RC, Krakauer J 1973. Distribution of gene frequency as a test of the theory of the selective neutrality of polymorphisms. Genetics 74, 175–195. [PubMed]
49. Beaumont MA, Nichols RA 1996. Evaluating loci for use in the genetic analysis of population structure. Proc. R. Soc. Lond. B 263, 1619–1626. (doi:10.1098/rspb.1996.0237)
50. Kauer MO, Dieringer D, Schlötterer C 2003. A microsatellite variability screen for positive selection associated with the ‘out of Africa’ habitat expansion of Drosophila melanogaster. Genetics 165, 1137–1148. [PubMed]
51. Beaumont MA, Balding DJ 2004. Identifying adaptive genetic divergence among populations from genome scans. Mol. Ecol. 13, 969–980. (doi:10.1111/j.1365-294X.2004.02125.x) [PubMed]
52. Filatov DA. 2009. Processing and population genetic analysis of multigenic datasets with ProSeq3 software. Bioinformatics 25, 3189–3190. (doi:10.1093/bioinformatics/btp572) [PMC free article] [PubMed]
53. Tettelin H, Riley D, Cattuto C, Medini D 2008. Comparative genomics: the bacterial pan-genome. Curr. Opin. Microbiol. 11, 472–477. (doi:10.1016/j.mib.2008.09.006) [PubMed]
54. Mira A, Martín-Cuadrado AB, D'Auria G, Rodríguez-Valera F 2010. The bacterial pan-genome: a new paradigm in microbiology. Int. Microbiol. 13, 45–57. (doi:10.2436/20.1501.01.110) [PubMed]
55. Whitaker RJ, Grogan DW, Taylor JW 2003. Geographic barriers isolate endemic populations of hyperthermophilic archaea. Science 301, 976–978. (doi:10.1126/science.1086909) [PubMed]
56. Vos M, Velicer GJ 2008. Isolation by distance in the spore-forming soil bacterium Myxococcus xanthus. Curr. Biol. 18, 386–391. (doi:10.1016/j.cub.2008.02.050) [PubMed]
57. Beaumont MA. 2005. Adaptation and speciation: what can Fst tell us? Trends Ecol. Evol. 20, 435–440. (doi:10.1016/j.tree.2005.05.017) [PubMed]
58. Storz JF. 2005. Using genome scans of DNA polymorphism to infer adaptive population divergence. Mol. Ecol. 14, 671–688. (doi:10.1111/j.1365-294X.2004.02437.x) [PubMed]
59. Oke V, Long SR 1999. Bacteroid formation in the rhizobium–legume symbiosis. Curr. Opin. Microbiol. 2, 641–646. (doi:10.1016/S1369-5274(99)00035-1) [PubMed]
60. Kaneko T, et al. 2002. Complete genomic sequence of nitrogen-fixing symbiotic bacterium Bradyrhizobium japonicum USDA110. DNA Res. 9, 189–197. (doi:10.1093/dnares/9.6.189) [PubMed]
61. Long SR. 2001. Genes and signals in the rhizobium–legume symbiosis. Plant Physiol. 125, 69–72. (doi:10.1104/pp.125.1.69) [PubMed]
62. Udvardi M, Poole PS 2013. Transport and metabolism in legume–rhizobia symbioses. Annu. Rev. Plant Biol. 64, 781–805. (doi:10.1146/annurev-arplant-050312-120235) [PubMed]
63. Barnett MJ, Long SR 2015. The Sinorhizobium meliloti SyrM regulon: effects on global gene expression are mediated by syrA and nodD3. J. Bacteriol. 197, 1792–1806. (doi:10.1128/JB.02626-14) [PMC free article] [PubMed]
64. Becker A, et al. 2009. A portal for rhizobial genomes: RhizoGATE integrates a Sinorhizobium meliloti genome annotation update with postgenome data. J. Biotechnol. 140, 45–50. (doi:10.1016/j.jbiotec.2008.11.006) [PMC free article] [PubMed]
65. Rausher MD, Delph LF 2015. When does understanding phenotypic evolution require identification of the underlying genes? Evolution 69, 1655–1664. (doi:10.1111/evo.12687) [PubMed]
66. Huberty LE, Gross KL, Miller CJ 1998. Effects of nitrogen addition on successional dynamics and species diversity in Michigan old-fields. J. Ecol. 86, 794–803. (doi:10.1046/j.1365-2745.1998.8650794.x)
67. Ratcliff WC, Underbakke K, Denison RF 2011. Measuring the fitness of symbiotic rhizobia. Symbiosis 55, 85–90. (doi:10.1007/s13199-011-0150-2)
68. Cevallos MA, Encarnación S, Leija A, Mora Y, Mora J 1996. Genetic and physiological characterization of a Rhizobium etli mutant strain unable to synthesize poly-beta-hydroxybutyrate. J. Bacteriol. 178, 1646–1654. [PMC free article] [PubMed]
69. Marroquí S, Zorreguieta A, Santamaría C, Temprano F, Soberón M, Megías M, Downie JA 2001. Enhanced symbiotic performance by Rhizobium tropici glycogen synthase mutants. J. Bacteriol. 183, 854–864. (doi:10.1128/JB.183.3.854-864.2001) [PMC free article] [PubMed]
70. Sachs JL, Russell JE, Hollowell AC 2011. Evolutionary instability of symbiotic function in Bradyrhizobium japonicum. PLoS ONE 6, e26370 (doi:10.1371/journal.pone.0026370) [PMC free article] [PubMed]
71. Fischer H-M. 1994. Genetic regulation of nitrogen fixation in rhizobia. Microbiol. Rev. 58, 352–386. [PMC free article] [PubMed]
72. Raymond J, Siefert JL, Staples CR, Blankenship RE 2004. The natural history of nitrogen fixation. Mol. Biol. Evol. 21, 541–554. (doi:10.1093/molbev/msh047) [PubMed]
73. Hu Y, Ribbe MW 2013. Nitrogenase assembly. Biochim. Biophys. Acta Bioenerg. 1827, 1112–1122. (doi:10.1016/j.bbabio.2012.12.001) [PMC free article] [PubMed]
74. Dixon R, Kahn D 2004. Genetic regulation of biological nitrogen fixation. Nat. Rev. Microbiol. 2, 621–631. (doi:10.1038/nrmicro954) [PubMed]
75. Edgren T, Nordlund S 2004. The fixABCX genes in Rhodospirillum rubrum encode a putative membrane complex participating in electron transfer to nitrogenase. J. Bacteriol. 186, 2052–2060. (doi:10.1128/JB.186.7.2052) [PMC free article] [PubMed]
76. Sperotto RA, Gross J, Vedoy C, Passaglia LMP, Schrank IS 2004. The electron transfer flavoprotein fixABCX gene products from Azospirillum brasilense show a NifA-dependent promoter regulation. Curr. Microbiol. 49, 267–273. (doi:10.1007/s00284-004-4318-3) [PubMed]
77. Lomovatskaya LA, Romanenko AS, Rykun OV 2015. Transmembrane adenylate cyclase controls the virulence factors of plant pathogenic Pseudomonas siringae and mutualistic Rhizobium leguminosarum. Microbiology 84, 473–478. (doi:10.1134/S0026261715040116)
78. Tian CF, Garnerone A-M, Mathieu-Demaziere C, Masson-Boivin C, Batut J 2012. Plant-activated bacterial receptor adenylate cyclases modulate epidermal infection in the Sinorhizobium meliloti-Medicago symbiosis. Proc. Natl Acad. Sci. USA 109, 6751–6756. (doi:10.1073/pnas.1120260109) [PubMed]
79. Brear EM, Day DA, Smith PMC 2013. Iron: an essential micronutrient for the legume-rhizobium symbiosis. Front. Plant Sci. 4, 1–15. (doi:10.3389/fpls.2013.00359) [PMC free article] [PubMed]
80. Enz S, Mahren S, Stroeher UH, Braun V 2000. Surface signaling in ferric citrate transport gene induction: interaction of the FecA, FecR, and FecI regulatory proteins. J. Bacteriol. 182, 637–646. (doi:10.1128/JB.182.3.637-646.2000) [PMC free article] [PubMed]
81. Falush D, Bowden R 2006. Genome-wide association mapping in bacteria? Trends Microbiol. 14, 353–355. (doi:10.1016/j.tim.2006.06.003) [PubMed]
82. Sheppard SK, et al. 2013. Genome-wide association study identifies vitamin B5 biosynthesis as a host specificity factor in Campylobacter. Proc. Natl Acad. Sci. USA 110, 11 923–11 927. (doi:10.5061/dryad.28n35) [PubMed]
83. Bäckhed F, Ley RE, Sonnenburg JL, Peterson DA, Gordon JI 2005. Host-bacterial mutualism in the human intestine. Science 307, 1915–1920. (doi:10.1126/science.1104816) [PubMed]
84. Sachs JL, Essenberg CJ, Turcotte MM 2011. New paradigms for the evolution of beneficial infections. Trends Ecol. Evol. 26, 202–209. (doi:10.1016/j.tree.2011.01.010) [PubMed]
85. Sachs JL, Skophammer RG, Bansal N, Stajich JE 2014. Evolutionary origins and diversification of proteobacterial mutualists. Proc. R. Soc. B 281, 20132146 (doi:10.1098/rspb.2013.2146) [PMC free article] [PubMed]
86. Cordero OX, Polz MF 2014. Explaining microbial genomic diversity in light of evolutionary ecology. Nat. Rev. Microbiol. 12, 263–273. (doi:10.1038/nrmicro3218) [PubMed]
87. Cohan FM. 2006. Towards a conceptual and operational union of bacterial systematics, ecology, and evolution. Phil. Trans. R. Soc. B 361, 1985–1996. (doi:10.1098/rstb.2006.1918) [PMC free article] [PubMed]
88. Guttman DS, Dykhuizen DE 1994. Detecting selective sweeps in naturally occurring Escherichia coli. Genetics 138, 993–1003. [PubMed]
89. Denef VJ, Mueller RS, Banfield JF 2010. AMD biofilms: using model communities to study microbial evolution and ecological complexity in nature. ISME J. 4, 599–610. (doi:10.1038/ismej.2009.158) [PubMed]
90. Polz MF, Alm EJ, Hanage WP 2013. Horizontal gene transfer and the evolution of bacterial and archaeal population structure. Trends Genet. 29, 170–175. (doi:10.1016/j.tig.2012.12.006) [PMC free article] [PubMed]
91. Krause DJ, Didelot X, Cadillo-Quiroz H, Whitaker RJ 2014. Recombination shapes genome architecture in an organism from the archaeal domain. Genome Biol. Evol. 6, 170–178. (doi:10.1093/gbe/evu003) [PMC free article] [PubMed]
92. Barrett R, Schluter D 2008. Adaptation from standing genetic variation. Trends Ecol. Evol. 23, 38–44. (doi:10.1016/j.tree.2007.09.008) [PubMed]
93. Innan H, Kim Y 2004. Pattern of polymorphism after strong artificial selection in a domestication event. Proc. Natl Acad. Sci. USA 101, 10 667–10 672. (doi:10.1073/pnas.0401720101) [PubMed]
94. Przeworski M, Coop G, Wall JD 2005. The signature of positive selection on standing genetic variation. Evolution 59, 2312–2323. (doi:10.1111/j.0014-3820.2005.tb00941.x) [PubMed]

Articles from Proceedings of the Royal Society B: Biological Sciences are provided here courtesy of The Royal Society