Identification of a SNP and development of a PCR-based allele-specific marker of the sporulation-deficient (sporeless) trait of the Tamogitake 108Y2D mutant using next-generation sequencing.

The mass scattering of basidiospores during the cultivation of edible mushrooms causes serious problems, such as allergic reactions in workers. Sporulation-deficient (sporeless) cultivars would be very useful for preventing these issues. We aimed to identify the single-nucleotide polymorphism (SNP) that is responsible for the single dominant sporeless mutation of the Tamogitake 108Y2D mutant using next-generation sequencing (NGS) and TILLING technology and to develop an allele-specific PCR marker for sporeless breeding. By comparing the sequences of the wild-type and its mutant genomes, we identified 685 mutation loci in gene regions and pinpointed one SNP only consistent with sporeless phenotype for 105 segregants, i.e., a C to T located at position 1,950 of the exonic region of a putative fungal transcription factor that generated a stop codon. We developed an allele-specific marker based on the identified SNP, and its high practicality was validated using tests against progenies from several hybrids and wild isolates from different geographical origins. Thus, the allele-specific PCR marker developed here will be useful for marker-assisted selection in the breeding of the sporeless trait of this mushroom. Furthermore, the technical success of SNP identification and marker development based on NGS genome data can help achieve efficient mutation breeding in mushrooms.

Introduction

In 2018, the financial value of mushroom production in Japan reached 246 billion JPY (https://www.maff.go.jp/j/tokei/kouhyou/tokuyo_rinsan/index.html), accounting for over 85% of the forestry output value. The current mushroom production has been expanded and diversified from hilly and mountainous areas to city suburbs because of the development of cultivation technologies in air-conditioned facilities. Conversely, a large number of basidiospores derived from the developing fruiting bodies, which are essential for the sexual reproduction of mushrooms, have led to various adverse effects in production facilities, as well as in the natural environment. As one of the most serious examples, mushroom growers may develop an allergic pulmonary disease (commonly known as mushroom worker’s lungs) after inhaling scattered spores (Hausen , Olsen 1987). In addition, the abundant spores that are released give rise to various other problems, such as damage to production facilities, reduced commercial value (because of the attachment of spores to fruiting bodies), spread of mushroom disease such as mycovirus infection (Grogan ), and depletion of genetic diversity in the natural population of the mushroom species that are cultivated (Hibbett and Donoghue 1996, Obatake ). These problems are particularly serious in the cultivation of Pleurotus species, such as Pleurotus ostreatus (Olsen 1987), Pleurotus pulmonarius (Ohira 1979), Pleurotus eryngii (Obatake ), Pleurotus florida (Hausen ), and Tamogitake (Pleurotus cornucopiae var. citrinopileatus) (Gisusi , Yoneyama ), because their spore release begins at a very early stage of the development of fruiting bodies and persists throughout maturation.

Among the approaches used to solve the problems caused by spore dispersion, perhaps the most effective method is the breeding of sporulation-deficient (sporeless) cultivars using sporeless mutants. Such mutants have attracted attention as they can be used as breeding material. The spontaneous sporeless mutants have been found to date in Coprinopsis cinerea (Day 1959, Gibbins and Lu 1982), Schizophyllum commune (Bromberg and Schwalb 1977), Lentinula edodes (Hasebe ), Pleurotus ostreatus (Eger ), and Pleurotus pulmonarius (Ohira 1979). The sporeless mutants have also been produced by artificial mutation processes, such as ultraviolet (UV) irradiation or the application of chemical mutagens such as nitrosoguanidine, include those of C. cinerea (Inage , Takemaru and Kamada 1972), Agrocybe cylindracea (Kiuchi 1998, Murakami 1998), P. ostreatus (Imbernon and Labarére 1989, Pandey and Ravishankar 2010), P. pulmonarius (Imbernon and Labarére 1989), Pleurotus eryngii (Obatake ), and Pleurotus sajor-caju (Pandey and Ravishankar 2010). In addition, by cytological study, Kanda and Murakami (1998) suggested that gene mutations in various different sporulation processes lead to sporeless mutation. Although many studies on the sporeless mutation have been done as described above, breeding of sporeless mutant cultivars with usefulness for commercial production has only been achieved in three mushroom species, i.e., A cylindracea (Murakami 1993), P. ostreatus (Baars , Baars and Hesen 2008) and P. eryngii (Obatake ).

In the breeding of mutants, substantial amounts of labor and time are spent in selecting mutant traits as well as isolating effective mutants. Thus, marker-assisted selection (MAS) using DNA markers is often adopted as a strategy that enables more efficient plant and animal breeding (Al-Samarai and Al-Kazaz 2015, Eathington , Liu and Cordes 2004, Snowdon and Friedt 2008). Also in mushrooms, a crossing is main strategy in breeding, and the production of fruiting bodies through the cultivation is necessary to evaluate the hybrids obtained by crossing. Harvesting fruiting bodies in the cultivation takes from one to several months. So, efficiency of MAS using DNA markers can also be applied to breeding mushrooms that are sporeless (Okuda , 2011, 2012). To determine mutation carriers, it is necessary to generate actual fruiting bodies; a substantial amount of labor and time can be saved if discrimination of sporeless traits becomes possible in the mycelium. In this regard, the development of DNA markers based on the gene that carries the sporeless mutation would allow the highly accurate detection of breeding mushrooms harboring the mutation. However, reports of genes involved in sporulation or in the sporeless phenotype are very limited (Okuda ). Therefore, to the best of our knowledge, the DNA marker developed for P. pulmonarius (Okuda ) is the only marker based on a causative gene, stpp1, which belongs to the MutS homolog gene family associated with mismatch repair in DNA synthesis or recombination in meiosis.

The annual production of Tamogitake is approximately 350 tons (Yoneyama ), which is not very large among the cultivated mushrooms of Japan. However, accumulating research, including human intervention trials and epidemiological studies of its beneficial effects, has led to an increasing demand for this mushroom as a health-promoting food (Chen , Fujino and He 1998, Hagiwara , Meng , Misaki , Tanaka ). In recent years, cultivation efforts have spread from Hokkaido, the primary production area, to the lesser-farmed regions of Chugoku, Shikoku, and Kyushu in Japan. Together, the problems on the expansion of spores described above have been occurring (Gisusi , Yoneyama , 2017). Thus, Gisusi identified one spontaneous sporeless mutant of Tamogitake that was estimated to be controlled by a single recessive gene; the authors bred the sporeless cultivar. In addition, Okuda developed linkage markers to detect the mutation mentioned above. However, this mutant had some adverse traits for commercial cultivation, such as the unfavorable morphology, or low yield of fruiting bodies and it’s not easy to improve these disadvantages because of recessive gene even if it was crossed; thus, its cultivation has not spread. Therefore, a new cultivar carrying a novel sporeless mutation not linked to traits that are adverse for cultivation is in strong demand. To address this issue, Yoneyama isolated a novel sporeless mutant (108Y2D, Fig. 1) with an extremely minor effect on cultivation traits by performing mutagenesis using UV irradiation of protoplasts prepared from dikaryon of Tamogitake. The number of spores falling from the mature form of this strain was less than 0.1% of that of the wild-type (Yoneyama ). Furthermore, the inheritable character of this mutant was revealed to be dominant and controlled by a single gene mutation (Yoneyama ). Although the function of the genes causing this mutation has not been clarified, cytological studies of basidia of this sporeless mutant have shown that they are impaired in post-meiotic stage of sporulation process (Yoneyama ).

Fig. 1.

Representative fruiting bodies of the breeding sporeless, Tamogitake No. 108Y2D (Pleurotus cornucopiae var. citrinopileatus) (A). Wild-type phenotype (white deposits (spore print) on the black paper after a few hours to overnight) (cultivar Pc05-1, B) and sporeless phenotype (no spore prints are detected) (No. 108Y2D, C).

This novel sporeless mutant is expected to be effective as a breeding material to resolve the fundamental problem of spore dispersion in commercial cultivation facilities of Tamogitake mushrooms, as described above (Fig. 1).

Recent developments in DNA sequencing technologies, such as next-generation sequencing (NGS), have provided techniques for reading multiple genomic sequences and comparing them, for the detection of polymorphisms (Davey ). Such advances in NGS methods are extremely useful in the mutation analysis of mushrooms. The reason is that the mutant genome and its original wild-type type genome can be isolated experimentally (Matsumoto ) and compared whole genome sequences of them directly. TILLING (Targeted Induced Local Lesions in Genomes) was developed as a method to detect the presence of mutation of a target gene among a mutant population using PCR, and has become available for a variety of uses (Comai and Henikoff 2006, Till ).

This study aimed to introduce the molecular techniques such as MAS in the breeding of Tamogitake sporeless cultivar using the novel mutant (Yoneyama ). We employed NGS and TILLING technology to identify the SNP with matching sporeless phenotype of the test-cross progenies. Then, we developed an allele-specific PCR marker based this polymorphism. Furthermore, we discussed about the gene that harbors the identified SNP in relation with this sporeless mutation.

Materials and Methods

Fungal strains

Yoneyama isolated a sporeless dominant mutant, 3B108, of the Tamogitake cultivar HfpriPc05-1 (Hfpri is omitted below, composition of the haplophase nucleus: Y1 + Y2) via UV irradiation. The authors developed the 108Y2D mutant (Yoneyama ) as a start material for following sporeless breeding by di-mon mating between 3B108 (composition of the haplophase nucleus: LP17 + Y2) and Y2. An LP17 was named for Y1 in which a sporeless mutation occurred by UV irradiation (Yoneyama ). The monokaryotic isolates Y1, Y2, and LP17 were prepared from Pc05-1 and 108Y2D using the protoplasting method (Matsumoto ), respectively. In the genome-wide mutation analysis described below, genomic DNA prepared from Y1 and LP17 monokaryons was used as the wild-type and mutant genome, respectively. The sporulation phenotype (Fig. 1) and mating type were determined for 105 progenies isolated from the 108Y2D sporeless mutant. Total of 110 monokaryons and 22 dikaryons were used for the validation test for allele-specific marker (Table 1, Supplemental Table 1). Out of monokaryons, 101 isolates and 9 isolates were randomly isolated from fruiting bodies of hybrids produced by the crossing between sporeless mutant-type and wild-type and the crossing between wild-types, respectively (Table 1). Dikaryons were wild isolates collected from fields (Supplemental Table 1). These isolates are maintained at the Forest Products Research Institute, Hokkaido Research Organization (Hokkaido, Japan).

Table 1.

Validation of the allele-specific marker for the detection of the sporeless trait using monokaryons from hybrids and dikaryons from wild isolates

Materials	Phenotype observed		Discrimination by real time PCR		Rate of correct
	spo–	spo+	spo–	spo+	(%)
F4 generation (monokaryons)
J06182a	25	15	25	15	100
J31035a	10	13	10	13	100
J58068a	21	17	21	17	100
F1 generation (monokaryons)
Y7Y12b	0	9	0	9	100
Wild isolatesc (dikaryons)	0	22	0	22	100

The phenotypes of each material isolate were judged based on the presence or absence of the spore print shown in Fig. 1. Observations of the spore print were carried out using the fruiting body obtained in the cultivation of dikaryons which were crossed between monokaryons from each hybrid and an optional tester monokaryon (wild-type). RT-qPCR was performed as described in the text.

F4 generation, monokaryons from the 4th cross between sporeless mutant and wild-type (spo− × spo+).

F1 generation, monokaryons from the 1st cross between wild-types (spo+ × spo+).

Wild isolates were dikaryon kept at Hokkaido Research Organization, Forest Products Research Institute (Supplemental Table 1).

DNA extraction for NGS and allele identification of the sporeless locus by PCR

For NGS, the Y1 and LP17 monokaryotic isolates were cultured on MYG medium (2% malt extract, 2% glucose, and 0.2% yeast extract) and lyophilized mycelia were prepared as described by Okuda . Their genomic DNA was extracted using the modified cetyl trimethyl ammonium bromide (CTAB) method (Yoder and Scheffer 1969). For allele identification among monokaryotic isolates using PCR, the Y1, LP17, and 105 progenies were cultured at 25°C for 2 weeks on MA medium (2% malt extract and 2% agar) overlaid with cellophane, and then mycelia were harvested and lyophilized for DNA extraction. Genomic DNAs of each isolate mycelia were extracted using DNeasy plant mini kit (Qiagen, Hilden, Germany).

NGS, genome assembly, and genome-wide mutation analysis

For reference sequence analysis, 10 μg of total DNA from Y1 monokaryon was used as the starting material. Sequencing libraries were constructed and sequenced on a PacBio RS II apparatus using DNA fragmented to a length of 20 kb, followed by selection, sequencing, and de novo assembly via SMRT analysis (Van Buren ). The Y1 genome was used as a reference sequence and for the development of markers. For re-sequencing of the LP17 genome, to detect SNPs and INDELs, sequencing libraries were constructed and sequenced using a HiSeq 2500 instrument (Illumina, San Diego, Ca, USA) and 100-bp paired-end sequencing. SNP and INDEL discovery were performed using the Genedata Expressionist Genomic Profiling 9.0.10 software (Genedata, Lexington, MA, USA). The analyzed sequences were purified using Trimmomatic 0.32 (Bolger ) to remove low-quality regions and sequence-adapter-derived results. Subsequently, sequence reads were mapped to the reference sequence from the Y1 genome using BWA-MEM (Li 2013). Low-quality and duplicate reads were removed from the mapping results and sequences were realigned. The differences between the reference sequence and the read sequence were detected as genetic mutations resulting from UV irradiation (Yoneyama ). The above-mentioned NGS, sequence assembly, and genome-wide mutation analysis were performed by Takara Bio Inc. (Otsu, Shiga, Japan).

Annotation and estimation of the mode of the sporeless mutation in silico

Gene prediction and annotation analysis of the Y1 genome sequence were conducted in silico. First of all, genes and coding sequences from the assembled contigs, including the Y1 genome, were estimated using AUGUSTUS (http://bioinf.uni-greifswald.de/augustus/) with default parameters and Laccaria bicolor as the model species. The amino acid sequences of the estimated genes were searched against the National Center for Biotechnology Information (NCBI) non-redundant protein database using the basic local alignment search tool (blastp). Gene ontology (GO) annotations of the genes were determined using the Blast2GO program (https://www.blast2go.com). Next, mutation patterns were estimated by comparing the relationship between the mapped mutation sites and the gene structure. This information analysis was performed by Takara Bio Inc.

Linkage analysis via TILLING and sequence analysis

To identify the contig that harbored the sporeless mutation site, a linkage analysis was performed between the sporeless phenotype and one to two mutation regions per contig using the TILLING method. Forty-eight progeny DNAs and their parental DNAs (from Y1 and LP17) were used in the TILLING analysis. Bulks (Supplemental Table 2) comprised four individuals each that possessed the sporeless trait or the wild-type trait based on sporulation phenotype data from Yoneyama . An equal amount of DNA from every individual in the bulk was mixed for the TILLING analysis.

Nucleotide sequences including mutation site in contig 12 for the test cross progenies and their parents were amplified by PCR using primer pair, 129FW (5ʹ-GACAGGACTCTCTAGTGATC-3ʹ) and 129RV (5ʹ-TCGCAGTATGGCTCTAACTG-3ʹ). PCR condition used were same as that written in TILLING experiment described above. The PCR product was electrophoresed and purified using by MinElute Gel Extraction Kit (QIAGEN, Hilden, Germany), and then subjected to sequencing. DNA sequencing was performed using a BigDye Terminator Cycle Sequencing Ready Reaction Kit with a Prism 3100 genetic analyzer (Applied Biosystems, Foster City, USA) with amplification primers.

PCR amplification was performed in a 30 μL volume containing 10 ng of pooled DNA, 10× ExTaq buffer (50 mM Tris-HCl, pH 9.2, 16 mM ammonium sulfate, 0.1% Tween 20, and 2.5 mM MgCl2), 0.1 mM dNTPs, 0.2 μM primers, and 0.25 U of ExTaq DNA polymerase (Takara Bio Inc.). Amplification primers including the mutation region within a length of ~1.2 kb were designed using Primer 3 Plus (https://primer3plus.com). PCR was conducted using an iCycler thermal cycler (Bio Rad medical, Waltham, MA, USA). The PCR conditions used for the amplification of all mutation sites were as follows: heat denaturation at 94°C for 2 min; followed by 30 cycles of 94°C for 10 s, 60°C (or the optimal annealing temperature of the different primers) for 10 s, and 72°C for 1 min; and an additional extension for 10 min at 72°C. The amplification products were heated and cooled to form heteroduplexes according to the Surveyor Mutation Detection kit protocol (Transgenomic, Omaha, NE, USA) and incubated at 42°C for 20 min in 30 μL of an enzyme solution containing 10 mM HEPES, pH 7.5, 10 mM MgCl2, 0.002% (w/v) Triton X-100, 0.2 μg/ml of BSA, and 0.4 μL of Surveyor nuclease (Transgenomic).

Subsequently, the digested PCR samples were dissolved in 40 μL of 0.1× TE buffer and electrophoresed using high-efficiency genome scanning (HEGS) electrophoresis with a procedure modified from Kawasaki and Murakami (2000) and Kawasaki . The detection of cleavage patterns was performed using the SYBER Gold Nucleic Acid Gel Stain (Invitrogen, Carlsbad, CA, USA).

SNP-specific primer design, PCR amplification, and real-time quantitative PCR (RT-qPCR)

To introduce MAS into the breeding of sporeless mutants, we developed an allele-specific marker, “Tamogitake-spo”, that specifically amplified mutants using the following DNAs: Y1, LP17, Pc05-1, 108Y2D, and 105 isolates of the test crossing population. The primer design was based on the method of Hayashi , in which 2 mismatched bases against the wild-type sequence were added at the 3ʹ end and adjacent to the 3ʹ end of the forward primer, respectively. This enabled that only the mutant allele would be amplified. The sequences of “Tamogitake-spo” marker were as follows: 5ʹ-TCGATGATGAGCAGAGGATGA-3ʹ (forward) and 5ʹ-TTCGCAGTATGGCTCTAACTGT-3ʹ (reverse). PCR was performed using ExTaq DNA Polymerase (Takara Bio Inc). The reaction mixture usually contained ~100 ng of genomic DNA, 10× ExTaq buffer, 2 μL of 0.1 mM dNTPs, 1 μL of 1.0 μM primers, and 0.25 U of ExTaq DNA polymerase (5 U/μL). The PCR conditions used for the amplification of mutation sites were as follows: heat denaturation at 94°C for 2 min; followed by 30 cycles of 94°C for 10 s, 58°C for 10 s, and 72°C for 1 min; and an additional extension for 10 min at 72°C. The amplified products were analyzed on 1.5% agarose gels in 1× TBE buffer at 100 V for 30 min and visualized using ethidium bromide staining under UV light. RT-qPCR was performed on an iCycler iQTM machine (Bio-Rad Laboratories, Inc., Hercules, CA, USA) in a 20 μL reaction mixture containing 0.5 μL of the above-mentioned extracted sample or genomic DNA (5 ng/μL), 10 μM of each primer, and the THUNDERBIRD SYBR qPCR Mix (TOYOBO, Osaka, Japan), according to the manufacturer’s protocol. The specificity of the amplifications was verified by melting-curve analysis. The PCR procedure was initially started with denaturation at 95°C for 1 min; followed by 40 cycles of 95°C for 15 s, 58°C for 15 s, and 72°C for 30 s. After the amplification was completed, the melting-curve analysis was performed from 60°C to 95°C at a rate of 0.5°C/10 s with continuous fluorescence monitoring. The melting temperature (Tm) value was defined as the peak of the curve. The Tm value was calculated using the iCycler iQTM Real-Time Detection System software, version 3.0 (Bio-Rad Laboratories). Finally, the sporulation phenotype was determined based on the Tm value and Cq value.

Results

NGS and assembly of the Tamogitake Y1 genome

In this study, a total of 38.58 Mb of reads of the Tamogitake Y1 genome sequence data were generated using the PacBio analytical system (Table 2). The results of the total nucleotide composition analysis indicated that the GC content of the genome sequences of Tamogitake Y1 was 49%. After de novo assembly, we identified a total of 88 contigs. The minimum length of the assembled contigs was 1,315 bp and the maximum was 4,706 kb. The average length of the assembled contigs was 438 kb, whereas that of N50 was 2,853 kb (Table 2). The genomic sequence data obtained in this study were deposited in the DDBJ Sequence Database (BHFW01000001-BHFW01000088) and were used as the reference sequence.

Table 2.

Summary of the results of next-generation sequencing (NGS) and de novo assembly of Tamogitake (Pleurotus cornucopiae var. citrinopileatus) Y1 genome sequences

NGS and de novo assembly result summary	Number
Contig numbers	88
Reads sequences (bp)	38,586,867
GC content (%)	49.46
Mean length of contig (bp)	438,487
Shortest contig (bp)	1,315
Longest contig (bp)	4,706,332
N50	2,853,238

Mutations caused by UV-based mutagenesis in the LP17 genome

By mapping the LP17 genome of the 108Y2D mutant to the Y1 reference genome, which was the original genome of LP17 (Materials and Methods), a total of 685 mutation sites were detected in 30 contigs with an average depth of 166 (Table 3). Among them, single-nucleotide substitutions and insertions were the main mutation types, as they accounted for 348 (51%) and 275 (40%) of the mutations, respectively (Table 3).

Table 3.

Mapping LP17 sequence to reference genome Y1 and detection of DNA sequence variation

	value
Reads number	64,748,548
Percent of mappinga	96.64
Average depth	166.71
Total mutation number	685
One base substitution	348
Insert	275
Deletion	62
Contig containing mutation	30

%, Percent of mapping read to total read number.

Linkage analysis via TILLING to select the candidate causative locus of the sporeless trait

To estimate the candidate causative locus of the sporeless mutation, we extracted 11 mutational sites arbitrary that occurred via serious mutation, such as frameshift and nonsense mutations (Supplemental Table 3). Subsequently, we performed a linkage analysis between them and the sporulation phenotype via TILLING analysis using the bulked DNAs from about half of the isolated populations (Yoneyama ), as described in the Materials and Methods. As a result, out of the 11 mutation sites investigated (Supplemental Table 3), only one mutation site (located in contig 12 (2,729,931 bp)) was linked to the sporulation phenotype (Fig. 2). This result suggested that this mutation site on contig 12 is located very close to the target mutation locus. Furthermore, there were no other mutations in the near region of this site (Supplemental Table 4). With expecting the high possibility as be a target locus of this study, we sequenced 300 bp around this mutation site (1,229,199 bp in contig 12) in all of the 105 isolates of the test cross population (Yoneyama ). The results of their sequence type were consistent with the each of sporulation-phenotype (wild-type/sporeless-type) perfectly.

Fig. 2.

Representative electrophoretic gel in TILLING analysis using primer pairs, unitig12-F and unitig12-R (Supplemental Table 3). Numbers: Bulk number (Supplemental Table 2). M: 100 bp ladder marker. Putative mutation in the sporeless bulks (no. 1–6) and mixed bulk (no. 9) could be identified by the presence of two bands (103 bp and 366 bp), with sizes adding up to the full-length PCR product (469 bp). Wild-type bulks showed only the 469 bp band.

Sporeless mutation estimated to be a nonsense mutation in silico

An AUGUSTUS-based analysis was used to estimate the region of the mutation site described above from the deduced CDS sequence (total length of the gene, 3,518 bp; 11 exons separated by 10 introns; encoding a protein of 990 amino acids) (Fig. 3, Supplemental Fig. 1). A blastp search using the encoded amino acid sequence showed high homology with “fungi-specific transcription factors”. In addition, a single-base substitution mutation (C to T) located at position 1,229,199 in contig 12, which was located in the fifth exon region, led to a nonsense mutation (Fig. 3, Supplemental Fig. 1); this corresponded to the amino acid at position 650 (1,950th bp in nucleotide sequences) of the wild-type putative protein (hypothetical protein PC_1229, accession no. LC528953 in DDBJ), which was changed from arginine to a stop codon in the mutant genome. This result did not change when the genetic code of AUGUSTUS was changed from L. bicolor to C. cinerea.

Fig. 3.

Position of the mutation in hypothetical protein PC_1229 and its structure. GAL4 (GAL4-like Zn(II)2Cys6 binuclear cluster DNA-binding domain) and Fungal_TF_MHR (fungal transcription factor regulatory middle homology region) were domains predicted by the blast search. Descriptions based on both coding sequences and amino acid sequences (accession no. LC528953 in DDBJ) are shown in Supplemental Fig. 1.

Validation of a primer pair for the sporeless mutation

Regular PCR was performed individually using these primers, and the results were consistent with the genotype showing the PCR amplification of a band of 172 bp in all mutants (Fig. 4A) and no amplification in the wild-type sequence. Furthermore, RT-qPCR was performed using the same primer pair, “Tamogitake-spo”, which confirmed that the mutant type and the wild type could be distinguished reproducibly using this marker, as shown in Fig. 4B.

Fig. 4.

Validation of the “Tamogitake-spo” marker for the detection of the sporeless trait using genomic DNA of sporeless mutant monokaryon LP17 (1) and dikaryon 108Y2D (3), and wild-type monokaryon Y1 (2) and dikaryon Pc05-1(4), respectively. A. Electrophoretic analysis of samples from A. PCR products were fractionated on a 1.5% agarose gel, followed by visualization using ethidium bromide staining. M: 100 bp ladder. B. Melting curve analysis using PCR products of the “Tamogitake-spo” marker.

To verify the performance of this allele-specific marker in the breeding process, RT-qPCR was performed on culture extracts from mycelia, which were prepared according to Izumitsu . The results obtained for all of the tested samples, which had various genetic backgrounds, such as 110 single-spore isolates from hybrids of breeding materials, and 22 wild isolates were completely consistent with each of the phenotypes (Table 1, Supplemental Table 1).

Discussion

In this study, we successfully performed whole-genome sequencing of Y1 (38.58 Mb) using the PacBio analytical system and carried out resequencing of the genome of its LP17 mutant (Yoneyama , 2017). We also estimated the mutation site that was responsible for the novel sporeless trait of Tamogitake 108Y2D (Yoneyama ). Furthermore, for breeding of the sporeless Tamogitake cultivars, we developed an allele-specific DNA marker based on the sequence of the mutation site in the LP17 genome.

In this study, we assembled 38.58 Mb of the Tamogitake genome from 88 contigs (Table 2). The genomes of two species of the genus Pleurotus, P.ostreatus (PC9 and PC15 strains) and P. eryngii (ATCC90797 strain), to which Tamogitake belongs taxonomically, had been analyzed and published previously (https://mycocosm.jgi.doe.gov). They were assembled to sizes of 34.3, 35.6, and 44.6 Mb, respectively. In addition, a size range of 35–45 Mb is normal for the assembled genome of basidiomycetes, based on published data (https://mycocosm.jgi.doe.gov). Therefore, the size of the Tamogitake genome obtained in this study seems to cover almost the entire genome of this fungus, and the convergence into only 88 contigs afforded by the NGS analysis may ease the identification of the mutated region. In addition, the genome size of mushroom species is smaller than that of crops, and haplophase nuclei exist independently in cells (dikaryotic cell, n + n) corresponding to diploid cells (2n). Furthermore, two constituent nuclei can be isolated separately using the protoplasting technique (Matsumoto ). Thus, the mutated haploid-phase nucleus, LP17, can be compared with its original haplophase nucleus, Y1, to detect the mutation sites. In this study, we mapped 96.64% of the sequence of the mutant LP17 genome against that of the original Y1 genome with an average depth of 166.71 (Table 3). Thus, this study provided an accurate genome-wide comparative analysis.

Based on the genetic analysis of the sporeless mutation targeted in this study, Yoneyama showed that this phenotype is caused by a single dominant gene that controls the mutation generated on the Y1 side, which is the constituent nucleus of the original strain, cultivar Pc05-1. A total of 685 mutation sites, which were distributed 30 contigs, were detected by genome comparison analysis between Y1 and LP17 (Table 3). In this study, to narrow down the site of the sporeless mutation out of the 685 mutation sites, for saving labor and time as much as possible, first of all, the annotation of the genome regions of the mutation sites was performed using the gene prediction software AUGUSTUS, followed by the estimation of each mutation based on the annotation results described above. Subsequently, we selected strong mutation sites, such as frameshift or nonsense mutations, preferentially and examined whether they were linked to the mutation phenotype. A linkage analysis between the sporeless phenotype and the selected mutation sites was performed via the TILLING method using 48 isolates (i.e., phenotypic bulks of four isolates) (Fig. 2, Supplemental Table 2). We identified the contig harboring the mutation that was linked to the sporeless phenotype at only 11 randomly selected locations. The TILLING results were validated by a sequence analysis of the relevant regions in all 105 isolates. Therefore, the experimental scheme indicated in this study, which combined genome-wide comparative analysis via NGS and the execution of TILLING, was considered effective in the absence of gene information on the targeted phenotype for narrowing down the mutation sites.

It is necessary to obtain fruiting bodies by cultivation to detect the presence of the sporeless trait based on the sporulation phenotype, which is time consuming and labor intensive. Therefore, the establishment of a MAS approach that can be determined at the mycelium level is desired. However, in general, it is difficult to design marker primers to differentiate the mutant gene from various wild-type genes based on gene SNPs. In this study, according to the protocol of Hayashi , we designed a forward primer with a 3ʹ end of ‘T’ to match the mutant gene from 108Y2D exclusively, and a mismatch nucleotide (the second nucleotide on the 3ʹ end) was introduced in the primer at the same time. Thus, the two nucleotide mismatches with the DNA harboring only the wild-type gene should not amplify the fragment, while the presence of only one nucleotide mismatch with the DNA harboring the mutant gene results in the amplification of the fragment. The “Tamogitake-spo” marker designed in this study distinguished the wild-type from the mutant-type genes without being affected by differences in the genetic background (Fig. 4, Table 1). Furthermore, amplification and detection using this marker primer were performed successfully, not only in the PCR–agarose-electrophoresis system, but also in the RT-qPCR (without electrophoresis) system (Fig. 4). Therefore, the “Tamogitake-spo” marker based on a SNP derived from the candidate causative gene of the sporeless mutation can be used to detect the presence or absence of the mutation very accurately. In addition, the use of RT-qPCR is expected to increase the efficiency of the selection of this trait greatly, suggesting that it is extremely advantageous for dealing with a large number of isolated generations in the breeding. To breed and utilize a recessive sporeless mutation in Tamogitake (Gisusi ), Okuda previously designed mutation-detection markers. However, the linkage marker reported by them was not mutation-site-sequence specific and not efficient for breeding because of recessive gene. Considering the above results described for the “Tamogitake-spo” marker, it is considered that this marker can be used more easily at the practical site of breeding with much greater precision than that of the marker developed by Okuda ; moreover, its use in the breeding of sporeless Tamogitake in the future is expected. At present, we are promoting the development of sporeless Tamogitake mushrooms carrying this mutant trait by practicing programs using the “Tamogitake-spo” marker. We could have currently bred the new cultivars and we applied to the plant variety protection system on the Ministry of Agriculture, Forestry and Fisheries for registering new cultivars by improving mutant 108Y2D (Fig. 1A) in a few years, since we developed “Tamogitake-spo” marker.

Cytological studies of basidia of this sporeless mutant have shown that impaired post-meiotic sporulation results in an extreme decrease in sporulation (Yoneyama ). As described in the Introduction, sporeless mutants resulting from gene mutations in various different sporulation processes have been isolated in basidiomycetes. Furthermore, several genes related to chromosome behavior, including involvement in sporulation, have been identified using C. cinerea as an experimental model organism for research on the sexual development of fungi (Kanda ). However, to our knowledge, only a few causative genes have been investigated in cultivated mushrooms, for example, dmc1 and rad51 of P. ostreatus (Mikosch , 2002) and stpp1 of P. pulmonarius (Okuda ). The mutation site responsible for the sporeless trait investigated in this study was estimated to result from a SNP located in the coding region of a characteristic transcription factor of fungi (Fig. 3, Supplemental Fig. 1). An in silico analysis performed using the AUGUSTUS software (http://bioinf.uni-greifswald.de/augustus/) showed that the 1,950th base C in the total coding sequence was replaced with a T. We have not determined the 5ʹ and 3ʹ ends of this gene using the RACE method. However, the exon sequences about 1.5 kb upstream and 0.4 kb downstream of the mutation site determined by cDNA sequence analysis were consistent with AUGUSTUS predictions (data not shown). Therefore, as predicted by the in silico analysis described above, the arginine residue at position 650 in the wild type was presumed to generate a stop codon in the mutant via a C < T substitution (Fig. 3, Supplemental Fig. 1). This candidate gene region shares significant homology with the DNA-binding domains found in fungi C6 zinc finger transcription factors, such as the GAL4 protein from yeast. It has been reported that the pro1+ gene encodes a C6 zinc finger transcription factor that required for fruiting body development in Sordaria macrospora (Masloff ). Also, recently, it has been shown that the ada-6 gene encodes a same domain that regulates conidiation and sexual development in Neurospora crassa (Sun ). The reproduction stages in which the ada-6 gene functions in N. crassa are presumed to be similar to the sporulation stage of mushrooms biologically, which is very interesting considering the relationship between mutations in this gene and the sporeless phenotype. In mushrooms, although direct evidence is lacking, our preliminary analysis of the expression activity of this gene by RT-qPCR revealed that it was expressed at the highest level in the gill of the mature fruit body, where basidia are formed (Supplemental Fig. 2). In addition, the expressions of several fungal transcription factors containing a Zn(2)-C6 fungal-type DNA-binding domain and a fungal-specific transcription factor domain was detected in the gills of mature fruiting body of L. edodes (Sakamoto ). Thus, it is considered that these findings described above support the involvement of this gene in sporulation in mushrooms. The fact that the mutation of this gene showed dominant heritability in the current work will be extremely useful for the breeding of the sporeless cultivar of Tamogitake.

Furthermore, if the function of this gene on the sporulation of Tamogitake will be verified by the recombinant DNA experiment, it is expected to lead to the development of sporeless mushrooms based on the mutation of this gene in other mushroom species, via genome editing and other recombination techniques.

The results of the current study demonstrated the success of NGS combined with TILLING for mutation-site identification of a novel sporeless trait (Yoneyama ). In addition, an allele-specific SNP marker was developed for the breeding of the sporeless Tamogitake cultivar, which could be used via both regular PCR detection or RT-qPCR detection. The developed SNP marker could be introduced successfully for the MAS of the sporeless trait in Tamogitake breeding, to verify its usefulness for the detection of this mutation. Furthermore, the candidate causative gene (hypothetical protein PC_1229) of the sporeless mutation estimated in this study has the potential to be used as a target gene in the sporeless breeding of other mushroom species in the future.

Author Contribution Statement

SY designed and executed the study, prepared all tables and figures, and wrote the manuscript. NS, NA, TA and MT performed experiments and analyzed data. TM designed and supervised the experiments, and corrected the manuscript. All authors read the manuscript and agreed to its submission.

Supplemental Figures

Supplemental Tables