Fine-scale empirical data on niche divergence and homeolog expression patterns in an allopolyploid and its diploid progenitor species.

Polyploidization is pervasive in plants, but little is known about the niche divergence of wild allopolyploids (species that harbor polyploid genomes originating from different diploid species) relative to their diploid progenitor species and the gene expression patterns that may underlie such ecological divergence. We conducted a fine-scale empirical study on habitat and gene expression of an allopolyploid and its diploid progenitors. We quantified soil properties and light availability of habitats of an allotetraploid Cardamine flexuosa and its diploid progenitors Cardamine amara and Cardamine hirsuta in two seasons. We analyzed expression patterns of genes and homeologs (homeologous gene copies in allopolyploids) using RNA sequencing. We detected niche divergence between the allopolyploid and its diploid progenitors along water availability gradient at a fine scale: the diploids in opposite extremes and the allopolyploid in a broader range between diploids, with limited overlap with diploids at both ends. Most of the genes whose homeolog expression ratio changed among habitats in C. flexuosa varied spatially and temporally. These findings provide empirical evidence for niche divergence between an allopolyploid and its diploid progenitor species at a fine scale and suggest that divergent expression patterns of homeologs in an allopolyploid may underlie its persistence in diverse habitats.

Positional bias in variant calls against draft reference assemblies.

BACKGROUND: Whole genome resequencing projects may implement variant calling using draft reference genomes assembled de novo from short-read libraries. Despite lower quality of such assemblies, they allowed researchers to extend a wide range of population genetic and genome-wide association analyses to non-model species. As the variant calling pipelines are complex and involve many software packages, it is important to understand inherent biases and limitations at each step of the analysis. RESULTS: In this article, we report a positional bias present in variant calling performed against draft reference assemblies constructed from de Bruijn or string overlap graphs. We assessed how frequently variants appeared at each position counted from ends of a contig or scaffold sequence, and discovered unexpectedly high number of variants at the positions related to the length of either k-mers or reads used for the assembly. We detected the bias in both publicly available draft assemblies from Assemblathon 2 competition as well as in the assemblies we generated from our simulated short-read data. Simulations confirmed that the bias causing variants are predominantly false positives induced by reads from spatially distant repeated sequences. The bias is particularly strong in contig assemblies. Scaffolding does not eliminate the bias but tends to mitigate it because of the changes in variants' relative positions and alterations in read alignments. The bias can be effectively reduced by filtering out the variants that reside in repetitive elements. CONCLUSIONS: Draft genome sequences generated by several popular assemblers appear to be susceptible to the positional bias potentially affecting many resequencing projects in non-model species. The bias is inherent to the assembly algorithms and arises from their particular handling of repeated sequences. It is recommended to reduce the bias by filtering especially if higher-quality genome assembly cannot be achieved. Our findings can help other researchers to improve the quality of their variant data sets and reduce artefactual findings in downstream analyses.

Co-expression network analysis of duplicate genes in maize (Zea mays L.) reveals no subgenome bias.

BACKGROUND: Gene duplication is prevalent in many species and can result in coding and regulatory divergence. Gene duplications can be classified as whole genome duplication (WGD), tandem and inserted (non-syntenic). In maize, WGD resulted in the subgenomes maize1 and maize2, of which maize1 is considered the dominant subgenome. However, the landscape of co-expression network divergence of duplicate genes in maize is still largely uncharacterized. RESULTS: To address the consequence of gene duplication on co-expression network divergence, we developed a gene co-expression network from RNA-seq data derived from 64 different tissues/stages of the maize reference inbred-B73. WGD, tandem and inserted gene duplications exhibited distinct regulatory divergence. Inserted duplicate genes were more likely to be singletons in the co-expression networks, while WGD duplicate genes were likely to be co-expressed with other genes. Tandem duplicate genes were enriched in the co-expression pattern where co-expressed genes were nearly identical for the duplicates in the network. Older gene duplications exhibit more extensive co-expression variation than younger duplications. Overall, non-syntenic genes primarily from inserted duplications show more co-expression divergence. Also, such enlarged co-expression divergence is significantly related to duplication age. Moreover, subgenome dominance was not observed in the co-expression networks - maize1 and maize2 exhibit similar levels of intra subgenome correlations. Intriguingly, the level of inter subgenome co-expression was similar to the level of intra subgenome correlations, and genes from specific subgenomes were not likely to be the enriched in co-expression network modules and the hub genes were not predominantly from any specific subgenomes in maize. CONCLUSIONS: Our work provides a comprehensive analysis of maize co-expression network divergence for three different types of gene duplications and identifies potential relationships between duplication types, duplication ages and co-expression consequences.

Genomic distribution of maize facultative heterochromatin marked by trimethylation of H3K27.

Trimethylation of histone H3 Lys-27 (H3K27me3) plays a critical role in regulating gene expression during plant and animal development. We characterized the genome-wide distribution of H3K27me3 in five developmentally distinct tissues in maize (Zea mays) plants of two genetic backgrounds, B73 and Mo17. There were more substantial differences in the genome-wide profile of H3K27me3 between different tissues than between the two genotypes. The tissue-specific patterns of H3K27me3 were often associated with differences in gene expression among the tissues and most of the imprinted genes that are expressed solely from the paternal allele in endosperm are targets of H3K27me3. A comparison of the H3K27me3 targets in rice (Oryza sativa), maize, and Arabidopsis thaliana provided evidence for conservation of the H3K27me3 targets among plant species. However, there was limited evidence for conserved targeting of H3K27me3 in the two maize subgenomes derived from whole-genome duplication, suggesting the potential for subfunctionalization of chromatin regulation of paralogs. Genomic profiling of H3K27me3 in loss-of-function mutant lines for Maize Enhancer of zeste-like2 (Mez2) and Mez3, two of the three putative H3K27me3 methyltransferases present in the maize genome, suggested partial redundancy of this gene family for maintaining H3K27me3 patterns. Only a portion of the targets of H3K27me3 required Mez2 and/or Mez3, and there was limited evidence for functional consequences of H3K27me3 at these targets.

Epigenetic and genetic influences on DNA methylation variation in maize populations.

DNA methylation is a chromatin modification that is frequently associated with epigenetic regulation in plants and mammals. However, genetic changes such as transposon insertions can also lead to changes in DNA methylation. Genome-wide profiles of DNA methylation for 20 maize (Zea mays) inbred lines were used to discover differentially methylated regions (DMRs). The methylation level for each of these DMRs was also assayed in 31 additional maize or teosinte genotypes, resulting in the discovery of 1966 common DMRs and 1754 rare DMRs. Analysis of recombinant inbred lines provides evidence that the majority of DMRs are heritable. A local association scan found that nearly half of the DMRs with common variation are significantly associated with single nucleotide polymorphisms found within or near the DMR. Many of the DMRs that are significantly associated with local genetic variation are found near transposable elements that may contribute to the variation in DNA methylation. Analysis of gene expression in the same samples used for DNA methylation profiling identified over 300 genes with expression patterns that are significantly associated with DNA methylation variation. Collectively, our results suggest that DNA methylation variation is influenced by genetic and epigenetic changes that are often stably inherited and can influence the expression of nearby genes.

Population genomics of the facultatively mutualistic bacteria Sinorhizobium meliloti and S. medicae.

The symbiosis between rhizobial bacteria and legume plants has served as a model for investigating the genetics of nitrogen fixation and the evolution of facultative mutualism. We used deep sequence coverage (>100×) to characterize genomic diversity at the nucleotide level among 12 Sinorhizobium medicae and 32 S. meliloti strains. Although these species are closely related and share host plants, based on the ratio of shared polymorphisms to fixed differences we found that horizontal gene transfer (HGT) between these species was confined almost exclusively to plasmid genes. Three multi-genic regions that show the strongest evidence of HGT harbor genes directly involved in establishing or maintaining the mutualism with host plants. In both species, nucleotide diversity is 1.5-2.5 times greater on the plasmids than chromosomes. Interestingly, nucleotide diversity in S. meliloti but not S. medicae is highly structured along the chromosome - with mean diversity (θ(π)) on one half of the chromosome five times greater than mean diversity on the other half. Based on the ratio of plasmid to chromosome diversity, this appears to be due to severely reduced diversity on the chromosome half with less diversity, which is consistent with extensive hitchhiking along with a selective sweep. Frequency-spectrum based tests identified 82 genes with a signature of adaptive evolution in one species or another but none of the genes were identified in both species. Based upon available functional information, several genes identified as targets of selection are likely to alter the symbiosis with the host plant, making them attractive targets for further functional characterization.

Reshaping of the maize transcriptome by domestication.

Through domestication, humans have substantially altered the morphology of Zea mays ssp. parviglumis (teosinte) into the currently recognizable maize. This system serves as a model for studying adaptation, genome evolution, and the genetics and evolution of complex traits. To examine how domestication has reshaped the transcriptome of maize seedlings, we used expression profiling of 18,242 genes for 38 diverse maize genotypes and 24 teosinte genotypes. We detected evidence for more than 600 genes having significantly different expression levels in maize compared with teosinte. Moreover, more than 1,100 genes showed significantly altered coexpression profiles, reflective of substantial rewiring of the transcriptome since domestication. The genes with altered expression show a significant enrichment for genes previously identified through population genetic analyses as likely targets of selection during maize domestication and improvement; 46 genes previously identified as putative targets of selection also exhibit altered expression levels and coexpression relationships. We also identified 45 genes with altered, primarily higher, expression in inbred relative to outcrossed teosinte. These genes are enriched for functions related to biotic stress and may reflect responses to the effects of inbreeding. This study not only documents alterations in the maize transcriptome following domestication, identifying several genes that may have contributed to the evolution of maize, but highlights the complementary information that can be gained by combining gene expression with population genetic analyses.

High-density genome-wide association mapping implicates an F-box encoding gene in Medicago truncatula resistance to Aphanomyces euteiches.

• The use of quantitative disease resistance (QDR) is a promising strategy for promoting durable resistance to plant pathogens, but genes involved in QDR are largely unknown. To identify genetic components and accelerate improvement of QDR in legumes to the root pathogen Aphanomyces euteiches, we took advantage of both the recently generated massive genomic data for Medicago truncatula and natural variation of this model legume. • A high-density (≈5.1 million single nucleotide polymorphisms (SNPs)) genome-wide association study (GWAS) was performed with both in vitro and glasshouse phenotyping data collected for 179 lines. • GWAS identified several candidate genes and pinpointed two independent major loci on the top of chromosome 3 that were detected in both phenotyping methods. Candidate SNPs in the most significant locus (σ(A)²= 23%) were in the promoter and coding regions of an F-box protein coding gene. Subsequent qRT-PCR and bioinformatic analyses performed on 20 lines demonstrated that resistance is associated with mutations directly affecting the interaction domain of the F-box protein rather than gene expression. • These results refine the position of previously identified QTL to specific candidate genes, suggest potential molecular mechanisms, and identify new loci explaining QDR against A. euteiches.

Phylogenetic signal variation in the genomes of Medicago (Fabaceae).

Genome-scale data offer the opportunity to clarify phylogenetic relationships that are difficult to resolve with few loci, but they can also identify genomic regions with evolutionary history distinct from that of the species history. We collected whole-genome sequence data from 29 taxa in the legume genus Medicago, then aligned these sequences to the Medicago truncatula reference genome to confidently identify 87 596 variable homologous sites. We used this data set to estimate phylogenetic relationships among Medicago species, to investigate the number of sites needed to provide robust phylogenetic estimates and to identify specific genomic regions supporting topologies in conflict with the genome-wide phylogeny. Our full genomic data set resolves relationships within the genus that were previously intractable. Subsampling the data reveals considerable variation in phylogenetic signal and power in smaller subsets of the data. Even when sampling 5000 sites, no random sample of the data supports a topology identical to that of the genome-wide phylogeny. Phylogenetic relationships estimated from 500-site sliding windows revealed genome regions supporting several alternative species relationships among recently diverged taxa, consistent with the expected effects of deep coalescence or introgression in the recent history of Medicago.

Whole-genome nucleotide diversity, recombination, and linkage disequilibrium in the model legume Medicago truncatula.

Medicago truncatula is a model for investigating legume genetics, including the genetics and evolution of legume-rhizobia symbiosis. We used whole-genome sequence data to identify and characterize sequence polymorphisms and linkage disequilibrium (LD) in a diverse collection of 26 M. truncatula accessions. Our analyses reveal that M. truncatula harbors both higher diversity and less LD than soybean (Glycine max) and exhibits patterns of LD and recombination similar to Arabidopsis thaliana. The population-scaled recombination rate is approximately one-third of the mutation rate, consistent with expectations for a species with a high selfing rate. Linkage disequilibrium, however, is not extensive, and therefore, the low recombination rate is likely not a major constraint to adaptation. Nucleotide diversity in 100-kb windows was negatively correlated with gene density, which is expected if diversity is shaped by selection acting against slightly deleterious mutations. Among putative coding regions, members of four gene families harbor significantly higher diversity than the genome-wide average. Three of these families are involved in resistance against pathogens; one of these families, the nodule-specific, cysteine-rich gene family, is specific to the galegoid legumes and is involved in control of rhizobial differentiation. The more than 3 million SNPs that we detected, approximately one-half of which are present in more than one accession, are a valuable resource for genome-wide association mapping of genes responsible for phenotypic diversity in legumes, especially traits associated with symbiosis and nodulation.

Selection, genome-wide fitness effects and evolutionary rates in the model legume Medicago truncatula.

Sequence data for >20 000 annotated genes from 56 accessions of Medicago truncatula were used to identify potential targets of positive selection, the determinants of evolutionary rate variation and the relative importance of positive and purifying selection in shaping nucleotide diversity. Based upon patterns of intraspecific diversity and interspecific divergence, c. 50-75% of nonsynonymous polymorphisms are subject to strong purifying selection and 1% of the sampled genes harbour a signature of positive selection. Combining polymorphism with expression data, we estimated the distribution of fitness effects and found that the proportion of deleterious mutations is significantly greater for expressed genes than for genes with undetected transcripts (nonexpressed) in a previous RNA-seq experiment and greater for broadly expressed genes than those expressed in only a single tissue. Expression level is the strongest correlate of evolutionary rates at nonsynonymous sites, and despite multiple genomic features being significantly correlated with evolutionary rates, they explain less than 20% of the variation in nonsynonymous rates (dN) and <15% of the variation in either synonymous rates (dS) or dN:dS. Among putative targets of selection were genes involved in defence against pathogens and herbivores, genes with roles in mediating the relationship with rhizobial symbionts and one-third of annotated histone-lysine methyltransferases. Adaptive evolution of the methyltransferases suggests that positive selection in gene expression may have occurred through evolution of enzymes involved in epigenetic modification.

Patterns of polymorphism and selection in the subgenomes of the allopolyploid Arabidopsis kamchatica.

Genome duplication is widespread in wild and crop plants. However, little is known about genome-wide selection in polyploids due to the complexity of duplicated genomes. In polyploids, the patterns of purifying selection and adaptive substitutions may be affected by masking owing to duplicated genes or homeologs as well as effective population size. Here, we resequence 25 accessions of the allotetraploid Arabidopsis kamchatica, which is derived from the diploid species A. halleri and A. lyrata. We observe a reduction in purifying selection compared with the parental species. Interestingly, proportions of adaptive non-synonymous substitutions are significantly positive in contrast to most plant species. A recurrent pattern observed in both frequency and divergence-diversity neutrality tests is that the genome-wide distributions of both subgenomes are similar, but the correlation between homeologous pairs is low. This may increase the opportunity of different evolutionary trajectories such as in the HMA4 gene involved in heavy metal hyperaccumulation.

Genome assembly and annotation of Arabidopsis halleri, a model for heavy metal hyperaccumulation and evolutionary ecology.

The self-incompatible species Arabidopsis halleri is a close relative of the self-compatible model plant Arabidopsis thaliana. The broad European and Asian distribution and heavy metal hyperaccumulation ability make A. halleri a useful model for ecological genomics studies. We used long-insert mate-pair libraries to improve the genome assembly of the A. halleri ssp. gemmifera Tada mine genotype (W302) collected from a site with high contamination by heavy metals in Japan. After five rounds of forced selfing, heterozygosity was reduced to 0.04%, which facilitated subsequent genome assembly. Our assembly now covers 196 Mb or 78% of the estimated genome size and achieved scaffold N50 length of 712 kb. To validate assembly and annotation, we used synteny of A. halleri Tada mine with a previously published high-quality reference assembly of a closely related species, Arabidopsis lyrata. Further validation of the assembly quality comes from synteny and phylogenetic analysis of the HEAVY METAL ATPASE4 (HMA4) and METAL TOLERANCE PROTEIN1 (MTP1) regions using published sequences from European A. halleri for comparison. Three tandemly duplicated copies of HMA4, key gene involved in cadmium and zinc hyperaccumulation, were assembled on a single scaffold. The assembly will enhance the genomewide studies of A. halleri as well as the allopolyploid Arabidopsis kamchatica derived from A. lyrata and A. halleri.

The Cardamine hirsuta genome offers insight into the evolution of morphological diversity.

Finding causal relationships between genotypic and phenotypic variation is a key focus of evolutionary biology, human genetics and plant breeding. To identify genome-wide patterns underlying trait diversity, we assembled a high-quality reference genome of Cardamine hirsuta, a close relative of the model plant Arabidopsis thaliana. We combined comparative genome and transcriptome analyses with the experimental tools available in C. hirsuta to investigate gene function and phenotypic diversification. Our findings highlight the prevalent role of transcription factors and tandem gene duplications in morphological evolution. We identified a specific role for the transcriptional regulators PLETHORA5/7 in shaping leaf diversity and link tandem gene duplication with differential gene expression in the explosive seed pod of C. hirsuta. Our work highlights the value of comparative approaches in genetically tractable species to understand the genetic basis for evolutionary change.

Discovering functional modules across diverse maize transcriptomes using COB, the Co-expression Browser.

Tools that provide improved ability to relate genotype to phenotype have the potential to accelerate breeding for desired traits and to improve our understanding of the molecular variants that underlie phenotypes. The availability of large-scale gene expression profiles in maize provides an opportunity to advance our understanding of complex traits in this agronomically important species. We built co-expression networks based on genome-wide expression data from a variety of maize accessions as well as an atlas of different tissues and developmental stages. We demonstrate that these networks reveal clusters of genes that are enriched for known biological function and contain extensive structure which has yet to be characterized. Furthermore, we found that co-expression networks derived from developmental or tissue atlases as compared to expression variation across diverse accessions capture unique functions. To provide convenient access to these networks, we developed a public, web-based Co-expression Browser (COB), which enables interactive queries of the genome-wide networks. We illustrate the utility of this system through two specific use cases: one in which gene-centric queries are used to provide functional context for previously characterized metabolic pathways, and a second where lists of genes produced by mapping studies are further resolved and validated using co-expression networks.

Genomic signature of adaptation to climate in Medicago truncatula.

Local adaptation and adaptive clines are pervasive in natural plant populations, yet the effects of these types of adaptation on genomic diversity are not well understood. With a data set of 202 accessions of Medicago truncatula genotyped at almost 2 million single nucleotide polymorphisms, we used mixed linear models to identify candidate loci responsible for adaptation to three climatic gradients-annual mean temperature (AMT), precipitation in the wettest month (PWM), and isothermality (ITH)-representing the major axes of climate variation across the species' range. Loci with the strongest association to these climate gradients tagged genome regions with high sequence similarity to genes with functional roles in thermal tolerance, drought tolerance, or resistance to herbivores of pathogens. Genotypes at these candidate loci also predicted the performance of an independent sample of plant accessions grown in climate-controlled conditions. Compared to a genome-wide sample of randomly drawn reference SNPs, candidates for two climate gradients, AMT and PWM, were significantly enriched for genic regions, and genome segments flanking genic AMT and PWM candidates harbored less nucleotide diversity, elevated differentiation between haplotypes carrying alternate alleles, and an overrepresentation of the most common haplotypes. These patterns of diversity are consistent with a history of soft selective sweeps acting on loci underlying adaptation to climate, but not with a history of long-term balancing selection.

Maize gene atlas developed by RNA sequencing and comparative evaluation of transcriptomes based on RNA sequencing and microarrays.

Transcriptome analysis is a valuable tool for identification and characterization of genes and pathways underlying plant growth and development. We previously published a microarray-based maize gene atlas from the analysis of 60 unique spatially and temporally separated tissues from 11 maize organs [1]. To enhance the coverage and resolution of the maize gene atlas, we have analyzed 18 selected tissues representing five organs using RNA sequencing (RNA-Seq). For a direct comparison of the two methodologies, the same RNA samples originally used for our microarray-based atlas were evaluated using RNA-Seq. Both technologies produced similar transcriptome profiles as evident from high Pearson's correlation statistics ranging from 0.70 to 0.83, and from nearly identical clustering of the tissues. RNA-Seq provided enhanced coverage of the transcriptome, with 82.1% of the filtered maize genes detected as expressed in at least one tissue by RNA-Seq compared to only 56.5% detected by microarrays. Further, from the set of 465 maize genes that have been historically well characterized by mutant analysis, 427 show significant expression in at least one tissue by RNA-Seq compared to 390 by microarray analysis. RNA-Seq provided higher resolution for identifying tissue-specific expression as well as for distinguishing the expression profiles of closely related paralogs as compared to microarray-derived profiles. Co-expression analysis derived from the microarray and RNA-Seq data revealed that broadly similar networks result from both platforms, and that co-expression estimates are stable even when constructed from mixed data including both RNA-Seq and microarray expression data. The RNA-Seq information provides a useful complement to the microarray-based maize gene atlas and helps to further understand the dynamics of transcription during maize development.

Candidate genes and genetic architecture of symbiotic and agronomic traits revealed by whole-genome, sequence-based association genetics in Medicago truncatula.

Genome-wide association study (GWAS) has revolutionized the search for the genetic basis of complex traits. To date, GWAS have generally relied on relatively sparse sampling of nucleotide diversity, which is likely to bias results by preferentially sampling high-frequency SNPs not in complete linkage disequilibrium (LD) with causative SNPs. To avoid these limitations we conducted GWAS with >6 million SNPs identified by sequencing the genomes of 226 accessions of the model legume Medicago truncatula. We used these data to identify candidate genes and the genetic architecture underlying phenotypic variation in plant height, trichome density, flowering time, and nodulation. The characteristics of candidate SNPs differed among traits, with candidates for flowering time and trichome density in distinct clusters of high linkage disequilibrium (LD) and the minor allele frequencies (MAF) of candidates underlying variation in flowering time and height significantly greater than MAF of candidates underlying variation in other traits. Candidate SNPs tagged several characterized genes including nodulation related genes SERK2, MtnodGRP3, MtMMPL1, NFP, CaML3, MtnodGRP3A and flowering time gene MtFD as well as uncharacterized genes that become candidates for further molecular characterization. By comparing sequence-based candidates to candidates identified by in silico 250K SNP arrays, we provide an empirical example of how reliance on even high-density reduced representation genomic makers can bias GWAS results. Depending on the trait, only 30-70% of the top 20 in silico array candidates were within 1 kb of sequence-based candidates. Moreover, the sequence-based candidates tagged by array candidates were heavily biased towards common variants; these comparisons underscore the need for caution when interpreting results from GWAS conducted with sparsely covered genomes.

Fine-scale population recombination rates, hotspots, and correlates of recombination in the Medicago truncatula genome.

Recombination rates vary across the genome and in many species show significant relationships with several genomic features, including distance to the centromere, gene density, and GC content. Studies of fine-scale recombination rates have also revealed that in several species, there are recombination hotspots, that is, short regions with recombination rates 10-100 greater than those in surrounding regions. In this study, we analyzed whole-genome resequence data from 26 accessions of the model legume Medicago truncatula to gain insight into the genomic features that are related to high- and low-recombination rates and recombination hotspots at 1 kb scales. We found that high-recombination regions (1-kb windows among those in the highest 5% of the distribution) on all three chromosomes were significantly closer to the centromere, had higher gene density, and lower GC content than low-recombination windows. High-recombination windows are also significantly overrepresented among some gene functional categories-most strongly NB-ARC and LRR genes, both of which are important in plant defense against pathogens. Similar to high-recombination windows, recombination hotspots (1-kb windows with significantly higher recombination than the surrounding region) are significantly nearer to the centromere than nonhotspot windows. By contrast, we detected no difference in gene density or GC content between hotspot and nonhotspot windows. Using linear model wavelet analysis to examine the relationship between recombination and genomic features across multiple spatial scales, we find a significant negative correlation with distance to the centromere across scales up to 512 kb, whereas gene density and GC content show significantly positive and negative correlations, respectively, only up to 64 kb. Correlations between recombination and genomic features, particularly gene density and polymorphism, suggest that they are scale dependent and need to be assessed at scales relevant to the evolution of those features.