- Research article
- Open Access
Adaptive genomic structural variation in the grape powdery mildew pathogen, Erysiphe necator
BMC Genomicsvolume 15, Article number: 1081 (2014)
Powdery mildew, caused by the obligate biotrophic fungus Erysiphe necator, is an economically important disease of grapevines worldwide. Large quantities of fungicides are used for its control, accelerating the incidence of fungicide-resistance. Copy number variations (CNVs) are unbalanced changes in the structure of the genome that have been associated with complex traits. In addition to providing the first description of the large and highly repetitive genome of E. necator, this study describes the impact of genomic structural variation on fungicide resistance in Erysiphe necator.
A shotgun approach was applied to sequence and assemble the genome of five E. necator isolates, and RNA-seq and comparative genomics were used to predict and annotate protein-coding genes. Our results show that the E. necator genome is exceptionally large and repetitive and suggest that transposable elements are responsible for genome expansion. Frequent structural variations were found between isolates and included copy number variation in EnCYP51, the target of the commonly used sterol demethylase inhibitor (DMI) fungicides. A panel of 89 additional E. necator isolates collected from diverse vineyard sites was screened for copy number variation in the EnCYP51 gene and for presence/absence of a point mutation (Y136F) known to result in higher fungicide tolerance. We show that an increase in EnCYP51 copy number is significantly more likely to be detected in isolates collected from fungicide-treated vineyards. Increased EnCYP51 copy numbers were detected with the Y136F allele, suggesting that an increase in copy number becomes advantageous only after the fungicide-tolerant allele is acquired. We also show that EnCYP51 copy number influences expression in a gene-dose dependent manner and correlates with fungal growth in the presence of a DMI fungicide.
Taken together our results show that CNV can be adaptive in the development of resistance to fungicides by providing increasing quantitative protection in a gene-dosage dependent manner. The results of this work not only demonstrate the effectiveness of using genomics to dissect complex traits in organisms with very limited molecular information, but also may have broader implications for understanding genomic dynamics in response to strong selective pressure in other pathogens with similar genome architectures.
Grapevine powdery mildew is one of the most widespread and devastating diseases of wine, table and raisin grapes, the vast majority of which are cultivars of Vitis vinifera. This disease is caused by the fungus Erysiphe necator Schw. [syn. Uncinula necator (Schw.) Burr.], an obligate biotroph that can infect all green tissues of a grapevine (Figure 1A-C ). Infected leaves exhibit reduced photosynthesis and often undergo premature senescence and abscission. Early berry infection causes berries to crack, and the overall impact on the crop includes decreased yields, increased acidity, and decreased anthocyanin and sugar content of mature fruit . Even low levels of powdery mildew infection on the berries can lead to ruined table grapes and wines with negative sensory attributes and decreased varietal character [2, 3].
Most cultivated varieties are susceptible to powdery mildew, and as a consequence growers are forced to apply fungicides to control the disease, often as frequently as every 7–10 days when disease pressure is high. It has been estimated that as much as 20% of the total costs associated with wine grape production in California goes to expenses related to powdery mildew control . Elemental sulfur was the first effective fungicide recommended for vineyards in 1848 to control powdery mildew, and it continues to be widely used, mainly due to its efficacy and low cost . Although its multi-site mode of action remains effective at controlling powdery mildew, the limitations to sulfur use include phytotoxicity at higher temperatures, the need for application as a protectant at frequent intervals, potential off-characters in wine, and the risk of unintended environmental consequences [9, 10]. Since the 1960s, new classes of fungicides have been developed and introduced, many with single-site modes of action and beneficial properties such as systemic effectiveness and longer times between applications . One important class of single-site fungicides is the sterol demethylase inhibitors (DMI), which includes the azole fungicides. DMIs inhibit fungal growth by targeting the cytochrome P450 lanosterol C-14α-demethylase (CYP51, also known as ERG11). CYP51 is a key enzyme in eukaryotic sterol biosynthesis, catalyzing the removal of the 14-α carbon of lanosterol during biosynthesis of ergosterol, the major fungal sterol and an important membrane component . Ergosterol depletion and accumulation of deleterious methylated sterols caused by CYP51 inhibition affect membrane integrity and functionality, reducing fungal growth .
DMIs are prone to resistance issues due to their single metabolic target. Triadimefon was the first DMI fungicide to be approved for powdery mildew control in vineyards in 1982, and resistance was documented as soon as 1986 [9, 13]. Resistance to other DMIs has since been reported, including cross-resistance between different compounds . DMI resistance in fungi has been associated with several molecular mechanisms, such as CYP51 mutations  and increased expression of CYP51. A single point mutation in CYP51 causing a tyrosine for phenylalanine substitution at position 136 (Y136F) is associated with DMI resistance in E. necator, and similar mutations conferring DMI resistance have been identified in other fungi including Blumeria graminis f. sp. hordei (Bgh, henceforth ), Candida albicans, and Mycosphaerella graminicola. The amino acid 136 is within the highly conserved CR2 domain involved in substrate recognition . Structural modeling has demonstrated that DMIs bind CYP51 at the sixth ligand of the heme group, and it has been suggested that the Y136F mutation confers resistance by obstructing the heme-binding site .
Pathogen populations, particularly species with mixed sexual and asexual reproduction systems like E. necator, rapidly evolve to circumvent host resistance genes or to counteract chemical control methods . Though new fungicides are being developed and breeding for genetic resistance is ongoing, little is known regarding how E. necator evolves to become resistant to fungicides or how it overcomes host immune responses. Whole-genome sequencing, transcriptomics, and comparative genomic approaches have been successfully applied to other fungal plant pathogens . However, publicly available genomic information for E. necator is scarce and limited to SSR markers [24–26], likely reflecting the difficulty of working with an obligate parasite. Genomic analyses of related powdery mildew species (Figure 1D), including the ascomycetes that cause powdery mildew in barley (Bgh; [27, 28]) and wheat (B. graminis f.sp. tritici; Bgt, henceforth; ) have provided insights into their evolution and mechanisms of pathogenicity.
In this study we applied a shotgun approach to sequence the genomes of five isolates of E. necator, and used RNA-seq and comparative genomics to predict and annotate protein-coding genes. Our results not only show that the E. necator genome is exceptionally large and repetitive compared to other fungal plant pathogens (as previously seen in other powdery mildews), but also that frequent structural variations occur between field isolates. Structural variation resulted in a wide range of copy numbers in the DMI fungicide target gene CYP51 (EnCYP51). A panel of 89 additional E. necator isolates collected from both synthetic fungicide-treated and non-treated vineyard sites was screened for both copy number variation and presence/absence of the Y136F mutation in the EnCYP51 gene. The vast majority (95%) of isolates with increased copy number were collected from fungicide-treated vineyards. Increased EnCYP51 copy numbers were almost always (94%) detected with the EnCYP51 mutant allele (Y136F), suggesting that an increase in copy number becomes advantageous only when the EnCYP51 gene is present in its fungicide-tolerant allelic form. We also show that EnCYP51 genomic copy number is correlated with both its mRNA transcript level and with fungal growth in the presence of DMI fungicide. Taken together, our results demonstrate that copy number variation can be adaptive in the development of resistance to DMI fungicides in E. necator.
Genome sequencing and assembly of five E. necatorisolates
Between 2012 and 2013, 94 E. necator isolates were collected from seven vineyard sites in California (Additional file 1: Table S1). Five isolates with diverse genetic and geographical backgrounds (Branching, C-strain, e1-101, Lodi, and Ranch9) were selected for genome sequencing based on microsatellite profiling results, collection site, and fungicide treatment history (Additional file 1: Table S1; Additional file 2: Figure S1; Additional file 3: Text S1). The Branching and e1-101 isolates were collected from vineyard sites without fungicide treatment, while the C-strain, Lodi, and Ranch9 isolates were all collected from sites regularly treated with synthetic fungicides to control powdery mildew.
We used an Illumina shotgun sequencing approach to sequence the haploid genomes of the five selected E. necator isolates. We generated both medium reads (∼250 bases) and short reads (∼150 bases) using the Illumina MiSeq and HiSeq2500 sequencers, respectively. The longer reads were used to assemble the draft genome sequence for the C-strain isolate, while the shorter reads were used for resequencing and de novo assembly of the other isolates. After quality trimming and filtering, most sequencing reads (98.5 ± 0.2%) were assembled into 6,475 ± 444 scaffolds with a total length of 50.5 ± 1.2 Mb (Table 1; Additional file 4: Table S2; mitochondrial DNA scaffolds: 188 kb ± 1.3 kb). Most of the reads of each isolate (97.4 ± 0.3%) aligned as pair-ends to the assembled contigs of the other isolates suggesting a similar genome completeness and representation across all five isolates. Unambiguously mapped reads were used to assess single nucleotide polymorphism (SNP) frequencies between the isolates following the GATK pipeline for SNP calling. The GATK analysis identified an average SNP frequency of 0.79 ± 0.08 SNPs/kb (Additional file 5: Table S3).
Based on the distribution of DNA k-mer counts , we estimated the genome size of E. necator at 126 ± 18 Mb, which is similar to the 120 Mb estimated genome size of Bgh[27, 28] and significantly smaller than the 180 Mb genome size estimated for Bgt. The discrepancy between assembled scaffold size and estimated genome size suggests that a large proportion of the E. necator genome is highly repetitive and collapsed into common contigs during assembly. This hypothesis is also supported by the high sequencing coverage observed in the highly fragmented and gene-poor fraction of the assemblies (Additional file 6: Figure S2A). In contrast, the longest and most gene-rich contigs exhibit the expected sequencing coverage based on estimated genome size and sequencing depth (Additional file 6: Figure S2A).
Depth of coverage was used to estimate the size of the interspersed repeats. Repeats, including transposable elements (TEs), microsatellites, and low-complexity regions, were estimated to represent 62.9 ± 3% of the E. necator genome (Additional file 6: Figure S2B), a proportion similar to the repeat content estimated for the Bgh genome (64%; ). TEs were annotated using a combination of ab initio discovery and comparative analysis. TEs accounted for 77.9 ± 3.6% of the total repeats, with LTR retrotransposons and LINEs being the most abundant classes of TEs in the genome (Additional file 6: Figure S2B). TEs were abundant among the sequenced transcripts, suggesting that in E. necator TEs are transcriptionally active (Additional file 6: Figure S2D).
To assess the quality and completeness of the assembled gene space, we applied the CEGMA pipeline  to map a set of 248 low copy Core Eukaryotic Orthologous groups (KOGs) which are conserved across higher eukaryotes. On average 96.2 ± 0.5% of the KOGs aligned as complete gene copies to the scaffolds, compared to 98.1 ± 0.4% that aligned as fragmented partial gene copies (Table 1; Additional file 4: Table S2). These results indicate that our shotgun sequencing approach generated a nearly complete assembly of the E. necator gene space, despite the obvious failure in resolving the structure of the repetitive fraction of the genome. For the analyses described below we used the genome assembly of the C-strain isolate as a reference, because it appeared to be the most complete based on assembly metrics and CEGMA evaluation (Table 1; Additional file 4: Table S2).
Gene prediction and annotation
To enhance gene discovery and gene structure prediction, we sequenced and assembled the E. necator transcripts expressed during the interaction with grape leaves at 0.5, 1, 3, and 6 days after conidia inoculation (Additional file 7: Table S4; Additional file 8: Figure S3A&B). Assembled transcripts were then used as evidence in a gene prediction pipeline that combined evidence-based, homology-based, and ab initio methods (Additional file 8: Figure S3C). After removing 4,511 genes associated with TEs, we estimated that the E. necator genome contains 6,533 protein-coding genes. A similar total number of protein-coding genes was identified in Bgh (6,470 ; Additional file 9: Table S5) and Bgt (6,525 ). Nearly 92% of the predicted genes have detectable similarity to proteins known from other ascomycetes (BLASTP e-value < 10−3; NCBI Ascomycete nr; Additional file 9: Table S5; Additional file 10: Table S6; Additional file 11: Figure S4). Percent BLASTP matches and alignment coverage between the peptides of E. necator and other fungi, support the overall quality of the gene prediction pipeline results (Additional file 12: Table S7). Prediction accuracy was also supported by the identification of higher frequencies in nucleotide variants in the intergenic space (0.89 ± 0.01 SNPs/kb) and intronic regions (0.57 ± 0.07 SNPs/kb) than in the exons (0.35 ± 0.04 SNPs/kb; Additional file 5: Table S3). An average of 1.7 introns per gene were predicted and, similarly to observations in other fungal genomes , intron lengths had a narrow distribution with mean size of 83.7 bp (exon mean size = 526.5 bp; Additional file 8: Figure S3D-F). Genes on longer scaffolds displayed shorter intergenic spaces than those on smaller scaffolds, suggestive of an uneven distribution of genes within the genome (Additional file 8: Figure S3G). As reported in previous analyses of powdery mildews [27, 28] and other obligate biotrophs [33–35], the E. necator genome showed a loss of enzymes involved in secondary metabolism, nitrate and sulfate metabolism, further supporting the convergent adaptation of obligate biotrophy (Additional file 3: Text S1 ). We did not find any of the genes known to be required for Repeat Induced Point mutation (RIP) suggesting that the mechanisms that usually control TEs in fungi are only partially functional in E. necator, as previously observed in other powdery mildew pathogens (Additional file 3: Text S1; ).
Secretome and candidate effector proteins
Microbial plant pathogens secrete effector proteins, some of which have been shown to subvert the plant innate immune response and enable infection . The predicted secretome of E. necator contained a total of 607 genes. GO term enrichment analysis of the secretome showed enrichment in hydrolase and peptidase molecular function (MF) terms (Additional file 13: Table S8), although many of these enzymes possess β-1,3-glucanase or chitinase activity and are likely involved in remodeling the pathogen’s own cell wall during early infection instead of breaking down the host tissue as in necrotrophy .
The E. necator secretome also included several classes of proteins that were previously identified as pathogen candidate effectors, such as EKA-like proteins  and ribonucleases-like proteins [39, 40]. These putative effector proteins were detected in transcriptome data from early infection time points suggesting that some might be involved in the very early stages of infection (Additional file 10: Table S6; Additional file 14: Table S9; Additional file 6: Figure S2D). Following the same approach described in , we identified 150 candidate secreted effector proteins (CSEPs; Additional file 3: Text S1), which displayed a significant enrichment in sequence motifs associated with secreted peptides of haustoria-forming pathogenic fungi (; Additional file 14: Table S9). These candidate effectors did not present any signature of positive selection unlike those described in other plant pathogens (; Additional file 3: Text S1). The absence of elevated levels of non-synonymous polymorphisms among the candidate effectors may hypothetically reflect the absence of antagonistic evolution in the host, since almost all cultivated grape varieties are fully susceptible to powdery mildew.
Structural variation between E. necatorisolates
Copy number variations (CNVs) are unbalanced changes in the structure of the genome and include deletions, insertions, and duplications of >1 kb in size . Although short sequencing read length and limited insert size do not allow reliable characterization of the type of structural variation underlying CNV, variation in sequencing coverage in genome assemblies can be used as an indicator of CNV between an assembled genome and sequencing reads from another genome (; Figure 2A&B). To detect structural variants between the five sequenced isolates, we applied the CNV-seq pipeline to the mapped reads . CNV-seq calculates the read-depth signal across all genomic coordinates using a sliding window and was shown to provide more accurate relative copy number estimation than other CNV-detection tools . Genome-wide CNV-seq results are summarized in Table 2. While the sizes of the detected CNVs were similar across the isolates, the total number of CNV sites had a broad range from 177 to 1,170, corresponding to 1.1% and 5.3% of the assemblies, respectively. The identified CNV loci not only spanned TEs (LTR: 364 ± 32 kb; non-LTR TEs: 99 ± 35 kb; Additional file 6: Figure S2C), but also 135 protein-coding genes (49.3 ± 4.0 CNV genes/genome; Additional file 10: Table S6). A broad diversity of genes was observed within CNV loci with no detectable enrichment for any particular function. CNV-seq estimations of copy number polymorphisms were validated by quantitative PCR (qPCR; Figure 2C).
Copy number variation of the EnCYP51locus in the five sequenced isolates
The EnCYP51 gene, which encodes the enzyme targeted by DMI fungicides, was among the protein-coding genes with a wide range of copy number variation across the five sequenced genomes. Sequencing coverage, confirmed by qPCR, estimated the presence of 3, 4, and 9 EnCYP51 copies in Ranch 9, C-strain, and Lodi, respectively, and single copies in e1-101 and Branching (Figure 3). Remarkably, the three isolates collected from fungicide-treated vineyards showed multiple copies of EnCYP51 with all copies being the Y136F mutated allele, while the two isolates from vineyards not treated with fungicides both had single copies of the EnCYP51 gene in its fungicide-susceptible allelic form (Additional file 15: Figure S5A&B). The Y136F substitution was the only non-synonymous polymorphism detected between the isolates in the EnCYP51 coding region (Additional file 15: Figure S5C). All EnCYP51 alleles were confirmed by Sanger sequencing, which also confirmed that all EnCYP51 duplicated copies were in the same allelic form (Additional file 15: Figure S5B).
Duplication boundaries, approximately 8.0 and 0.6 kb upstream and downstream of the EnCYP51 coding region, respectively, were common in all isolates (Figure 3). Only in Lodi did an additional duplication event occur involving a shorter region of 7.4 kb (Figure 3). In all isolates, the duplication events did not involve the two flanking genes, En-g4918 and En-g3071, which were confirmed by qPCR to be single copy in all genomes. In all isolates the downstream duplication boundary coincided with a fragment of a Gypsy element (Figure 3). Although numerous other TEs were found throughout the sequence flanking the duplicated region, we did not detect repeated homologous sequences that could provide the substrate for non-allelic homologous recombination and allow us to hypothesize the mechanism underlying the copy number polymorphisms in the locus . Microsynteny of the EnCYP51 locus between E. necator, Bgh, and Bgt (Figure 4) both confirms the accuracy of the sequence assembly and also shows the structural conservation in this region of the genome across species despite the frequent structural rearrangements observed in the EnCYP51 locus.
CNVs can modify the expression of genes that localize within the CNV in a gene dose-dependent manner . To determine the impact of CNV on EnCYP51 expression, we measured mRNA transcript expression levels at 1 and 5 days post conidia inoculation using qRT-PCR. Our results show that EnCYP51 was differentially expressed in the different isolates at both time points (Figure 5A) with a strong correlation between gene expression levels and EnCYP51 copy number (Figure 5B).
Copy number variation of the EnCYP51 locus occurs frequently in E. necatorand is correlated with fungal growth under fungicide treatment
To determine the extent of EnCYP51 CNV in E. necator beyond the five sequenced isolates, we screened an additional 89 isolates collected from multiple vineyard locations (Additional file 1: Table S1). All isolates tested for CNV were also genotyped to determine the allelic forms of EnCYP51. Out of the 89 total isolates, 48% carried the wild-type allele, while 52% carried the Y136F mutation (Figure 6). The majority (88%) of the isolates with the wild type allele were collected from non-sprayed vineyards, while almost all (95%) of the isolates with the Y136F mutation came from sprayed vineyards (Figure 6B). Quantitative PCR analysis showed the presence of multiple copies of EnCYP51 in 47 isolates with a wide range of copy numbers, from as few as 2 copies and up to 14 copies (Figure 6C). The two genes that flank EnCYP51 were mostly single copy confirming the conservation of the duplication boundaries (Figure 6D). Remarkably, 94% of the isolates with multiple copies of EnCYP51 also carried the mutant Y136F allele, while 96% of the isolates with single copies of EnCYP51 were wild type, further supporting the hypothesis that the acquisition of multiple copies of the fungicide-tolerant allele contributes to the evolutionary fitness. Hypergeometric tests supported the significant association between fungicide treatment and both presence of the Y136F allele (P = 4.8 × 10−4) and multiple copies of EnCYP51 (P = 4.3 × 10−5). The EnCYP51 CNV and its allelic form did not correlate with SSR-based clustering (Figure 6A; Additional file 3: Text S1) suggesting that the Y136F mutation and the acquisition of multiple of copies of the EnCYP51 gene occurred multiple times in the E. necator populations screened.
To test the hypothesis that an increase in EnCYP51 copy number leads to stronger resistance to DMI fungicides, we tested the capability of isolates with different copy numbers of EnCYP51 to grow in the presence of Rally® 40WSP, a DMI fungicide. Eight isolates with a range of CNVs in EnCYP51 were tested, including the five isolates we sequenced. A qPCR analysis was performed to quantify the accumulation of E. necator biomass eight days after conidia inoculation. As expected, under fungicide treatment we did not detect any growth of Branching and e1-101, which are the two isolates that carry the wild type and single copy EnCYP51, while a positive correlation between fungal growth in presence of the DMI fungicide and the estimated number of EnCYP51 copies was observed in the isolates containing the Y136F mutation (Figure 6E; Additional file 16: Figure S6B). While we know the number of EnCYP51 copies and the EnCYP51 allelic form in each isolate, we were unable to control for genetic differences in other genes in the genome that may mask or confound the role of EnCYP51. For example it has been suggested that ABC transporters contribute to DMI resistance , and although we did not see sequence polymorphisms and CNVs that may suggest a contribution of ABC transporters to the differential DMI susceptibility of the isolates we tested, we could not exclude the possibility that other loci in the genome are responsible for the differential growth in response to the fungicide treatment.
This study describes for the first time the genome of one of the most widespread and devastating diseases of grapevines. With an estimated size of ~125 Mb, the E. necator genome is among the largest sequenced ascomycete genomes. The E. necator genome also possessed a characteristic reduction in protein-coding genes, a feature also seen in the genomes of other biotrophic fungi . The expansion of genome size in powdery mildews [27–29], including E. necator, is associated with the proliferation of repetitive DNA, particularly TEs, which account for about 63% of the E. necator genome. TEs were also detected in the E. necator transcriptome, indicating that they are transcriptionally active at least during host infection. Active TEs are a major source of mutations in the genome. Transposition of TEs not only contributes to the expansion of genome size , but can also cause chromosome breaks and rearrangements , gene deletions , gene duplications , illegitimate recombination , and changes in gene expression [56, 57]. TE activity has also been associated with the evolution of pathogen virulence factors .
Frequent structural variations in the genome also point to a highly unstable genomic architecture in E. necator, although a direct link between TE proliferation and CNV cannot yet be established based on our results. In this study, the analysis of read depth showed that up to 5% of the genome assemblies, corresponding to approximately 2.8 Mb, are CNV, which indicates that at least among the five isolates we sequenced CNV involves larger genomic regions than DNA sequence polymorphisms (SNP: 41,657 ± 4,151 SNPs/genome). Studies with model organisms, including fungi, plants, and animals, have already shown that the extent of genetic differences between individual organisms due to alterations in the number of gene copies can often be greater than differences in the nucleotide sequences . CNVs can result from different types of structural variations, such as deletions, translocations, inversions, tandem duplications, and novel insertions. Because read depth analysis cannot identify with certainty the type of rearrangement that underlies the CNV, characterization of the type of medium- and large-scale structural genomic variations would require longer reads, larger inserts between paired-end reads, and longer scaffolds than those used in this study. Nonetheless, qPCR validation of CNVs (Figure 2C & Figure 3 inset) confirmed that read depth analysis can be used to reliably detect CNV regions in the genome.
CNVs can influence phenotypic variation and adaptation by modification of gene dosage, gene expression, or disruption of genes that span boundaries of structural rearrangements. CNVs have been associated with phenotypic differences and complex diseases in humans , adaptation to different soil types in Arabidopsis lyrata, higher aluminum tolerance in maize , nematode resistance in soybean , drug resistance in Plasmodium falciparum, and DDT resistance in Drosophila melanogaster. A role of CNV in contributing to effector protein evolution was described for the oomycete plant pathogen, Phytophtora sojae.
Here we show that (i) the DMI fungicide target EnCYP51 is CNV in 95% of isolates collected from vineyards treated with fungicides, (ii) increasing EnCYP51 copy number leads to higher expression of EnCYP51, and (iii) increasing EnCYP51 copy number correlates with increased ability to withstand fungicide treatment. Furthermore, Sanger sequencing confirmed that the vast majority (94%) of the isolates with multiple CYP51 copies carried the mutant Y136F allele in every gene copy. Our results also show that isolates with identical or very similar SSR haplotypes can have contrasting EnCYP51 alleles and copy numbers, which suggests that structural rearrangement in the locus and subsequent selection under fungicide pressure happened frequently and on multiple occasions in the populations we screened. This observation is compatible with findings of high CNV rates over family generational scales .
While the Y136F mutation reduces the affinity between CYP51 and the inhibitor [12, 21], sensitivity to DMIs has also been shown to depend on the expression levels of CYP51. In other plant pathogens the overexpression of CYP51 has resulted in decreased susceptibility to DMIs [15, 67–69]. Segmental duplication of a 126-bp fraction of the CYP51 promoter led to increased expression of CYP51 in a fungicide resistant strain of Penicillium digitatum, while duplication of the entire chromosome bearing the CYP51 gene resulted in increased fungicide resistance in the yeast Candida glabrata. Although polymorphisms have been found in resistant CYP51 genes [71, 72], no clear association between a specific mutation and an increase in CYP51 expression has been reported in fungi.
Based on the strong association between CNV, the Y136F allele, and fungicide treatment, we hypothesize that in E. necator CNV can provide an additional layer of genetic diversity that contributes to pathogen fitness once a favorable allele is acquired. The development of increasing quantitative fungicide resistance could happen in two steps. First, occasional E. necator individuals with the Y136F mutation in EnCYP51 survive DMI fungicide treatments and, as a consequence, Y136F allelic frequency increases in vineyards treated with DMIs. Second, structural rearrangements that increase the number of EnCYP51 copies occur in individuals carrying the Y136F, which are then favored by further DMI treatments. We do not rule out the possibility that structural rearrangements in the EnCYP51 locus may occur with the same frequency in the presence of the wild type allele, but based on the low frequency of CNV in the wild type allele we observed, the acquisition of multiple copies of the wild type alleles does not seem to provide evolutionary fitness in fungicide treated vineyards. Alternatively, the strong connection between the mutation and multiple copies in E. necator could also be a result of a fitness cost due to the Y136F mutation. The mutation occurs in a highly conserved substrate recognition site , and perhaps the substitution leads to a less efficient sterol biosynthesis in E. necator and the increased copy number is needed to compensate the loss of affinity to the substrate. The compensatory effect provided by the acquisition of multiple CNV copies may explain why out of the 94 isolates we studied none carried the Y136F allele as a single copy gene. However, we did not observe any obvious impact of the Y136F mutation on fungal growth when growth rates of wild type and mutant strains in the absence of fungicide were compared (Additional file 16: Figure S6A).
Taken together, our results show that CNV of the EnCYP51 gene contributes to E. necator fitness by providing increasing quantitative protection against DMI fungicide treatments in a gene-dosage dependent manner. Monitoring powdery mildew evolution for the development of fungicide resistance should include measuring the number of the EnCYP51 gene copies. This knowledge will help in managing fungicide programs and maintaining effective control of powdery mildew while minimizing the development of populations increasingly resistant to DMI treatments.
As part of this study we generated valuable genomic information for the grape powdery mildew pathogen, including: (i) whole-shotgun sequence assemblies of five isolates of E. necator, (ii) an estimate of the distribution of genetic diversity as both sequence polymorphisms and structural variants among E. necator isolates, (iii) characterization of the expression of E. necator genes in infected tissues using RNAseq analysis, (iv) gene models of 6,533 protein-coding genes identified with an a hybrid approach that involved both ab initio and reference based predictions using genome-guided de novo assembled transcripts. In addition to providing the first description of the large and highly repetitive genome, in this work we show that frequent structural variations occur between field isolates, including copy number variation in EnCYP51, the target of sterol demethylase inhibitor fungicides. The EnCYP51 gene in isolates collected from fungicide-treated vineyards contained both a mutation known to result in higher fungicide tolerance and an increased copy number that correlated with gene expression and fungal growth in the presence of fungicides. Our results also suggest that CNV can provide an additional layer of genetic diversity that contributes to pathogen fitness once a favorable allele is acquired.
The results of this work not only demonstrate the effectiveness of using genomics to dissect complex traits in organisms with very limited molecular information, but also may have broader implications for understanding genomic dynamics in response to strong selective pressure in other pathogens with similar genome architectures.
A total of 94 E. necator isolates were collected between 2012 and 2013 from seven vineyard sites in California. The DNA of four additional isolates from the eastern United States, kindly provided by Dr. L. Cadle-Davidson (USDA-ARS, Geneva, NY), were included in our analysis to allow comparisons of SSR fingerprint data among different studies (; Additional file 1: Table S1). Using the procedure described in , isolates were single-spore purified and maintained on detached leaves of V. vinifera cv. Carignan. Leaves were surface sterilized for 3 minutes with 3.0% sodium hypochlorite, rinsed with sterile water, dried between layers of sterile paper towels, and placed onto a 0.8% water agar medium in 100 mm × 15 mm Petri dishes with the adaxial side up and the cut petiole inserted into the medium. Four subsequent sub-transfers of a single chain of conidia were made to new leaves to ensure a pure colony was maintained in culture.
Genome sequencing and assembly
Genomic DNA was extracted from mycelia of each isolate collected from colonies growing on detached Carignan leaves approximately 20 days post inoculation. A cyclone separator tube adapter attached to the vacuum port was used to collect the E. necator tissue from the surface of the leaf to minimize plant tissue contamination. DNA was extracted using a modified CTAB protocol , and then fragmented using a Diogenode Bioruptor NGS. Sequencing libraries were prepared using the KAPA LTP Library Sequencing Preparation Kit. Insert size selection was carried out using the Life E-gel SizeSelect to achieve a mean library fragment size of approximately 600 bp for the MiSeq and 550 bp for the HiSeq2500. Insert size and library quality was confirmed before sequencing using the Agilent 2100 Bioanalyzer (Agilent Technologies). Sequencing was carried out on Illumina MiSeq and HiSeq2500 machines at the DNA Technologies Core at UC Davis. Sequenced reads were quality trimmed (Q > 20) using sickle (v1.210; ) and adapters were removed with scythe (version 0.991; ). High quality reads were assembled de novo using CLC Workbench v6.1. Assembly parameters were optimized to achieve maximal assembly completeness of the gene space estimated using the Core Eukaryotic Genes Mapping Approach (CEGMA) analysis . All assemblies were generated using a word size of 64 and a bubble size of 244. Contigs with homology only to non-fungal sequences in the complete NCBI nt collection were considered contaminant and discarded. Basic assembly metrics were extracted using the assemblathon_stats.pl script  and are shown in Additional file 4: Table S2. Genome sizes were estimated based on k-mer count distribution .
Transposable elements annotation
RepeatModeler (v1.0.7; ) was used to perform ab initio repeat family prediction. The scaffolds of the C-strain isolate were used as the only input for RepeatModeler. With the exception of TEs marked as “unknown”, the repeat sequences identified by RepeatModeler were combined with the RepeatMasker database and used as reference for RepeatMasker (v4.0.2; ) to screen the DNA assemblies for interspersed repeats and low complexity DNA sequences.
Transcriptome sequencing, de novoassembly, and gene prediction
Third and fourth fully expanded Carignan leaves were surface sterilized as described above and placed in Petri dishes containing 0.8% water agar medium. A portable paint sprayer was used to inoculate leaves with a 106 conidia/mL suspension of the C-strain isolate in 0.001% Tween. The leaves were maintained under ambient conditions. Time points were collected at 0.5, 1, 3, and 6 days post inoculation. Accumulation of E. necator transcripts was confirmed by qRT-PCR (Additional file 8: Figure S3B; Additional file 17: Table S10). For each time point, three independent pools each containing two leaves were collected. RNA was extracted from combined leaf and fungal tissue using a CTAB method . cDNA libraries were prepared using the Illumina TruSeq RNA Sample Preparation Kit. Library quality was confirmed with the Agilent 2100 Bioanalyzer, and sequencing was carried out on the Illumina HiSeq 2500. Sequence reads were quality filtered and trimmed using scythe and sickle as described above. For gene prediction we applied a hybrid approach: (i) transcripts were de novo assembled using RNA-seq reads and a genome-guided assembly strategy; (ii) the assembled transcripts were employed both as training set for ab initio prediction and as evidence for homology-based prediction together with other fungal proteomes. To identify transcribed loci in the genome for genome-guided transcript assembly, reads were mapped onto the C-strain scaffolds with the spliced (intron-aware) aligner TopHat (version 2.0.9; ) with the following parameters: --r 100 --mate-std-dev 50 --min-intron-length 20. Reads mapped on each transcribed partition were then de novo assembled independently using Trinity (version r20131110; with jaccard_clip option on ). ORFs were extracted using the perl script transcripts_to_best_scoring_ORFs.pl, part of Trinity. De novo assembled transcripts together with the gene structures of the CEGs identified in the C-strain genome by CEGMA were used to train Augustus  for ab initio gene discovery on repeat-masked scaffolds. MAKER (version 2.8; ) was then used to integrate the ab intio prediction with homology based gene prediction using Bgh, Bgt, and B. cinerea peptide sequences as templates (Additional file 8: Figure S3C). Predicted peptides matching known TE associated proteins (BLASTP, e-value < 10−3) were removed. Predicted transcripts not matching any peptide in NCBI nr (BLASTX, e-value < 10−3) were removed unless RNA-seq data mapping provided evidence of their expression. The transcripts of most of these predicted protein-coding genes were present in the RNA-seq data, ranging from 49.1% at 0.5 dpi to 80.6% at 6 dpi of the total 6,533 genes (Additional file 8: Figure S3H-J).
Functional annotation was performed using BLASTP (NCBI nr Asomycota database; e-value < 10−3, −v 10 -b 10, −F F) followed by Blast2GO . Pfam domains were annotated using the Pfam batch search server (; p-value < 10−3), while proteins were assigned to Carbohydrate Active enZYme (CAZy) families based on similarity searches against the non-redundant sequences of the CAZy database using the CAZymes Analysis Toolkit (CAT; ). Proteins were clustered in tribes based on sequence similarity using BLASTP (e-value < 10−6) followed by Markov clustering with TribeMCL . The presence of secretion signal peptides was evaluated using SignalP v.4.1 . Perl scripts were used to scan the first 110 amino acids of each protein for putative effector motifs as described in  (Additional file 10: Table S6).
Gene expression profiling during an infection time course
The quality filtered and trimmed pair-end RNA-seq reads were aligned to a combined reference of the 6,533 predicted protein-coding transcripts and the V. vinifera var. PN40024 genomic scaffolds to profile the E. necator genes during susceptible grape leaf infection. Because RNA-seq samples include both E. necator and Carignan RNAs, a combined reference was used to minimize the non-specific mapping of V. vinifera reads onto E. necator transcripts . The aligner Bowtie2 (version 2.1.0; ) was used for read mapping with parameters --end-to-end --sensitive. SAM output files were parsed to extract read counts per transcript with the python script sam2counts.py (version 0.91; ). Non-specific mappings were discarded by filtering with sam2counts.py (−q 30). DEseq (version 1.12.1, ) was used to normalize raw transcript counts and to compare infection time points.
Assessing genetic diversity
Illumina pair-end genomic sequence reads from each isolate were mapped onto the C-strain scaffold using Bowtie2 (−−end-to-end --sensitive; ). Mapping metrics are reported in Additional file 4: Table S2. Non-uniquely mapped reads and optical duplicates were removed with samtools (version 0.1.19; ) and Picardtools (version 1.104), respectively. The GATK (version 2.4.3) RealignerTargetCreator and IndelRealigner programs were applied to realign the reads mapped on indel sites. The GATK UnifiedGenotyper was then used to identify SNPs with parameters --glm SNP --min_base_quality_score 20 --ploidy 1 --outputmode EMIT_VARIANT_ONLY. VCF files generated by GATK were then parsed with custom perl scripts to extract counts of intergenic, intronic, and exonic variants. To identify genes that are polymorphic and under positive selection we applied the same procedure as described in . Synthetic sequences incorporating the GATK-detected SNPs were generated using the GATK FastaAlternateReferenceMaker tools. Orthologous transcripts were then aligned using ClustalW and subjected to pair-wise polymorphism and positive selection analysis using the program Yn00 . Any pair-wise comparisons that yielded a dN value > 0 were classified as polymorphic and those with dN/dS values > 1 were classified as under positive selection. SSR fingerprinting was carried out as described in Additional file 3: Text S1 (see also Additional files 18, 19 and 20 for additional information). Missing core-eukaryotic genes in the 5 isolates (Additional file 21: Table S14) were identified as described in Additional file 3: Text S1.
Estimation of genome-wide copy number variation
Sequence depth of coverage was used to identify CNV loci in the five E. necator isolates using CNV-seq . Because CNV-seq requires a control sample for comparison, we performed all possible pair-wise comparisons of the 5 genomes. BAM files filtered to exclude ambiguous mapping and optical duplications as described above were used as input for CNV-seq analysis, which was performed with the parameters --global-normalization --genome-size 120000000 --annotate --p-value 0.001. Gene copy numbers were measured using quantitative PCR. DNA was extracted as described above. Additional file 17: Table S10 lists all qPCR amplified genes and their corresponding primers. qPCR was performed on an Applied Biosystems 7500 PCR System using Power SYBR Green Master Mix. All qRT-PCR reactions were run in triplicates and performed with the following cycling conditions: 50°C for 2 min, 95°C for 10 min, followed by 40 cycles of 95°C for 15 s and 60°C for 60 s. Primer efficiencies were calculated using 4-fold DNA dilutions (1:1, 1:4, 1:16, 1:64, and 1:256) in duplicate. The efficiencies of the primer sets used in this study were all above 90% (Additional file 17: Table S10). Specificity of the primers was checked by analyzing dissociation curves ranging from 60 to 95°C. The 2−ΔΔCT method was used to normalize and calibrate CYP51 copy numbers as described in . Data were normalized to the E. necator single copy elongation factor gene EnEF1 (En-g1817; Figure 2A) by subtracting the EnEF1 Ct value and using Branching DNA (1 CYP51 copy) as a calibrator to adjust for any minor variance among the DNA dilution quantities of the different samples as follows: ΔΔCT = (Ct, CYP51 – Ct, EnEF1)x – (Ct, CYP51 – Ct, EnEF1)y, where x = test sample and y = Branching isolate. CYP51 alleles were determined by PCR as described in  and confirmed by Sanger sequencing of PCR amplicons generated using primers listed in Additional file 17: Table S10, and the following PCR conditions: 95°C for 1 min, 46 cycles of 95°C for 15 s, 63°C for 30 s, 68°C for 30 s, and a final extension of 68°C for 7 min.
Measurement of CYP51expression by qRT-PCR
Detached Carignan leaves were inoculated with E. necator conidia and RNA was collected 1 and 5 days after inoculation as described above. cDNA was prepared from the isolated RNA using M-MLV Reverse Transcriptase (Promega). All qRT-PCR reactions were performed as follows: 50°C for 2 min, 95°C for 10 min, followed by 40 cycles of 95°C for 15 s and 60°C for 60 s. The E. necator elongation factor gene (En-g1817 or EnEF1; Additional file 17: Table S10), whose expression was confirmed to be constant across inoculation time points (Additional file 10: Table S6), was used as the reference gene and processed in parallel with the genes of interest. All expression values are presented as fold-EnEF1 (number of target molecules/number of EF1 molecules), calculated using the ΔCt method as described in .
Fungicide resistance testing of E. necatorisolates
Carignan leaves were surface sterilized as described above. A cork borer was used to cut 16 mm diameter leaf disks, and the disks were placed in 6-well tissue culture plates containing 0.8% water agar medium. Rally® 40WSP fungicide (active ingredient: 40% myclobutanil) was applied to the surface of each leaf disk at 7 different concentrations in triplicate with concentrations ranging from 0 to 50 mg/L of the active ingredient. 50 µL of each fungicide dilution was applied to the surface of the leaf disk, and the residual liquid was allowed to evaporate. After 24 hours, the leaf disks were inoculated with one of the five E. necator isolates. Each inoculum was prepared as a suspension of 105 conidia/mL in 0.001% Tween, and 30 µL was pipetted onto each leaf disk. The liquid was allowed to evaporate, and the plates were stored at ambient conditions. Eight days post inoculation, combined leaf and fungal DNA was extracted from each disk using the CTAB method and fungal biomass accumulation was determined by qPCR amplification of the E. necator elongation factor EnEF1. DNA accumulation levels were linearized with the formula 2 – (EnEF1 Ct – VvACT Ct) using the grape actin gene as reference (Additional file 17: Table S10), which was confirmed not to change with the progression of infection, as reference. To account for growth rate differences among isolates in absence of the fungicide, growth results in the presence of the fungicide were calculated as the percentage of the fungal biomass accumulation in the absence of the fungicide. Relative growth values across all fungicide concentrations are shown in Figure 6E and Additional file 16: Figure S6B.
The whole genome shotgun projects have been deposited in the GenBank database [accession nos. JNVN00000000 (isolate C-strain), JNUS00000000 (isolate Branching), JNUT00000000 (isolate Ranch9), JOKO00000000 (isolate e1-101), and JNUU00000000 (isolate Lodi)]. RNA-sequencing data used in this study have been deposited in the National Center for Biotechnology Information Gene Expression Omnibus (GEO) database, http://www.ncbi.nlm.nih.gov/geo (no. GSE58958). The raw reads are available via Sequence Read Archive project no. SRP043708 (accession nos. SRX642773 - SRX642784).
Gadoury DM, Cadle-Davidson L, Wilcox WF, Dry IB, Seem RC, Milgroom MG: Grapevine powdery mildew (Erysiphe necator): a fascinating system for the study of the biology, ecology and epidemiology of an obligate biotroph. Mol Plant Pathol. 2012, 13 (1): 1-16. 10.1111/j.1364-3703.2011.00728.x.
Calonnec A, Cartolaro P, Poupot C, Dubourdieu D, Darriet P: Effects of Uncinula necator on the yield and quality of grapes (Vitis vinifera) and wine. Plant Pathol. 2004, 53 (4): 434-445. 10.1111/j.0032-0862.2004.01016.x.
Stummer BE, Francis IL, Zanker T, Lattey KA, Scott ES: Effects of powdery mildew on the sensory properties and composition of Chardonnay juice and wine when grape sugar ripeness is standardised. Aust J Grape Wine Res. 2005, 11 (1): 66-76. 10.1111/j.1755-0238.2005.tb00280.x.
Tamura K, Peterson D, Peterson N, Stecher G, Nei M, Kumar S: MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Mol Biol Evol. 2011, 28 (10): 2731-2739. 10.1093/molbev/msr121.
Felsenstein J: Confidence limits on phylogenies: an approach using the bootstrap. Evolution. 1985, 39: 783-791. 10.2307/2408678.
Tamura K, Nei M, Kumar S: Prospects for inferring very large phylogenies by using the neighbor-joining method. Proc Natl Acad Sci U S A. 2004, 101 (30): 11030-11035. 10.1073/pnas.0404206101.
Fuller KB, Alston JM, Sambucci O: The value of powdery mildew resistance in grapes: evidence from California. Robert Mondavi Institute Center for Wine Economics. 2014, Davis: UC, 1-28.
Russell PE: A century of fungicide evolution. J Agric Sci. 2005, 143: 11-25. 10.1017/S0021859605004971.
Gubler WD, Ypema HL, Ouimette DE, Bettiga LJ: Occurrence of resistance in Uncinula necator to triadimefon, myclobutanil, and fenarimol in California grapevines. Plant Dis. 1996, 80 (8): 902-909. 10.1094/PD-80-0902.
Hinckley EL, Matson PA: Transformations, transport, and potential unintended consequences of high sulfur inputs to Napa Valley vineyards. Proc Natl Acad Sci U S A. 2011, 108 (34): 14005-14010. 10.1073/pnas.1110741108.
Gisi U, Chin KM, Knapova G, Küng Färber R, Mohr U, Parisi S, Sierotzki H, Steinfeld U: Recent developments in elucidating modes of resistance to phenylamide, DMI and strobilurin fungicides. Crop Prot. 2000, 19 (8–10): 863-872.
Becher R, Wirsel SG: Fungal cytochrome P450 sterol 14alpha-demethylase (CYP51) and azole resistance in plant and human pathogens. Appl Microbiol Biotechnol. 2012, 95 (4): 825-840. 10.1007/s00253-012-4195-9.
Miller TC, Gubler WD: Sensitivity of California isolates of Uncinula necator to trifloxystrobin and spiroxamine, and update on triadimefon sensitivity. Plant Dis. 2004, 88 (11): 1205-1212. 10.1094/PDIS.2004.88.11.1205.
Ypema HL, Ypema M, Gubler WD: Sensitivity of Uncinula necator to benomyl, triadimefon, myclobutanil, and fenarimol in California. Plant Dis. 1997, 81 (3): 293-297. 10.1094/PDIS.19184.108.40.2063.
Ma Z, Proffer TJ, Jacobs JL, Sundin GW: Overexpression of the 14αDemethylase target gene (CYP51) Mediates Fungicide Resistance in Blumeriella jaapii. Appl Environ Microbiol. 2006, 72 (4): 2581-2585. 10.1128/AEM.72.4.2581-2585.2006.
Délye C, Laigret F, Corio-Costet MF: A mutation in the 14α-demethylase gene of Uncinula necator that correlates with resistance to a sterol biosynthesis inhibitor. Appl Environ Microbiol. 1997, 63 (8): 2966-2970.
Délye C, Bousset L, Corio-Costet MF: PCR cloning and detection of point mutations in the eburicol 14alpha-demethylase (CYP51) gene from Erysiphe graminis f. sp. hordei, a “recalcitrant” fungus. Curr Genet. 1998, 34 (5): 399-403. 10.1007/s002940050413.
Sanglard D, Ischer F, Koymans L, Bille J: Amino acid substitutions in the cytochrome from azole-resistant Candida albicans clinical isolates contribute to resistance to azole antifungal agents. Antimicrob Agents Chemother. 1998, 42 (2): 241-253. 10.1093/jac/42.2.241.
Leroux P, Albertini C, Gredt M, Walker A-S: Mutations in the CYP51 gene correlated with changes in sensitivity to sterol 14α-demethylation inhibitors in field isolates of Mycosphaerella graminicola. Pest Manag Sci. 2007, 63: 688-698. 10.1002/ps.1390.
Aoyama Y, Noshiro M, Gotoh O, Imaoka S, Funae Y, Kurosawa N, Horiuchi T, Yoshida Y: Sterol 14-demethylase P450 (P45014DM*) is one of the most ancient and conserved P450 species. J Biochem. 1996, 119 (5): 926-933. 10.1093/oxfordjournals.jbchem.a021331.
Mullins JGL, Parker JE, Cools HJ, Togawa RC, Lucas JA, Fraaije BA, Kelly DE, Kelly SL: Molecular modelling of the emergence of azole resistance in Mycosphaerella graminicola. PLoS One. 2011, 6 (6): e20973-10.1371/journal.pone.0020973.
McDonald BA, Linde C: Pathogen population genetics, evolutionary potential, and durable resistance. Annu Rev Phytopathol. 2002, 40 (1): 349-379. 10.1146/annurev.phyto.40.120501.101443.
Spanu PD: The genomics of obligate (and nonobligate) biotrophs. Annu Rev Phytopathol. 2012, 50: 91-109. 10.1146/annurev-phyto-081211-173024.
Wakefield L, Gadoury DM, Seem RC, Milgroom MG, Sun Q, Cadle-Davidson L: Differential gene expression during conidiation in the grape powdery mildew pathogen, Erysiphe Necator. Phytopathol. 2011, 101 (7): 839-846. 10.1094/PHYTO-11-10-0295.
Frenkel O, Portillo I, Brewer MT, Péros JP, Cadle-Davidson L, Milgroom MG: Development of microsatellite markers from the transcriptome of Erysiphe necator for analysing population structure in North America and Europe. Plant Pathol. 2012, 61 (1): 106-119. 10.1111/j.1365-3059.2011.02502.x.
Brewer MT, Milgroom MG: Phylogeography and population structure of the grape powdery mildew fungus, Erysiphe necator, from diverse Vitis species. BMC Evol Biol. 2010, 10: 268-10.1186/1471-2148-10-268.
Spanu PD, Abbott JC, Amselem J, Burgis TA, Soanes DM, Stuber K, Ver Loren van Themaat E, Brown JK, Butcher SA, Gurr SJ, Lebrun MH, Ridout CJ, Schulze-Lefert P, Talbot NJ, Ahmadinejad N, Ametz C, Barton GR, Benjdia M, Bidzinski P, Bindschedler LV, Both M, Brewer MT, Cadle-Davidson L, Cadle-Davidson MM, Collemare J, Cramer R, Frenkel O, Godfrey D, Harriman J, Hoede C, et al: Genome expansion and gene loss in powdery mildew fungi reveal tradeoffs in extreme parasitism. Science. 2010, 330 (6010): 1543-1546. 10.1126/science.1194573.
Hacquard S, Kracher B, Maekawa T, Vernaldi S, Schulze-Lefert P, Ver Loren van Themaat E: Mosaic genome structure of the barley powdery mildew pathogen and conservation of transcriptional programs in divergent hosts. Proc Natl Acad Sci U S A. 2013, 110 (24): E2219-E2228. 10.1073/pnas.1306807110.
Wicker T, Oberhaensli S, Parlange F, Buchmann JP, Shatalina M, Roffler S, Ben-David R, Dolezel J, Simkova H, Schulze-Lefert P, Spanu PD, Bruggmann R, Amselem J, Quesneville H, Ver Loren van Themaat E, Paape T, Shimizu KK, Keller B: The wheat powdery mildew genome shows the unique evolution of an obligate biotroph. Nature Genet. 2013, 45 (9): 1092-1096. 10.1038/ng.2704.
Simpson JT: Exploring genome characteristics and sequence quality without a reference. Bioinformatics. 2014, 30 (9): 1228-1235. 10.1093/bioinformatics/btu023.
Parra G, Bradnam K, Ning Z, Keane T, Korf I: Assessing the gene space in draft genomes. Nucleic Acids Res. 2009, 37 (1): 289-297. 10.1093/nar/gkn916.
Kupfer DM, Drabenstot SD, Buchanan KL, Lai H, Zhu H, Dyer DW, Roe BA, Murphy JW: Introns and splicing elements of five diverse fungi. Eukaryot Cell. 2004, 3 (5): 1088-1100. 10.1128/EC.3.5.1088-1100.2004.
Baxter L, Tripathy S, Ishaque N, Boot N, Cabral A, Kemen E, Thines M, Ah-Fong A, Anderson R, Badejoko W, Bittner-Eddy P, Boore JL, Chibucos MC, Coates M, Dehal P, Delehaunty K, Dong S, Downton P, Dumas B, Fabro G, Fronick C, Fuerstenberg SI, Fulton L, Gaulin E, Govers F, Hughes L, Humphray S, Jiang RH, Judelson H, Kamoun S, et al: Signatures of adaptation to obligate biotrophy in the Hyaloperonospora arabidopsidis genome. Science. 2010, 330 (6010): 1549-1551. 10.1126/science.1195203.
Kemen E, Gardiner A, Schultz-Larsen T, Kemen AC, Balmuth AL, Robert-Seilaniantz A, Bailey K, Holub E, Studholme DJ, Maclean D, Jones JD: Gene gain and loss during evolution of obligate parasitism in the white rust pathogen of Arabidopsis thaliana. PLoS Biol. 2011, 9 (7): e1001094-10.1371/journal.pbio.1001094.
Duplessis S, Cuomo CA, Lin Y-c, Aerts A, Tisserant E, Grabherr MG, Kodira CD, Kohler A, Kües U, Lindquist EA, Lucas SM, Grigoriev IV, Szabo LJ, Martin F: Obligate biotrophy features unraveled by the genomic analysis of rust fungi. Proc Natl Acad Sci U S A. 2011, 108 (22): 9166-9171. 10.1073/pnas.1019315108.
Kamoun S: A catalogue of the effector secretome of plant pathogenic oomycetes. Annu Rev Phytopathol. 2006, 44: 41-60. 10.1146/annurev.phyto.44.070505.143436.
Spanu P, Kämper J: Genomics of biotrophy in fungi and oomycetes–emerging patterns. Curr Opin Plant Biol. 2010, 13 (4): 409-414. 10.1016/j.pbi.2010.03.004.
Sacristán S, Vigouroux M, Pedersen C, Skamnioti P, Thordal-Christensen H, Micali C, Brown JK, Ridout CJ: Coevolution between a family of parasite virulence effectors and a class of LINE-1 retrotransposons. PLoS One. 2009, 4 (10): e7463-10.1371/journal.pone.0007463.
Pliego C, Nowara D, Bonciani G, Gheorghe DM, Xu R, Surana P, Whigham E, Nettleton D, Bogdanove AJ, Wise RP, Schweizer P, Bindschedler LV, Spanu PD: Host-induced gene silencing in barley powdery mildew reveals a class of ribonuclease-like effectors. Mol Plant Microbe Interact. 2013, 26 (6): 633-642. 10.1094/MPMI-01-13-0005-R.
Pedersen C, van Themaat EVL, McGuffin L, Abbott J, Burgis T, Barton G, Bindschedler L, Lu X, Maekawa T, Weßling R, Cramer R, Thordal-Christensen H, Panstruga R, Spanu P: Structure and evolution of barley powdery mildew effector candidates. BMC Genomics. 2012, 13 (1): 694-10.1186/1471-2164-13-694.
Godfrey D, Bohlenius H, Pedersen C, Zhang Z, Emmersen J, Thordal-Christensen H: Powdery mildew fungal effector candidates share N-terminal Y/F/WxC-motif. BMC Genomics. 2010, 11: 317-10.1186/1471-2164-11-317.
Ma W, Guttman DS: Evolution of prokaryotic and eukaryotic virulence effectors. Curr Opin Plant Biol. 2008, 11 (4): 412-419. 10.1016/j.pbi.2008.05.001.
Girirajan S, Campbell CD, Eichler EE: Human copy number variation and complex genetic disease. Annu Rev Genet. 2011, 45: 203-226. 10.1146/annurev-genet-102209-163544.
Alkan C, Coe BP, Eichler EE: Genome structural variation discovery and genotyping. Nat Rev Genet. 2011, 12 (5): 363-376. 10.1038/nrg2958.
Xie C, Tammi MT: CNV-seq, a new method to detect copy number variation using high-throughput sequencing. BMC Bioinformatics. 2009, 10: 80-10.1186/1471-2105-10-80.
Duan J, Zhang JG, Deng HW, Wang YP: Comparative studies of copy number variation detection methods for next-generation sequencing technologies. PLoS One. 2013, 8 (3): e59128-10.1371/journal.pone.0059128.
Hastings PJ, Lupski JR, Rosenberg SM, Ira G: Mechanisms of change in gene copy number. Nat Rev Genet. 2009, 10 (8): 551-564. 10.1038/nrg2593.
Chaignat E, Yahya-Graison EA, Henrichsen CN, Chrast J, Schütz F, Pradervand S, Reymond A: Copy number variation modifies expression time courses. Genome Res. 2011, 21 (1): 106-113. 10.1101/gr.112748.110.
Chen A, Dubcovsky J: Wheat TILLING mutants show that the vernalization gene VRN1 down-regulates the flowering repressor VRN2 in leaves but is not essential for flowering. PLoS Genet. 2012, 8 (12): e1003134-10.1371/journal.pgen.1003134.
Hacquard S: Chapter Four - The Genomics of Powdery Mildew Fungi: Past Achievements, Present Status, and Future Prospects. Advances in Botanical Research. Edited by: Martin FM. 2014, Oxford, UK: Elsevier, 70: 109-142.
SanMiguel P, Gaut BS, Tikhonov A, Nakajima Y, Bennetzen JL: The paleontology of intergene retrotransposons of maize. Nature Genet. 1998, 20 (1): 43-45. 10.1038/1695.
Mcclintock B: The origin and behavior of mutable Loci in Maize. Proc Natl Acad Sci U S A. 1950, 36 (6): 344-355. 10.1073/pnas.36.6.344.
Harberd NP, Flavell RB, Thompson RD: Identification of a transposon-like insertion in a Glu-1 allele of wheat. Mol Gen Genet. 1987, 209 (2): 326-332. 10.1007/BF00329661.
Akhunov ED, Akhunova AR, Dvorak J: Mechanisms and rates of birth and death of dispersed duplicated genes during the evolution of a multigene family in diploid and tetraploid wheats. Mol Biol Evol. 2007, 24 (2): 539-550.
Lönnig W-E, Saedler H: Chromosome rearrangements and transposable elements. Annu Rev Genet. 2002, 36: 389-410. 10.1146/annurev.genet.36.040202.092802.
Coen ES, Carpenter R, Martin C: Transposable elements generate novel spatial patterns of gene expression in Antirrhinum majus. Cell. 1986, 47 (2): 285-296. 10.1016/0092-8674(86)90451-4.
Zhang Z, Saier MHJ: A novel mechanism of transposon-mediated gene activation. PLoS Genet. 2009, 5 (10): e1000689-10.1371/journal.pgen.1000689.
Raffaele S, Kamoun S: Genome evolution in filamentous plant pathogens: why bigger can be better. Nat Rev Microbiol. 2012, 10 (6): 417-430.
Pang AW, MacDonald JR, Pinto D, Wei J, Rafiq MA, Conrad DF, Park H, Hurles ME, Lee C, Venter JC, Kirkness EF, Levy S, Feuk L, Scherer SW: Towards a comprehensive structural variation map of an individual human genome. Genome Biol. 2010, 11 (5): R52-10.1186/gb-2010-11-5-r52.
Turner TL, Bourne EC, Von Wettberg EJ, Hu TT, Nuzhdin SV: Population resequencing reveals local adaptation of Arabidopsis lyrata to serpentine soils. Nat Genet. 2010, 42 (3): 260-263. 10.1038/ng.515.
Maron LG, Guimaraes CT, Kirst M, Albert PS, Birchler JA, Bradbury PJ, Buckler ES, Coluccio AE, Danilova TV, Kudrna D, Magalhaes JV, Pineros MA, Schatz MC, Wing RA, Kochian LV: Aluminum tolerance in maize is associated with higher MATE1 gene copy number. Proc Natl Acad Sci U S A. 2013, 110 (13): 5241-5246. 10.1073/pnas.1220766110.
Cook DE, Lee TG, Guo X, Melito S, Wang K, Bayless AM, Wang J, Hughes TJ, Willis DK, Clemente TE, Diers BW, Jiang J, Hudson ME, Bent AF: Copy number variation of multiple genes at Rhg1 mediates nematode resistance in soybean. Science. 2012, 338 (6111): 1206-1209. 10.1126/science.1228746.
Nair S, Nash D, Sudimack D, Jaidee A, Barends M, Uhlemann AC, Krishna S, Nosten F, Anderson TJ: Recurrent gene amplification and soft selective sweeps during evolution of multidrug resistance in malaria parasites. Mol Biol Evol. 2007, 24 (2): 562-573.
Schmidt JM, Good RT, Appleton B, Sherrard J, Raymant GC, Bogwitz MR, Martin J, Daborn PJ, Goddard ME, Batterham P, Robin C: Copy number variation and transposable elements feature in recent, ongoing adaptation at the Cyp6g1 locus. PLoS Genet. 2010, 6 (6): e1000998-10.1371/journal.pgen.1000998.
Qutob D, Tedman-Jones J, Dong S, Kuflu K, Pham H, Wang Y, Dou D, Kale SD, Arredondo FD, Tyler BM, Gijzen M: Copy number variation and transcriptional polymorphisms of Phytophthora sojae RXLR effector genes Avr1a and Avr3a. PLoS One. 2009, 4 (4): e5066-10.1371/journal.pone.0005066.
DeBolt S: Copy number variation shapes genome diversity in Arabidopsis over immediate family generational scales. Genome Biol Evol. 2010, 2: 441-453. 10.1093/gbe/evq033.
Hamamoto H, Hasegawa K, Nakaune R, Lee YJ, Makizumi Y, Akutsu K, Hibi T: Tandem repeat of a transcriptional enhancer upstream of the sterol 14α-demethylase gene (CYP51) in Penicillium digitatum. Appl Environ Microbiol. 2000, 66 (8): 3421-3426. 10.1128/AEM.66.8.3421-3426.2000.
Schnabel G, Jones AL: The 14α-demethylase (CYP51A1) gene is overexpressed in Venturia inaequalis strains resistant to myclobutanil. Phytopathol. 2001, 91 (1): 102-110. 10.1094/PHYTO.2001.91.1.102.
Luo C-X, Schnabel G: The cytochrome P450 lanosterol 14α-demethylase gene is a demethylation inhibitor fungicide resistance determinant in Monilinia fructicola field isolates from Georgia. Appl Environ Microbiol. 2008, 74 (2): 359-366. 10.1128/AEM.02159-07.
Marichal P, Vanden Bossche H, Odds FC, Nobels G, Warnock DW, Timmerman V, Van Broeckhoven C, Fay S, Mose-Larsen P: Molecular biological characterization of an azole-resistant Candida glabrata isolate. Antimicrob Agents Chemother. 1997, 41 (10): 2229-2237.
Hulvey J, Popko JT, Sang H, Berg A, Jung G: Overexpression of ShCYP51B and ShatrD in Sclerotinia homoeocarpa isolates exhibiting practical field resistance to a demethylation inhibitor fungicide. Appl Environ Microbiol. 2012, 78 (18): 6674-6682. 10.1128/AEM.00417-12.
Ma B, Tredway LP: Induced overexpression of cytochrome P450 sterol 14alpha-demethylase gene (CYP51) correlates with sensitivity to demethylation inhibitors (DMIs) in Sclerotinia homoeocarpa. Pest Manag Sci. 2013, 69 (12): 1369-1378. 10.1002/ps.3513.
Riaz SBU, Lejkina I, Gubler WD, Walker MA: Report of a new grape powdery mildew morphotype with branched conidiophores. Plant Pathol Quar. 2013, 3 (1): 19-27.
Möller EM, Bahnweg G, Sandermann H, Geiger HH: A simple and efficient protocol for isolation of high molecular weight DNA from filamentous fungi, fruit bodies, and infected plant tissues. Nuc Acids Res. 1992, 20 (22): 6115-6116. 10.1093/nar/20.22.6115.
Joshi N, Fass J: Sickle: A sliding-window, adaptive, quality-based trimming tool for FastQ files. [https://github.com/najoshi/sickle]
Buffalo V: Scythe - a Bayesian adapter trimmer. [https://github.com/vsbuffalo/scythe]
Bradnam K, Fass J, Alexandrov A, Baranay P, Bechner M, Birol I, Boisvert S, Chapman J, Chapuis G, Chikhi R, Chitsaz H, Chou W-C, Corbeil J, Del Fabbro C, Docking T, Durbin R, Earl D, Emrich S, Fedotov P, Fonseca N, Ganapathy G, Gibbs R, Gnerre S, Godzaridis E, Goldstein S, Haimel M, Hall G, Haussler D, Hiatt J, Ho I: Assemblathon 2: evaluating de novo methods of genome assembly in three vertebrate species. GigaScience. 2013, 2: 10-10.1186/2047-217X-2-10.
Smit AFA, Hubley R: RepeatModeler Open-1.0. [http://www.repeatmasker.org]
Smit AFA, Hubley R, Green P: RepeatMasker Open-3.0. [http://www.repeatmasker.org]
Blanco-Ulate B, Vincenti E, Powell AL, Cantu D: Tomato transcriptome and mutant analyses suggest a role for plant stress hormones in the interaction between fruit and Botrytis cinerea. Front Plant Sci. 2013, 4: 142-
Trapnell C, Pachter L, Salzberg SL: TopHat: discovering splice junctions with RNA-Seq. Bioinformatics. 2009, 25 (9): 1105-1111. 10.1093/bioinformatics/btp120.
Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, Adiconis X, Fan L, Raychowdhury R, Zeng Q, Chen Z, Mauceli E, Hacohen N, Gnirke A, Rhind N, di Palma F, Birren BW, Nusbaum C, Lindblad-Toh K, Friedman N, Regev A: Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nature Biotechnol. 2011, 29 (7): 644-652. 10.1038/nbt.1883.
Stanke M, Steinkamp R, Waack S, Morgenstern B: AUGUSTUS: a web server for gene finding in eukaryotes. Nucleic Acids Res. 2004, 32 (Web Server issue): W309-W312.
Cantarel BL, Korf I, Robb SM, Parra G, Ross E, Moore B, Holt C, Sanchez Alvarado A, Yandell M: MAKER: an easy-to-use annotation pipeline designed for emerging model organism genomes. Genome Res. 2008, 18 (1): 188-196.
Conesa A, Götz S, García-Gómez J, Terol J, Talón M, Robles M: Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics. 2005, 21 (18): 3674-3676. 10.1093/bioinformatics/bti610.
Finn RD, Mistry J, Tate J, Coggill P, Heger A, Pollington JE, Gavin OL, Gunasekaran P, Ceric G, Forslund K, Holm L, Sonnhammer EL, Eddy SR, Bateman A: The Pfam protein families database. Nucleic Acids Res. 2010, 38 (Database issue): D211-222-
Park BH, Karpinets TV, Syed MH, Leuze MR, Uberbacher EC: CAZymes Analysis Toolkit (CAT): web service for searching and analyzing carbohydrate-active enzymes in a newly sequenced organism using CAZy database. Glycobiol. 2010, 20 (12): 1574-1584. 10.1093/glycob/cwq106.
Enright AJ, Van Dongen S, Ouzounis CA: An efficient algorithm for large-scale detection of protein families. Nucleic Acids Res. 2002, 30 (7): 1575-1584. 10.1093/nar/30.7.1575.
Petersen TN, Brunak S, Von Heijne G, Nielsen H: SignalP 4.0: discriminating signal peptides from transmembrane regions. Nat Methods. 2011, 8 (10): 785-786. 10.1038/nmeth.1701.
Saunders DG, Win J, Cano LM, Szabo LJ, Kamoun S, Raffaele S: Using hierarchical clustering of secreted protein families to classify and rank candidate effectors of rust fungi. PLoS One. 2012, 7 (1): e29847-10.1371/journal.pone.0029847.
Westermann AJ, Gorski SA, Vogel J: Dual RNA-seq of pathogen and host. Nat Rev Microbiol. 2012, 10 (9): 618-630. 10.1038/nrmicro2852.
Langmead B, Salzberg SL: Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012, 9 (4): 357-359. 10.1038/nmeth.1923.
Buffalo V: sam2counts.py - convert SAM mapping results to reference sequence counts. [https://github.com/vsbuffalo/sam2counts]
Anders S, Huber W: Differential expression analysis for sequence count data. Genome Biol. 2010, 11 (10): R106-10.1186/gb-2010-11-10-r106.
Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, Marth G, Abecasis G, Durbin R, Genome Project Data Processing S: The Sequence Alignment/Map format and SAMtools. Bioinformatics. 2009, 25 (16): 2078-2079. 10.1093/bioinformatics/btp352.
Cantu D, Segovia V, MacLean D, Bayles R, Chen X, Kamoun S, Dubcovsky J, Saunders DG, Uauy C: Genome analyses of the wheat yellow (stripe) rust pathogen Puccinia striiformis f. sp. tritici reveal polymorphic and haustorial expressed secreted proteins as candidate effectors. BMC Genomics. 2013, 14: 270-10.1186/1471-2164-14-270.
Yang Z: PAML 4: phylogenetic analysis by maximum likelihood. Mol Biol Evol. 2007, 24 (8): 1586-1591. 10.1093/molbev/msm088.
Ferreira ID, Rosario VE, Cravo PV: Real-time quantitative PCR with SYBR Green I detection for estimating copy numbers of nine drug resistance candidate genes in Plasmodium falciparum. Malaria J. 2006, 5: 1-10.1186/1475-2875-5-1.
This work was supported by funding to DC from the American Vineyard Foundation (Grants: 2013–1510, 2014–1657) and from the College of Agricultural and Environmental Sciences (UC Davis) and from the Louis P. Martini Endowment to MAW. LJ was supported by the Adolf L. & Richie C. Heck Research Fellowship, the Frank H. Bartholomew Scholarship, the John Ferrington Award, the Robert L. Balzer Scholarship, the Southeast Pennsylvania Region of the American Wine Society Scholarship, the Stewart Good Scholarship, the Wine Spectator Scholarship, and the Women For WineSense Scholarship. We thank Drs. Stephen Pearce (UC Davis) and Cristobal Uauy (John Innes Centre) for insightful comments and Irina Lejkina, Zirou Ye, and Alexandre Carriot for technical assistance.
The authors declare that they have no competing interests
Conceived and designed the experiments: DC, LJ, SR, WDG, MAW. Performed the experiments: LJ, SR, AMC, KCHA, BM. Performed bioinformatic analysis: DC, KCHA. Population structure analysis was carried out by SR. All authors contributed to the analysis of the data and to the review of the paper. Wrote the manuscript: DC, LJ. All authors read and approved the final manuscript.