The European race of Gremmeniella abietina hosts a single species of Gammapartitivirus showing a global distribution and possible recombinant events in its history.

The population genetics of the family Partitiviridae was studied within the European race of the conifer pathogen Gremmeniella abietina. One hundred sixty-two isolates were collected from different countries, including Canada, the Czech Republic, Finland, Italy, Montenegro, Serbia, Spain, Switzerland, Turkey and the United States. A unique species of G. abietina RNA virus-MS1 (GaRV-MS1) appears to occur indistinctly in G. abietina biotypes A and B, without a particular geographical distribution pattern. Forty-six isolates were shown to host GaRV-MS1 according to direct specific RT-PCR screening, and the virus was more common in biotype A than B. Phylogenetic analysis based on 46 partial coat protein (CP) cDNA sequences divided the GaRV-MS1 population into two closely related clades, while RNA-dependent RNA polymerase (RdRp) sequences revealed only one clade. The evolution of the virus appears to mainly occur through purifying selection but also through recombination. Recombination events were detected within alignments of the three complete CP and RdRp sequences of GaRV-MS1. This is the first time that recombination events have been directly identified in fungal partitiviruses and in G. abietina in particular. The results suggest that the population dynamics of GaRV-MS1 do not have a direct impact on the genetic structure of its host, G. abietina, though they might have had an innocuous ancestral relationship.

Introduction

The European race of Gremmeniella abietina (Lagerberg) Morelet (anamorph Brunchorstia pinea) is the causal agent of stem canker and shoot blight in many conifer species in Europe and North America (Donaubauer 1972). It leads to epidemics depending on annual weather conditions (Uotila 1988). The most recent serious outbreak of G. abietina took place in Sweden in 2001 (Wulff & Walheim 2003), when more than 50 000 ha of pine forest were severely affected. The up-to-date taxonomy of this species comprises three biotypes. Biotype B, or the Small Tree type (STT), is primarily found in northern European countries, but it appears to also be present in high mountainous areas of the Lake District of Turkey (Botella et al., 2012a, Botella et al., 2012b, Dogmus-Lehtijarvi et al., 2012). The Alpine biotype has been described in the Alpine region. These biotypes both grow in harsh conditions of snow and extreme temperatures and produce abundant pycnidia (4–8 celled) and apothecia in the field (Hellgren & Högberg 1995). The third currently known biotype is the biotype A, or the Large Tree Type (LTT). It is the most pathogenic biotype and shows the widest distribution, ranging from the Italian Apennines to northern Sweden and North America, where it was introduced in the 1970s and designated the European race (EU) (Skilling ). Conversely, it produces apothecia only rarely in forests (Uotila, 1992, Kaitera and Jalkanen, 1996), and its conidia are 4 celled. Population genetic studies performed through analyses of SCAR markers have consistently found highly similar alleles in the B and Alpine biotypes (Botella ), while Biotype A displays fixed nucleotide differences from the other two (Botella ). Biotype A exhibits low overall genetic variability (Hamelin ), with the exception of the Spanish population, which is a clonal population in a warmer location, likely as a consequence of a founder effect (Santamaria et al., 2005, Botella et al., 2010). Artificial pairings between the A and B biotypes revealed low rates of successful germination and growth, and accordingly, no viable hybrids have been found in nature (Uotila ), indicating that these biotypes should be separated into two different species (Uotila et al., 2000, Uotila and Hantula, 2012).

Tuomivirta et al. 2002 described the occurrence of two independent double-stranded (ds) RNA banding patterns in a single G. abietina biotype A isolate, which were later shown to represent putative members of the families Totiviriridae and Partitiviridae (Tuomivirta and Hantula, 2003a, Tuomivirta and Hantula, 2005). Only biotype B appears to harbour a taxonomically unassigned species with certain endornavirus-like characteristics (Tuomivirta ), and both biotypes A and B host mitoviruses (Tuomivirta and Hantula, 2003b, Botella et al., 2012b). Mitovirus-infected isolates appear to present greater mycelial growth than mitovirus-free isolates at the optimal growing temperature for G. abietina of 15 °C (Romeralo ). The viral communities appear to be separated into biotypes A and B (Hantula & Tuomivirta 2003), except in Spain, where two species of mitoviruses occur in the genetically unique population, which was most likely derived from biotype A.

Members of the family Partitiviridae are doubled-stranded RNA (dsRNA) viruses that are hypothesised to have evolved from the genus Totivirus (Ghabrial, 1998, Sabanadzovic et al., 2009). The family Partitiviridae includes dsRNA viruses with bisegmented genomes, which encode separately encapsidated RdRP and CP segments as well as satellite RNA molecules. Partitiviruses infect fungal, plant, and protozoan hosts. The classification of partitiviruses has been recently revised, and they now comprise five families: Alphapartitivirus, Betapartitivirus, Cryspoviruses, Deltapartitivirus, and Gammapartitivirus. The genera Alphapartitivirus and Betapartitivirus include viruses hosted by fungi or plants, while the genera Gammapartitivirus and Cryspovirus comprise only mycoviruses (Nibert ).

Transmission of partitiviruses, similar to other mycoviruses, may occur vertically through sexual and asexual spores or horizontally via anastomosis between related fungal genotypes (Ihrmark et al., 2002, Tavantzis, 2010). However, a number of studies have demonstrated the capacity of partitiviruses to be transferred through somatically incompatible fungi under laboratory conditions (Ihrmark et al., 2002, Vainio et al., 2010, Vainio et al., 2011, Vainio et al., 2012a) and in nature (Vainio ). Fungal and plant partitiviruses appear to show many similarities in their virion structure and close evolutionary relationships. Despite the lack of information about plant partitiviruses, a recent study based on RdRp amino acid sequences suggested that horizontal transfer of partitiviruses has occurred during the evolutionary history of fungi and plants (Gabrial ).

We previously analysed the European race of G. abietina (Botella ) and its community of mitoviruses (Botella et al., 2012a, Botella et al., 2012b). The present study builds on our prior objective of determining the evolutionary history of mycoviruses and their influence on the phylogeny of G. abietina, this time focussing on the family Partitiviridae. More specifically, this study has the following three main objectives: (1) to measure the prevalence of partitiviruses in the European race of G. abietina; (2) to analyse their genetic diversity and population structure; and (3) to assess evolutionary processes, such as recombination and selection, to better understand possible host-virus coevolution.

Materials and methods

Fungal collection

A total of 162 isolates belonging to the European race of Gremmeniella abietina were screened for the presence of possible partitiviruses. Overall, 68 biotype A isolates were examined: four isolates from the Czech Republic, 52 isolates from Finland, six isolates from Italy, five isolates from Serbia, and one isolate from Montenegro. For biotype B, there were 30 isolates from Finland, two isolates from Turkey (biotype B/Alpine), and one isolate from Switzerland (Alpine). From the Spanish population, 50 isolates were studied, and finally, from the European race in North America, there were 10 isolates from Canada and one isolate from the United States. The data collected for all of the strains are presented in Supplementary Table S1, while only virus-containing isolates are listed in Table 1 .

Table 1

Fungal isolates and corresponding GaRV-MS1 strains hosted.

Gremmeniellaisolates	Taxonomicalgroup	Location/Country	Tree host	GaRV-MS1 strain	CPd	RdRpd
Ka05-T2_2	Biotype A	Karhula/Finland	P. sylvestris	GaRV-MS1-33	KJ786397	KJ786353
Ke05-14	Biotype A	Keräkankare/Finland	P. sylvestris	GaRV-MS1-44	KJ786408	KJ786364
SurS4	Biotype A	Nummi-Pusula, Finland	P. sylvestris	GaRV MS1-2	Tuomivirta & Hantula 2005	Tuomivirta & Hantula 2005
C5	Biotype A	Nuttio/Finland	P. sylvestris	GaRV-MS1-1	Tuomivirta & Hantula, 2003a	Tuomivirta & Hantula, 2003a
IT12-1-1	Biotype A	Speziac, Italy	P. nigra	GaRV-MS1-24	KJ786388	KJ786344
IT12-3-7	Biotype A	Speziac, Italy	P. nigra	GaRV-MS1-19	KJ786383	KJ786339
Ko12-2	Biotype A	Kopaonik NPe, Serbia	P. sylvestris	GaRV-MS1-25	KJ786389	KJ786345
Ko12-3	Biotype A	Kopaonik NPe, Serbia	P. sylvestris	GaRV-MS1-23	KJ786385	KJ786341
Ko12-5	Biotype A	Kopaonik NPe, Serbia	P. sylvestris	GaRV-MS1-20	KJ786384	KJ786340
Ko12-8	Biotype A	Kopaonik NP4, Serbia	P. sylvestris	GaRV-MS1-23	KJ786387	KJ786343
DU01	Biotype A	Durmitor, Montenegro	P. sylvestris	GaRV-MS1-22	KJ786386	KJ786342
H1-4	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-3	KJ786411	KJ786412
P1-1	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-12	KJ786376	KJ786332
P1-2	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-4	KJ786367	KJ786323
P1-5	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-13	KJ786377	KJ786333
P1-8	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-15	KJ786379	KJ786335
P1-10	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-17	KJ786381	KJ786337
P1-11	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-7	KJ786371	KJ786327
P1-13	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-16	KJ786380	KJ786336
P1-14	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-9	KJ786373	KJ786329
P1-16	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-10	KJ786374	KJ786330
P3-3	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-40	KJ786404	KJ786360
P3-4	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-36	KJ786400	KJ786356
P3-6	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-37	KJ786401	KJ786357
P3-7	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-35	KJ786399	KJ786355
P3-9	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-38	KJ786402	KJ786358
P3-18	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-39	KJ786403	KJ786359
P4-1	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-5	KJ78636	KJ786324
P4-12	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-34	KJ786398	KJ786354
00P-3	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-26	KJ786390	KJ786346
00P-7	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-27	KJ786391	KJ786347
09P-1	SP	Astudillo, Spain	P. halepensis	GaRV-MS1-29	KJ786393	KJ786349
H4-4	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-8	KJ786372	KJ786328
Pha-16	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-30	KJ786394	KJ786350
Zo1	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-11	KJ786375	KJ786331
H4-2	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-31	KJ786395	KJ786351
H4-15	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-32	KJ786396	KJ786352
H4-33	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-6	KJ786370	KJ786326
H8-2	SP	Valle de Cerrato, Spain	P. halepensis	GaRV-MS1-28	KJ786392	KJ786348
NLC 508	Eu race	St.Jonh's, NFa, Canada	P. nigra	GaRV-MS1-43	KJ786407	KJ786363
NLC 509	Eu race	St.Jonh's, NFa, Canada	P. nigra	GaRV-MS1-42	KJ786406	KJ786362
NLC 510	Eu race	St.Jonh's, NFa, Canada	P. nigra	GaRV-MS1-41	KJ786405	KJ786361
QU-85-1074	Eu race	Qcb, Canada	P. resinosa	GaRV-MS1-14	KJ786378	KJ786334
QU-92-0584B	Eu race	Qcb, Canada	P. resinosa	GaRV-MS1-18	KJ786382	KJ786338
TU1	Biotype B	Yenişarbademli, Turkey	P. nigra ssp. pallasiana	GaRV-MS1-45	KJ786409	KJ786365
TU2	Biotype B	Yenişarbademli, Turkey	P. nigra ssp. pallasiana	GaRV-MS1-46	KJ786410	KJ786366

Newfoundland.

Quebec.

Passo di Centocroci Varese Ligure Spezia.

NCBI accession numbers.

Kopaonik National Park.

Total RNA extraction and cDNA synthesis for population analysis

All of the fungi were cultivated on malt extract agar medium (MEA) supplemented with cellophane and incubated at 19 °C. After two weeks, total RNA was extracted and purified using the commercial RNeasy® Plant Mini Kit (Qiagen), with a few modifications. Approximately 3 mg of freshly grown mycelium was scraped and subjected to freeze using liquid nitrogen. The mycelium was first ground in 2 ml RNase-free tubes with TissueLyser for 1 min at 30 Hz. This step was repeated after adding the RNeasy Lysis Buffer complemented with 10 μl of β-mercaptoethanol. The rest of the process was developed following the instructions of the manufacturer. RNA quality and quantity were tested through agarose gel electrophoresis and using a Qubit® 2.0 Fluorometer (Invitrogen), respectively.

A total of 100 ng of RNA was employed for the synthesis of first-strand cDNA using the recombinant M-MuLV Reverse Transcriptase (Thermo Scientific First Strand cDNA Synthesis kit). The sample aliquot was denatured at 65 °C for 5 min and immediately cooled in liquid nitrogen. The concentrations of pD(N)6 Random Hexamer primers, dNTP, and 5× buffer were as recommended by the manufacturer. Reverse transcription-polymerase chain reaction (RT-PCR) was successively carried out twice with Phusion High-Fidelity DNA polymerase (Thermo Fisher Scientific). In the case of both CP and RdRp RT-PCR, the final reaction volume was 25 μl, with the following components: 0.5 U of the enzyme, 1× Phusion HF buffer, 200 μM dNTPs, 0.5 μM each specific primer (PAR_CP_1F/PAR_CP_1R and SpRdRpgap1F/SpRdRpgap1R) (Table 2) and 2 μl of cDNA (<250 ng). The thermocycling conditions consisted of 30 s at 98 °C for initial denaturation, followed by 35 cycles of 10 s at 98 °C, 30 s at 60–65 °C (CP) or 55 °C (RdRp), and 30 s at 72 °C, with a final extension for 10 min at 72 °C. The amplified regions were 559 bp (CP) and 536 bp (RdRp) (Fig 1 ).

Table 2

Specific RT PCR primers designed for the study.

Primer code	Sequence 5′- 3′a	Genome region
PAR_CP_1F	GGCGTGGCTGCCCATGAGTT	CP
PAR_CP_1R	AGCCCGTTGRGACACGGGRA	CP
PARTI7_1R	GGCGTGGCTGCCCATGAGTT	CP
PARTI7_1F	CAGGAATCCGCCAACGCCGA	CP
PAR CP_4F	GCCAGCAGACGAGAAAAGAG	CP
PAR CP _4R	ATTTTCAAGCCTCCATCACG	CP
PAR 3	AATCCCAAAACACGCTTCTCG	CP
SpRdRp5′R	GACGAAACGCCTTCTTTGCTTT	RdRp
SpRdRpGap1F	ATGGACCAACCAGATACTCCGATG	RdRp
SpRdRpGap1R	GATCCTCTGTTTTCTCGCACTTCTC	RdRp
SpRdRpGap2R	GCT GTG AGC TGG CTC CG	RdRp
SpRdRpGap2F	TGT TCA ACG GAC GGC GAT G	RdRp
SpRdRpGap3R	GAGTTCAGCCTCGAAAAAGCA	RdRp

The following designations are used for degenerate sites: Y (C or T), R (A or G), K (G or T),W(A or T).

Fig 1

Representation of the genome organisation of GaRV-MS1-3 and recombination events occurred between strains. *The actual breakpoint position is undetermined (it was most likely overprinted by a subsequent recombination event). Minor Parent, parent contributing the smaller fraction of sequence. Major Parent, parent contributing the larger fraction of sequence. Unknown, only one parent, and a recombinant need be in the alignment for a recombination event to be detectable. The sequence listed as unknown was used to infer the existence of a missing parental sequence.

The RT-PCR products were separated via electrophoresis in 1 % agarose gels (SERVA) containing 1 % TBE (1 mM Tris-Boric acid at pH 8.0, 1 mM EDTA at pH 7.8, 1 mM Tris base at pH 8.0) and 10 μl of GelRed™ Nucleic Gel Acid Gel Stain (Biotium). The PCR products were visualised under UV light after 1 h and 30 min of running at 95 V/30 cm in 1 % TBE buffer solution. Gene Ruler™ 100 bp DNA Ladder Plus (Thermo Scientific) was used as a length marker. Finally, the GeneJET PCR Purification kit (Thermo Scientific) was employed to purify the PCR products, and the template DNA concentrations were determined with a Qubit® 2.0 Fluorometer (Invitrogen) following the instructions of the manufacturer.

DsRNA extraction and determination of the full-length sequence

The Spanish isolate of Gremmeniella abietina H1-4 was chosen for determination of the full-length sequence of a putative new strain of GaRV-MS1 (GenBank accession numbers for the CP, RdRp, and the unknown protein III: KJ786411–KJ786413). The double-stranded RNAs (dsRNAs) of the two corresponding linear molecules (1.4–3.0 kb) were purified, integrating some modifications into the cellulose affinity chromatography protocol of Morris & Dodds (1979) and Tuomivirta . Approximately 3 g of freshly growing mycelium was collected and frozen with liquid nitrogen in 50 ml RNase-free tubes. Once all of the nitrogen was evaporated, the mycelium was homogenised with 8 mm-diameter stainless steel beads (Retsch, Germany). The tube was vortexed on a standard vortex mixer at maximum velocity for 3 min. Then, 2.5 ml of lysis buffer (2.5 ml Tris–HCl pH 8.0 (Carl Roth GmbH + Co. KG, Germany), 0.5 ml of 1 % β-mercaptoethanol (Serva Electrophoresis, Germany), 2.5 ml of EDTA (Penta, Petr Švec, Czech Republic), 15 ml of SDS (pH 7,2; Carl Roth GmbH + Co. KG, Germany), and 25 ml of sterile distilled water (SDW) were added. The dsRNA-binding cellulose added was Sigmacell cellulose (C6288) (Sigma–Aldrich; Germany). The protocols for double-stranded (ds) cDNA synthesis and cloning were developed as described in related works (Tuomivirta and Hantula, 2003a, Tuomivirta and Hantula, 2003b, Botella et al., 2012b).

DNA sequencing

Sequencing of the cloned cDNA fragments and RT-PCR products for population analyses was conducted by Macrogen Europe. The universal primers M13 Forward (17 mer) and M13 Reverse (17 mer) were used. Moreover, it was necessary to design a set of specific primers (Table 2) to fill some gaps within the sequences of CP and RdRp as well as a T4 RNA-primer for sequencing the ends. Finally, to produce more reliable consensus sequences, the purified PCR products from two different reactions were sequenced in both directions. All of the sequences were submitted to GenBank (Table 1).

All of the obtained sequences were adjusted, compiled, and aligned using the Geneious Pro 5.5.9 software package (Biomatters Ltd). NCBI Protein Blast (Blastx) was employed to search for similar sequences and conserved domains. NCBI ORF Finder (http://www.ncbi.nlm.nih.gov/gorf/gorf.html) and Geneious Pro 5.5.6 were used to search for ORFs. The Vienna RNA web suite (Gruber ) was employed to predict secondary structures (http://rna.tbi.univie.ac.at).

Phylogenetic analyses and genetic diversity

MEGA 5.2.1 (Tamura ) was used to perform calculations to select the best substitution model for the construction of the phylogenetic tree based on the CP and RdRp amino acid sequences. The optimal substitution pattern was determined according to the Bayesian information criterion of the program, which specified the Jones-Taylor-Thornton (JTT) matrix-based method as the best model using a discrete Gamma distribution (+G) with five rate categories and assuming that a certain fraction of sites are evolutionarily invariable (+I). The tree topologies obtained through neighbour-joining (NJ) and maximum likelihood (ML) phylogenetic analyses were compared but were found to be highly similar. Therefore, only NJ trees are presented here (Fig 2 , Supplementary S1 and S2). For these analyses, gaps and missing data were completely deleted from the alignment. There were a total of 186 and 166 positions in the final CP and RdRp datasets, respectively.

Fig 2

Optimal unrooted NJ tree based on 46 partial nucleotide CP sequences of GaRV-MS1. Branches corresponding to partitions reproduced in less than 90 % bootstrap replicates are collapsed. GaRV-MS1 strains hosted in biotype A of G. abietina are represented by ■, those hosted in Biotype B are represented by ●, ▼ represents strains in EU race isolates, and finally, ♦ represents the strains in the Spanish population.

DnaSP v5 (Librado & Rozas 2009) was used to estimate genetic diversity parameters for the population. The average number of nucleotide differences per site between two sequences (Pi), number of haplotypes (h), haplotype diversity (Hd), number of segregating sites (S), number of total mutations (h) and the average number of differences per population were calculated.

Evolution and recombination

DnaSP v5 was employed for the analyses of synonymous and non-synonymous nucleotide substitutions, according to the Nei-Gojobori model using the Jukes & Cantor (1969) correction (Nei & Gojobori 1986). We calculated the rate of synonymous per synonymous or silent sites (dS), the rate non-synonymous per non-synonymous substitutions or non-silent sites (dN), and the ratio ω (dN dS −1). Sites (or codons) with alignment gaps or missing data were not used.

The occurrence of recombination was tested employing the default conditions of the suite of programs included in RDP3 (Martin ). Different detection methods are computed in this program: RDP, GENECONV, Bootscan, MaxChi, Chimaera, SiScan, PhylPro, LARD, and 3Seq. Significant recombination events were considered to be those detected by at least three different methods, after Bonferroni correction for multiple comparisons. RDP3 provides a detailed breakdown of recombination breakpoint locations and the identities of the recombinant and parental sequence.

Results

Genome organisation of GaRV-MS1-3

The presence of putative partitivirus molecules was previously demonstrated in the Spanish population of Gremmeniella abietina (Botella ). Among the positive isolates, H1-4 was chosen for the determination of a possible new strain. Following the nomenclature practices described by Tuomivirta & Hantula (2003a and 2005), it was designated G. abietina RNA virus MS1 strain 3 (GaRV-MS1-3). The complete genome of GaRV-MS1-3 was segmented into three putative particles (Fig 1). The largest segment showed a length of 1781 bp and a relatively high GC content (49.7 %). It exhibited a unique open reading frame (ORF) of 1620 bp, from position 63 to 1682. The 5′ and 3′ untranslated regions (UTR) were 62 bp and 99 bp long, respectively. Based on translation frame 3 and the standard code, the RdRp presented an ORF of 590 amino acids (aa) and a molecular weight of 62.155 kDa. The second longest segment presented a length of 1586 bp and a moderately higher GC content (50.9 %). There was only one ORF of 1332 bp identified, from position 100 to 1431 (Fig 1), which encodes the coat protein (CP). The length of the 5′ UTR was 100 bp, while that of the 3′ UTR was 155 bp. The CP presented a length of 444 aa and a molecular weight of 47.229 kDa. The smallest segment was 1187 bp long and exhibited a GC content of 44.1 %. The ORF, which encodes a hypothetical protein (unknown protein III), was 714 bp long, and the 5′ and 3′ UTRs were 347 and 126 bp, respectively. Based on frame 3 and the standard translation code, the unknown protein III consists of 237 aa, with a molecular weight of 26.564 kDa (Fig 1). Potential stem-loop and panhandle structures could be deduced from the RNA sequences of the three proteins when the entire molecules and their individual 5′ and 3′ UTRs were analysed in silico (data non shown).

Comparison of the CP, RdRp, and the unknown III amino acid sequences of GaRV-MS1-3 with the corresponding sequences of other partitiviruses in GenBank (Table 3 ) revealed high similarity to GaRV-MS1 and 2, which were previously described in Finnish isolates of G. abietina biotype A (Tuomivirta and Hantula, 2003b, Tuomivirta and Hantula, 2005), confirming that GaRV-MS1-3 belongs to the same species as GaRV-MS1-1 and GaRV-MS1-2 (Tuomivirta and Hantula, 2003b, Tuomivirta and Hantula, 2005).

Table 3

Relation of the most similar CP, RdRp, and unknown III protein sequences to GaRV-MS1-3 found in Genbank after Blastx search.

Genbank numbers	Description/corresponding virus	Identity
RdRp
NP_659027.1	Gremmeniella abietina RNA virus MS1-1	530/539 (98%)
YP_138540.1	Gremmeniella abietina RNA virus MS1-2	530/539 (98%)
YP_052856.2	Penicillium stoloniferum virus S	389/532 (73%)
YP_001686789.1	Botryotinia fuckeliana partitivirus 1	379/533 (71%)
ABV30675.1	Aspergillus ochraceous virus	364/534 (68%)
AGI52210.1	Verticillium dahliae partitivirus 1	365/532 (69%)
AGO04402.1|	Ustilaginoidea virens partitivirus	362/532 (68%)
NP_116716.1|	Discula destructiva virus 1	356/532 (67%)
NP_620301.1|	Discula destructiva virus 2	352/532 (66%)
CAY25801.2|	Aspergillus fumigatus partitivirus 1	354/535 (66%)
CAJ31886.1|	Ophiostoma partitivirus 1	352/535 (66%)
AFX73023.1|	Grapevine partitivirus	347/456 (76%)
NP_624350.1|	Mycovirus FusoV (Fusarium solani f. sp. robiniae)	320/516 (62%)
CP
YP_138541.1	Gremmeniella abietina RNA virus MS1-2	431/443 (97%)
NP_659028.1	Gremmeniella abietina RNA virus MS1-1	429/443 (97%)
YP_052857.1	Penicillium stoloniferum virus S	250/434 (58%)
CAJ01910.1	Penicillium stoloniferum virus S	250/434 (58%)
NP_620302.1	Discula destructiva virus 2	242/443 (55%)
YP_001686790.1	Botryotinia fuckeliana partitivirus 1	251/443 (57%)
ABV30676.1	Aspergillus ochraceous virus	247/434 (57%)
AGO04403.1	Ustilaginoidea virens partitivirus	222/443 (50%)
NP_116742.1	Discula destructiva virus 1	229/443 (52%)
CAZ61323.2	Aspergillus fumigatus partitivirus	222/444 (50%)
CAJ31887.1	Ophiostoma partitivirus 1
AGI52209.1	Verticillium dahliae partitivirus 1	204/409 (50%)
NP_624351.1	Mycovirus FusoV (Fusarium solani f. sp. robiniae)
YP_001686791.1	Botryotinia fuckeliana partitivirus 1	74/198 (37%)
Unknown III protein
NP_659029.1	hypothetical protein GaRvM1s3gp1 [Gremmeniella abietina RNA virus MS1-1]	232/237 (98%)
YP_138542.1	hypothetical protein GaRVMS2s3gp1 [Gremmeniella abietina RNA virus MS1-2]	229/237 (97%)
AGO04405.1	hypothetical protein [Ustilaginoidea virens partitivirus]	69/217 (32%)

Genetic variability of GaRV-MS1 population

The genetic parameters of CP and RdRp in the GaRV-MS1 population indicated a relatively low degree of genetic variation. Concerning the CP sequences, within a total of 46 strains presenting a sequence length of 559 bp (net sites), the total number of segregating sites (S) was 67, including 7 singleton and 60 parsimony-variable sites. The number of mutations occurring in the population (η) was 74; the number of haplotypes (h) was 21; the haplotype diversity (Hd) was 0.93; the nucleotide diversity estimated based on the average number of differences per site between two sequences (Pi) was 0.034; and the average number of differences was 19.19. For RdRp, the nucleotide variability was even lower; the 46-sequence alignment displayed 515 net sites (excluding sites with gaps/missing data). The total number of segregating sites (S) was 51; the total number of mutations (η) was 57; the number of haplotypes (h) was 18; the haplotype diversity (Hd) was 0.77; the nucleotide diversity estimated based on the average number of differences per site between two sequences (Pi) was 0.018; and the average number of differences was 9.34.

Taken together, these data confirm that GaRV-MS1 is highly conserved within the European race of Gremmeniella abietina.

Phylogenetic relationships of GaRV-MS1 within the European race of Gremmeniella abietina

The phylogenetic NJ trees based on the partial CP and RdRp sequences of 46 strains of GaRV-MS1 are shown in Fig 2 and Supplementary Figures S1 and S2. Analysis of RdRp clustering did not reveal any major clusters supported with high bootstrap values. However, in the case of CP sequences, two major clusters were observed. The reliability of the phylogram was supported by the bootstrap test, as the differentiated clusters were confirmed in more than 90 % of replicates. The first cluster was composed of eight strains belonging to Gremmeniella abietina isolates from three taxonomically distinct groups: the EU race from Newfoundland (Canada), biotype A from Finland, and biotype B from Turkey. The second, more numerous cluster was composed of 38 GaRV-MS1 strains hosted by biotype A and Spanish G. abietina isolates from Canada, Italy, Montenegro, Serbia, and Spain. The number of CP amino acid differences per site determined through calculations of the mean diversity between and within the two clusters in the phylogenetic tree were 0.012 and 0.022, respectively (in both cases, SE was 0.003).

According to this result, the taxonomic complexity of G. abietina and GaRV-MS1 is not correlated.

Evolutionary forces and recombination

To obtain a better understanding of the evolutionary dynamics of GaRV-MS1, the selection pressure on protein-coding regions was assessed based on the ω ratio (dN dS −1) for the 46 partial CP and RdRp sequences (strains) of GaRV-MS1 and the three complete CP, RdRp, and unknown III protein sequences. In the case of the partial CP strains, there were a total of 138.43 synonymous sites (SS) and 419.57 non-synonymous sites (NSS) on average. The RdRp partial sequences hosted a total average of 117.66 SS and 395.34 NSS. The average SS and NSS within the complete RdRp sequences were 177.17 and 533.83, respectively, and within the complete CP sequences, we observed 341.78 SS and 987.22 NSS. In the case of the three complete unknown III protein sequences from GaRV-MS1-1, 2, and 3, there were 341.78 SS and 987.22 NSS. The substitution rate (ω) was below one in all cases (Table 4 ), suggesting the existence of active purifying selection on these sequences.

Table 4

Selection pressure on protein-coding regions of GaRV-MS1 CP, RdRP, and unknown III protein.

Protein	Strains (sequences)	N	dS	dN	ω
RdRp	Partial	46	0.06	0.005	0.08
RdRp	Complete	3	0.11	0.012	0.11
CP	Partial	46	0.13	0.01	0.08
CP	Complete	3	0.14	0.01	0.07
Unknown III protein	Complete	3	0.11	0.01	0.09

No recombination events were detected within the GaRV-MS1 population when the partial CP and RdRp sequences were used for the analysis, nor were such events detected among the complete unknown III protein sequences of GaRV-MS1-1, 2 and 3. However, when the three complete RdRp sequences of GaRV-MS1-1, 2 and 3 were analysed, one recombination event was observed to have occurred between GaRV-MS1-2, as the ‘major parent’ (the sequence that contributed most to the recombinant), and GaRV-MS1-3, as the ‘minor parent’. In this analysis, the minor parent was noted as “unknown”, which means that only one parent and a recombinant are required in the alignment to detect a recombination event. GaRV-MS1-2 may be the actual recombinant (Fig 1A). Three detection methods resulted in significant P values: MaxChi (3.12 10−02), Chimaera (6.29 10−03) and, 3Seq (1.94 10−02). The analysis indicated possible misidentification of the initial breakpoint and recombinant.

In the case of the complete CP sequences, two recombination events appeared. In both cases, the recombinant sequence was GaRV-MS1-2 (Fig 1B). The first recombination event showed significant P values when using three detection methods: GENECONV (3.32 10−02), MaxChi (8.57 10−04), and 3Seq (3.50 10−02). The ‘major parent’ was GaRV-MS1-1, which was detected as ‘unknown’. The ‘minor parent’ was GaRV-MS1-3. The P values of the second recombination event were also significant for three detection methods RDP (0.038), Bootscan (0.02) and 3Seq (7.22 10−03). In this case, the ‘major’ and ‘minor’ parents were GaRV-MS1-3 and GaRV-MS1-1, respectively.

Occurrence of GaRV-MS1 in the European race of Gremmeniella abietina

Among a total of 162 Gremmeniella abietina isolates analysed for the occurrence of GaRV-MS1, 46 (or 28.75 %) of the fungal isolates hosted the virus. When the occurrence of GaRV-MS1 was analysed within each population/biotype, the highest incidence was found in the Spanish population (56 %), followed by the EU race in NA (45.45 %), biotype A (16.17 %) and, finally, biotype B (and the Alpine biotype; only 6.1 %).

The percentages varied slightly when the prevalence was analysed by geographical region considering the total studied isolates: region I (North of Europe: Finland), region II (Central Europe: Czech Republic and Switzerland), region III (Mediterranean Europe: Turkey, Serbia, Montenegro, Italy, and Spain) and region IV (North America: USA and Canada). The highest prevalence occurred within region III (22.9 %), followed by region IV (3.08 %) and region I (2.5 %). In region II, no occurrence of partitiviruses was detected.

Discussion

Most mycoviruses are ‘cryptoviruses’ (Tavantzis 2010); that is, they do not produce phenotypic changes in their fungal hosts (Ghabrial 1998). However, some mycoviruses, such as Cryphonectria hypovirus 1 (CHV-1), display parasite-like traits, such as hypovirulence (Milgoom & Cortesi 2004), and have typically evolved with their hosts (Bryner ). Other mycoviruses, such as Ustilago maydis virus H1, are involved in symbiotic interactions, which strongly affect the evolutionary dynamics and community structure of both symbiotic partners (Voth ). One of the basic relevant questions in evolutionary biology is the degree to which the diversification of parasites and/or symbionts is associated with the diversification of their hosts (Göker ). In this sense, the diversification of the European race of Gremmeniella abietina into different biotypes, and possible species (Uotila et al., 2000, Uotila and Hantula, 2012), is mostly the outcome of geographical separation and ecological adaptation. However, the virus communities of these fungi could provide some insights about their history of dispersion. According to our results, only one species of Gammapartitivirus (Nibert ) appears to occur in the European race of G. abietina (GaRV-MS1), and no host specialisation was observed. The NJ phylogenetic tree based on partial CP sequences showed two very close clades of GaRV-MS1 strains, and the nucleotide variability of the entire population was quite low (0.034). The RdRp alignment exhibited even lower variability. Highly similar strains appeared indistinctly in different countries and in both biotypes A and B (Fig 2). The absence of codivergence between GaRV-MS1 and its host, G. abietina, might be due to an apparent lack of either negative (parasitic) or positive (symbiotic) effects. Although a possible effect cannot be totally ruled out, Tuomivirta did not observe major differences between dsRNA-containing (partitivirus-like patterns) and dsRNA-free isolates. Thus, GaRV-MS1 would not be acting as an evolutionary force driving the species diversification of G. abietina. Nevertheless, combinations of events such as host switching, duplication, and parasite extinction may lead to topological incongruence and divergence (Göker ).

Virus transmission may have occurred from biotype A to B or vice versa, as hypothesised for Gremmeniella mitovirus 1 (GMV1) and 2 (GMV2) found in the Spanish population (specialised A population) and biotype B (Botella ). Transmission of partitiviruses between somatically incompatible fungi may also occur in nature (Vainio ). However, in the case of GaRV-MS1, it is difficult to assess where and when cell-to-cell contact might have taken place. The populations examined in our study are geographically separated, diminishing the possibility of recent transmission.

Purifying selection appears to be the main evolutionary force acting on GaRV-MS1, as amino acid substitutions are silenced (Kimura 1983) in the RdRp, CP and the unknown III protein. This pattern is assumed to occur in other fungal viruses as well, such as Ustilago maydis virus H1 (Voth ), Heterobasidion RNA virus 6 (HetRV6) (Vainio ), and GMV 1 and 2 (Botella et al., 2012a, Botella et al., 2012b), whereas fungal viruses such Cryphonectria hypovirus 1 (CHV1) are subjected to positive selection (Feau ). Additionally, recombination events were detected between the complete CP sequences of GaRV-MS1-1 and 3, suggesting that one or more of the ancestors of GaRV-MS1-2 was a recombinant. Likewise, recombination was identified in the full-length RdRp sequences of GaRV-MS1-1 and 2, with GaRV-MS1-3 indicated to be a recombinant. Because the donor sequences precisely substituted a homologous region of the receptor sequence without altering in its structure (Nuss 2005), we propose that homologous recombination occurred. RNA recombination events are common among plant ssRNA viruses (Chare, 2003, Alabi et al., 2010, Sztuba-Solińska et al., 2011) and human and animal dsRNA viruses such as rotaviruses (Suzuki et al., 1998, Phan et al., 2007), but there is less evidence of recombination in the dsRNAs of plant (Liu ) or fungal viruses (Carbone et al., 2004, Feau et al. 2014, Linder-Basso et al., 2005, Voth et al., 2006), and this phenomenon is particularly scarce in fungal partitiviruses. Nevertheless, some studies carried out on coronavirus have indicated that recombination might be more common than is often assumed but may go undetected because of the action of strong purifying selection, which will remove new, deleterious combinations of mutations (Banner & Lai 1991). Hence, in the case of GaRV-MS1, we note that recombination could have taken place under specific conditions, permitting the virus to cross the species border. However, it is not possible to determine when this event occurred. Limited information is available about possible factors enhancing viral recombination in fungi. Nevertheless, the RNA-silencing antiviral defence response has been shown to contribute to hypoviral RNA recombination in Cryphonectria parasitica (Sun ). The RNA silencing pathway contributes to viral RNA recombination by providing 5′ and 3′ fragments of viral RNA (Zhang & Nuss 2008).

The prevalence of GaRV-MS1 was found to be moderate (28.4 %) and lower than that of mitoviruses (74 %) (Botella ). However, it was relatively high in comparison with other pathogens such as Heterobasidion, where the frequency of partitivirus-infected isolates was shown to be ca. 5 % (Vainio ). The occurrence of GaRV-MS1 varied among the tested populations. Biotype A and the Spanish population in particular presented the highest rates, while biotype B and the Alpine type displayed a lower occurrence of partitiviruses. The same phenomenon was observed in GMV1 and GMV2 (Botella ), though these findings were most likely enhanced by the lower number of analysed isolates, especially from the Alpine region; nevertheless, a further explanation may be based on their mode of reproduction. Biotype A mainly reproduces via conidiospores, which is a known means of efficient transfer of mycoviruses (Polashock et al., 1997, Buck, 1998, Pearson et al., 2009), while biotype B and the Alpine type show sexual reproduction more frequently. The percentage of partitivirus-infected sexual and asexual spores can fluctuate greatly within the same species, and it is apparently determined by the combination of fungal host and virus genotypes (Tavantzis 2010). Furthermore, mycoviruses that flow through septal pores may also be constrained by fungal defence mechanisms (Leslie and Zeller, 1996, Cortesi and Milgroom, 1998, Nuss, 2011), such as the genetic self/nonself recognition system or RNA silencing (Nuss 2005). Evidence of such constraint was recently described after infection by Rosellinia necatrix partitivirus 2 (RnPV2) (Chiba ), although no phenotypic alterations were observed.

In conclusion, the partitivirus community of the European race of G. abietina appears to be composed of a single species (GaRV-MS1) with low genetic variability, which is seemingly stable within the different populations of the fungal host. This virus evolves not only through purifying selection but also to some extent via recombination. To the best of our knowledge, this is the first report of the occurrence of recombination in fungal partitiviruses.