Author Archives: admin

The contributors

Andrea Harper (CNAP, University of York) and Ian Bancroft (CNAP, University of York)

The analysis

The analysis used computed RPKM transcript abundance data for 130,978 principal isoforms in leaves, representing each of 184 trees from across Denmark and phenotyped for disease symptoms by our colleague Erik Dahl Kjaer’s group. A regression analysis was conducted, by which correlations between gene expression levels and dieback susceptibility were identified. Those with the most significant correlations, P<10-8, are listed.

Materials

Methods

Qualities of raw reads from the thirteen samples were assessed with FASTQC. Adapter- and quality- trimming was performed with Trim galore, a wrapper script using FASTQC and cutadapt (phred cutoff:20, error rate:0.1, adapter overlap: 1bp, min. length: 20bp, paired read length cut-off: 35bp). FASTQC was automatically run on all trimmed files to confirm trimmed read quality. Raw and trimmed metrics and GC content of raw reads can be seen in Table 1.

The variant detection software, VarScan2 pileup2snp (or mpileup2snp) was used ( –p-value 0.05 –min-coverage 10 –output-vcf -min-var-freq 0.95) to call SNPs (Koboldt, et al., 2012). The numbers of SNPs were normalised to the number of bases with >=10X coverage to take into account the different depths of sequencing.

Results

The Fera samples contain fewer sequencing reads compared to the fruiting body samples, mixed material samples and KW1 (Table 1). Sample KW1 covers 32% of the genome positions with a coverage of >=10X compared to between 0.17-0.27% for the Fera samples. The remaining samples range from 13% (AT1) to 32% (Upton) (Table 1). This suggests that limited SNP calling is possible in the FERA samples. The KW1 genome assembly has a size of 65Mb and transcriptome assembly has a size of ~35Mb suggesting that spliced transcripts represent 55% of the genome size.

Sample

Sample type

Raw reads x2

GC%

Filtered x2

% 10X KW1 genome positions covered

KW1

C. fraxinea

41,036,951

49

38,594,279

32%

FERA 105

C. fraxinea

1,442,015

48

1,440,015

0.26%

FERA 88

C. fraxinea

1,731,461

47

1,729,869

0.20%

FERA 93

C. fraxinea

1,355,421

47

1,353,849

0.24%

FERA 94

C. fraxinea

1,642,453

47

1,639,096

0.17%

FERA 232

C. fraxinea

1,919,713

48

1,918,416

0.27%

FERA 233

C. fraxinea

1,838,409

47

1,837,097

0.26%

Mature fruiting body (MFB1)

C. fraxinea

15,033,876

49

14,524,487

27%

Primordial fruiting body (PFB1)

C. fraxinea

7,298,316

49

7,187,891

16%

AT1

Mixed material

31,059,879

48

30,388,033

29%

AT2

Mixed material

23,096,021

44

21,924,824

13%

Holt

Mixed material

36,960,651

49

36.231.834

21%

Upton

Mixed material

43,099,662

48

41,367,071

32%

Table 1: Read number pre- and post-filtering, G+C content of raw data and percentage of KW1 genome positions covered at a depth of >=10X.

SNP calling using VarScan revealed a high level of genetic variation between the samples (Table 2). The mature fruiting body (MFB1) and primordial fruiting body (PFB1) samples showed a high level of variation (~16,000 SNPs) when aligned to the KW1 genome. Samples AT1 (~6,000), AT2 (~3,000), Holt (~14,000) and Upton (~12,000) also suggested a high level of genetic variation although this could be due to the alignment of non-Chalara reads present in the mixed sample, which have a high degree of similarity to regions of the C. fraxinea genome. This highlights the importance of the availability of multiple genomes due to such high levels of genetic variation.

We identified respectively 40, 20, 47, 22, 53, 35 single-nucleotide differences (10X coverage, 95% consensus) for samples FERA88, FERA93, FERA94, FERA105, FERA232, FERA233 when compared to the KW1 genome. In total, 204 unique SNPs were identified in the FERA samples when multi-sample calling was used. P-values were calculated using a Fisher’s Exact Test on the read counts supporting reference and variant alleles. Of the 204 SNPs, 165 were located within an annotated gene, 155 of these within annotated exons. These 155 SNPs were distributed across 70 genes which are we are currently annotating (BLAST & PFAM).

One of these genes (CHAFR746836.1.1_0031990) has a sequence that is similar (62% identity) to cerato-platanin from the basidiomycete Trametes versicolor (GenBank: EIW62259.1) and was first identified in the ascomycete, Ceratocystis fimbriata, the causal agent of “canker stain disease” in Platanus x acerifolia in Europe (Pazzagli, et al., 1999; Pazzagli, et al., 2006). It belongs to a family of cerato-platanin phytotoxic proteins which are found in the cell wall of the fungus (Boddi, et al., 2004) and are involved in the host-pathogen interaction (Pazzagli, et al., 1999). Six SNPs were identified in this protein. All the FERA samples except FERA88 and FERA105 differ from KW1 at all six SNPs; FERA 88 and FERA105 are identical to KW1 at all six SNPs (Fig 1). At these six SNPs, the primordial fruiting body sample resembles FERA 105 and FERA 88 and differs from FERA 90, FERA 232, FERA233, FERA 94 and the mature fruiting body sample. These two fruiting body samples originate from different sources. In the mature fruiting body sample, which is likely to contain both dikaryotic hyphae and haploid ascospores, all 6 SNPs appear to be heterokaryotic. In all other samples the SNPs are homokaryotic/homozygous.

Fig 1: (CLick for larger view) IGV view of putative cerato-platanin gene in C. fraxinea. The six SNPs called by VarScan are shown in the top track. The next four tracks show the homokaryotic/homozygous SNPs which are the same in FERA93, FERA232, FERA233 and FERA94, the sixth track shows possible heterokaryotic SNPs in the mature fruiting body (MFB1) and the last three tracks show that the homokaryotic/homozygous SNPs which resemble the KW1 genome sequence are also present in FERA88, FERA105 and the primordial fruiting body (PFB1).

Normalisation of the identified SNPs shows that the highest variation occurs in fruiting body samples, Holt and Upton when compared to the KW1 genome (Fig 1). Further analysis is required to understand the significance of this difference in variation and whether the Upton and Holt variation is due to the presence of other fungi, such as Phytothphora sp. and Togninia sp. in the mixed samples interfering with the SNP calling.

Fig 2: Normalised number of SNPs against bases with >=10X . Click for larger view

Contributors

Christine Sambles, David Studholme. University of Exeter, Devon.

Introduction

To check their species composition, we performed a metagenomic analysis on several datasets that had been generated from samples described as Fraxinus excelsior, Chalara fraxinea or as mixed infected material. We used MEGAN (v4, MEtaGenome Analyzer (Huson, et al., 2011)) and the assembled transcripts from F. excelsior and C. fraxinea to identify the taxonomic groups to which uninfected sample transcripts are allocated, to use as a reference database for binning. We identified many transcripts specific to the infected samples (not in uninfected Fraxinus) that fell outside of the Chalara and Fraxinus bins; these may indicate additional microbial species present during infection. These might include species acting synergistically with C. fraxinea during infection or opportunists present as a secondary consequence of infection. The uninfected F. excelsior sample identifies species that are part of the ‘normal’ or ‘healthy’ tree microbiota and could therefore be excluded from the list of infection-related species.

Material

Methods

We identified sequence similarity between assembled transcripts and GenBank protein sequences using BLASTX; we used as queries the transcripts from uninfected F. excelsior (ATU1), a C. fraxinea isolate (KW1) and four mixed material samples (AT1, AT2, Upton and Holt). We loaded the output from BLASTX into MEGAN and performed taxonomic binning using a minimum support value of 35, a minimum BLAST score of 50 and only retaining hits whose bit scores lie within 10% of the best score. The analyses were normalised, compared and rendered within MEGAN.

Results

For each of the six samples, we identified the numbers of transcripts assigned to Helotiales (the order in which C. fraxinea is classified) and to Viridiplantae (green plants). The results are summarised in Table 1. Transcripts from the (nominally) C. fraxinea isolate KW1, binned into the classes of Dothideomycetes, Eurotiomycetes, Leotiomycetes and Sordariomycetes, which all reside within the subphylum of Pezizomycotina. This result is consistent with the sequenced sample being pure Chalara. As expected, all of the bin-able transcripts from the F. excelsior ATU1 transcripts fell within the Viridiplantae kingdom, specifically within the group of flowering plants (Magnoliophyta).

Table 1: Reads binned to Helotiales and Viridiplantae in normalised comparison for each sample.

In the data from Upton mixed material, 74% of the transcripts were not binned within the Helotiales or Viridiplantae, which is where C. fraxinea and F. excelsior transcripts are expected to fall, based on the results from pure isolate of C. fraxinea and the uninfected F. excelsior. In the Upton data, 34% of the total number of transcripts was assigned to Oomycetes; specifically 33% to Phytophthora spp.. Additionally, 13% are not assigned to any taxon and a further 11% had no significant similarity to proteins in the GenBank database, detectable by BLAST. The presence of Phytophthora spp. might be attributed to cross-lane contamination during sequencing, since the Norwich laboratory handling the Upton data also work with Phytophthora infestans. A similar contamination had also been reported for Fera samples with Maize Chlorotic Mottle Virus (MCMV) and Sugarcane Mosaic Virus (SMV) sequences being present. This is a common problem in Illumina sequencing which has led to the incorporation of a taxonomic binning step into sequencing pipelines including at our own sequencing facility at the University of Exeter. These contamination issues highlight the importance of confirming the taxonomic distribution of sequence data in addition to quality checks before performing any downstream analyses. Once identified, the contaminant reads can be removed from the dataset.

The Holt mixed material sample analysis showed 1.5% transcripts binned to Togninia minima, an ascomycete in the order Calosphaeriales. T. minima is a pathogen of grapevines and Prunus spp., however, the closely related T. fraxinopennsylvanica (anamorph: Phaeoacremonium mortoniae) has been observed in dead vascular tissue of declining ash tree branches (Fraxinus latifolia) in California (Eskalen, et al., 2005a; Eskalen, et al., 2005b). It may be that T. fraxinopennsylvanica is present in the Holt material and that the reads were assigned to the species T. minima because that is the most closely related species for which extensive sequence data is available.

For the nominally pure sample of C. fraxinea isolate KW1, only 32% of reads were assigned to Helotiales. However, 66% of reads were assigned to Fungi with 16% not assigned to a taxa and 17% had no significant similarity to proteins in the GenBank BLAST database. This is likely to be due to insufficient sequence data in the GenBank database from Chalara and closely related species.

Further analysis using alignments to the F. excelsior and C. fraxinea genome will help interpret whether or not other taxa such as Togninia sp., indicated to be present by the MEGAN analysis, are present or whether they are mis-assigned reads due to the lack of Chalara– and Fraxinus-related proteins in the database.

The interpretation

The genome assembly contains numerous low GC, high coverage, repeat rich regions, with low overall gene density. This suggests an assembly with collapsed repeats of a genome that is overall rich in repeats.

Cite this:

The contributors

Dan MacLean

The analysis

To extend the analysis of Chalara’s ability to degrade wood and thus invade a tree directly I examined the lignin-decaying gene complement in the genome. The previous work on identifying decay related enzymes in [Floudas]: http://dx.doi.org/10.1126/science.1221748 was used as a base.

For each group in the Floudas decay-related list Chalara proteins were assayed to identify proteins in the with >50% sequence identity to the representative protein in at least 10 members. Counts of the number of Chalara proteins passing this threshold were used as the estimate of the number of members of each group in the Chalara genome.

Material

Analysis:

The AT1, AT2, ATU, Upton and Holt mixed-material RNASeq assemblies were searched against the Genbank nr database with blastx in order to investigate the virome of each sample. A number of novel viral sequences were detected. A blastx was then carried out against the Kenninghall wood (KW) mycelia assembly to determine whether any of these viral sequences were likely to be infecting Chalara fraxinea. While none of the novel viruses originally discovered in the other datasets were present in the KW assembly, a 2.4kb mitovirus was detected. A tblastx search of this mitovirus sequence against the mixed-material samples discovered shorter contigs in the Holt and Upton assemblies which had homology to the KW mitovirus but which had not been detected in the original blastx against Genbank.
Fera has been testing large numbers of samples as part of the response to the presence of Ash Dieback in the UK and so we sequenced six cultured isolates of Chalara fraxinea (RNASeq with MiSeq) to investigate the prevalence of mitoviruses in a small number of our cultures. Mitoviruses were detected in the Fera 88, 93 and 94 samples but not Fera 105, 232 and 233.

Interpretation:

The presence or absence of mitovirus does not appear to be geographically correlated but the Fera samples were not comprehensively sequenced so mitoviruses may be present in these samples but not sequenced. Further work is required to investigate whether these mitoviruses confer hypovirulence to Chalara fraxinea.

The contributors

The material

First flush leaf material was harvested from Tree 35 in Denmark in May 2013.

The analysis

mRNA was extracted and then a strand-specific, paired-end Illumina RNA-Seq library constructed. Around 193 million read pairs were obtained from a single HiSeq 2500 lane and these are now available from The Sainsbury Laboratory’s FTP server, with details in the github repository here.

A blastp search of the amino acid code using the FASTA sequence of Trypanosomal alternative oxidase (TAO) as a template found significant alignment with an E value of 6e-62.

MitoProt II identified a putative mitochondrial targeting sequence.

SWISS-MODEL was used to model the putative C. fraxinea AOX protein using trypanosomal alternative oxidase as a template (PDB 3VV9), with a QMEAN Z-score of -5.119, which is consistent with monotopic membrane proteins.

Read and assembly sequences were analysed to calculate GC content in reads and different genomic feature types. This figure shows that the gDNA reads and contigs/scaffolds have a skewed, almost bimodal GC content. The skew is absent from RNA sequence and gene feature sequences, indicating a difference in GC content in the non-genic regions of the genome.

were then analysed for the presence of small Compositionally Non-Homogenous Domains (GC-rich regions) using IsoPlotter 2.3. in Elhaik et al. This figure shows each of the scaffolds as black bars in which the white regions indicate statistically significant regions of non-homoegeneity, according to the method in in Elhaik et al.

These plots show in greater detail the changes in GC percent across these scaffolds. Many show strong drops in the GC content.

The pattern of GC poor patches observed here is reminiscent of the ‘isochore-like’ regions observed around the regions carrying effector genes in Leptosphaeria maculan’s that are believed to be involved in the plasticity of that organisms genome and may contribute to increased evolutionary potential and an ability to change hosts quickly see Raffaele and Kamoun, 2012.