Gs2019.10.465 3of14 2.Materials and Methods 2.1.Strains and Culture Condition All samples used for whole- ome bisulfite sequencing (WGBS)and RNA-se g were collected from identical batches of cultivated Pleurotus eryngii var.eryngii (Pe,strain ID:JKXB130DA)and Pleurotus tuoliensis (Pt,strain ID:JKBL130LB)at the three major developmental stages:mycelium, primordium and fruit body.Both monokaryotic mycelia were cultivated in PDA (Potato De extrose Agar)liquid medium in dark at 23 days witl naking culture(120 rpm).Mycelium fi by sterile gauze was read extraction.In or um,first v.For Pt cold siml k d ay. formed at 0C for 48 h.Th d P vere under hle light condition (3001x-10001y)at 14C for imordium initiation aftor 10 davs of cultivation.primordia with 2-3 cm were sampled from the bottles.finally.the remaining cultures were transferred to a light(15C)/dark(7C)photoperiod of 12h condition for about 15 days to induce fruit body formation.For RNA-seq,three biological replicates per condition were separately sampled.All samples were stored in liquid nitrogen until DNA and RNA extractions.Morphology of all sequencing samples were shown in Figure Sl 2.2.Whole-Genome Bisulfite Sequencing and Data Processing Genomic DNA (NA)from primordium and fruit body of P sis and P was extracte B.CA eared by son to th ng first Trimn nomatic [32]wa uality ads The the clean data were aligned to the corresr ponding draft genomes of Pt and Pe with 1-bp mismatch by Bismark program [33].Only uniquely mapped reads were used for further analysis. 2.3.Differential Methylation Analysis Cytosine sites in CG-contexts with more than four uniquely mapped reads were used to further analysis.To decide methylated site,bisulfite non-conversion rate of%evaluated from unmethylated ADNA was used as background methylation leve .Then,bin mial test was performed for each cyto ide if the observed methylation level signifi ntly higher than such ackground methylatio . above adjusted e metho ned a ger tha roach with lkb size with step size of 200-bp was used to identify DMRs.Only windows including at least 10 mCG in at least one tissue were for further DMRs identification.For each window,Fisher's exact test was performed to measure the difference of methylation level (reads count as agent)betweer adjacent developmental stages.P-values of above test were adjusted by Benjamini-Ho chberg methoc ues Windows with 0.01 and the changes of methylation level 0.15 were defined 5 rlap ng window with identical DMRs characteristics (hyper-or hypo-methylation) were am ikb of g enes)and g e bodies including at least 10 mcq/kb ir t loast ono (MP)and ne body (MGB) respectively,and were used for further identification of different methylation level described as above 2.4.RNA-seq and Data Process mRNA extracted from mycelium,primordium and fruit body of P.tuoliensis and P ervngii were sequenced by the Illumina Hiseq 2500 platform.Raw data were filtered by removing low-quality reads. Genes 2019, 10, 465 3 of 14 2. Materials and Methods 2.1. Strains and Culture Conditions All samples used for whole-genome bisulfite sequencing (WGBS) and RNA-seq were collected from identical batches of cultivated Pleurotus eryngii var. eryngii (Pe, strain ID: JKXB130DA) and Pleurotus tuoliensis (Pt, strain ID: JKBL130LB) at the three major developmental stages: mycelium, primordium and fruit body. Both monokaryotic mycelia were cultivated in PDA (Potato Dextrose Agar) liquid medium in dark at 23 ◦C for 7 days with shaking culture (120 rpm). Mycelium filtered by sterile gauze was ready for nucleic acid extraction. In order to acquire primordium, first, liquid cultures were transformed to cultivation bottles and grown in dark at 25 ◦C for 60 days with 70% humidity. For Pt, cold simulation was performed at 0 ◦C for 48 h. Then, cultivation bottles of Pt and Pe were under blue light condition (300lx–1000lx) at 14 ◦C for primordium initiation. After 10 days of cultivation, primordia with 2–3 cm were sampled from the bottles. Finally, the remaining cultures were transferred to a light (15 ◦C)/dark (7 ◦C) photoperiod of 12 h condition for about 15 days to induce fruit body formation. For RNA-seq, three biological replicates per condition were separately sampled. All samples were stored in liquid nitrogen until DNA and RNA extractions. Morphology of all sequencing samples were shown in Figure S1. 2.2. Whole-Genome Bisulfite Sequencing and Data Processing Genomic DNA (gDNA) from mycelium, primordium and fruit body of P. tuoliensis and P. eryngii was extracted by CTAB method and sheared by sonication to the size range of 200 to 300 bp. Library construction and Bisulfite treatment were processed as described in Hu et al. [31]. Sequencing was performed on the Illumina Hiseq 2500 platform (Illumina; San Diego, USA) with standard protocols. For data processing, first, Trimmomatic [32] was used to remove low-quality sequencing reads. Then, the clean data were aligned to the corresponding draft genomes of Pt and Pe with 1-bp mismatch by Bismark program [33]. Only uniquely mapped reads were used for further analysis. 2.3. Differential Methylation Analysis Cytosine sites in CG-contexts with more than four uniquely mapped reads were used to further analysis. To decide methylated site, bisulfite non-conversion rate of 0.3% evaluated from unmethylated λDNA was used as background methylation level. Then, binomial test was performed for each cytosine site to decide if the observed methylation level significantly higher than such background methylation level [34,35]. P-values of the above tests were then adjusted by Benjamini-Hochberg method to q-values. Cytosine sites with q-value lower than 0.01 and corresponding methylation level larger than 5% were defined as methylated CG context (mCG). Differentially methylated regions (DMRs) between adjacent developmental stages were identified as described in [14] with minor modification. A sliding-window approach with 1kb size with step size of 200-bp was used to identify DMRs. Only windows including at least 10 mCG in at least one tissue were for further DMRs identification. For each window, Fisher’s exact test was performed to measure the difference of methylation level (reads count as agent) between adjacent developmental stages. P-values of above test were adjusted by Benjamini-Hochberg method to q-values. Windows with a q-value < 0.01 and the changes of methylation level ≥ 0.15 were defined as DMRs. Overlapping window with identical DMRs characteristics (hyper- or hypo-methylation) were merged into a large region. Similarly, promoters (upstream 1kb of genes) and gene bodies including at least 10 mCG/kb in at least one tissue were defined as methylated promoter (MP) and methylated gene body (MGB), respectively, and were used for further identification of different methylation level described as above. 2.4. RNA-seq and Data Process mRNA extracted from mycelium, primordium and fruit body of P. tuoliensis and P. eryngii were sequenced by the Illumina Hiseq 2500 platform. Raw data were filtered by removing low-quality reads