Transcriptome Analysis of Light-Regulated Monoterpenes Biosynthesis in Leaves of Mentha canadensis L.

Light is a key environmental aspect that regulates secondary metabolic synthesis. The essential oil produced in mint (Mentha canadensis L.) leaves is used widely in the aromatics industry and in medicine. Under low-light treatment, significant reductions in peltate glandular trichome densities were observed. GC-MS analysis showed dramatically reduced essential oil and menthol contents. Light affected the peltate glandular trichomes’ development and essential oil yield production. However, the underlying mechanisms of this regulation were elusive. To identify the critical genes during light-regulated changes in oil content, following a 24 h darkness treatment and a 24 h recovery light treatment, leaves were collected for transcriptome analysis. A total of 95,579 unigenes were obtained, with an average length of 754 bp. About 56.58% of the unigenes were annotated using four public protein databases: 10,977 differentially expressed genes (DEGs) were found to be involved in the light signaling pathway and monoterpene synthesis pathway. Most of the TPs showed a similar expression pattern: downregulation after darkness treatment and upregulation after the return of light. In addition, the genes involved in the light signal transduction pathway were analyzed. A series of responsive transcription factors (TFs) were identified and could be used in metabolic engineering as an effective strategy for increasing essential oil yields.


Introduction
Light is one of the most important regulators of plant growth and development. It also has a significant impact on the differential accumulation of isoprenoid-derived metabolites [1]. Although significant research has been carried out on light-controlled seedling morphogenesis, seed germination, shade-avoidance, and flower induction [2][3][4][5], much less is known about how light regulates secondary metabolite biosynthesis, and about the regulatory mechanisms of monoterpenes biosynthesis. It has been reported that light intensity can alter terpenoids production via the activation of photosensitive enzymes involved in the methylerythritol phosphate (MEP) pathway [6]. Previous reports have shown that overexpressing Cryptochrome1 increased the amount of artemisinin in the medicinal plant Artemisia annua L. [7]. It has also been reported that the key transcription factor in the light signaling pathway, AaHY5, positively regulates the expression of genes involved in the artemisinin biosynthetic pathway [8].
Mentha L., belonging to the Labiatae family, is widely distributed and cultivated throughout the world. Because of their high economic value, Mentha species have become The organic compounds volatilized from mint leaves were collected via headspace sorption and analyzed via GC-MS. Fifteen volatile compounds were identified, accounting for 99.87% and 97.13% of the total volatile exudates in mint leaves under control and low light conditions, respectively. The changes in total volatile compounds and monoterpene compositions were investigated. The relative contents of volatile compounds were dramatically reduced compared to the internal standard substance, which was bornanone. The level of the main oxygenated monoterpene, menthol (constituting 32.02% of total volatile compounds), decreased to 12.85% following low light treatment. The proportion of pulegone increased from 0.95 to 8.70%, and that of piperitone increased from 1.23 to 28.31% ( Figure 2). The low light condition caused the accumulation of pulegone and piperitone, which are intermediate compounds in mint leaves, but had no reductive effect on the proportions of total monoterpene hydrocarbons or oxygenated monoterpenes (Table 1). The organic compounds volatilized from mint leaves were collected via headspace sorption and analyzed via GC-MS. Fifteen volatile compounds were identified, accounting for 99.87% and 97.13% of the total volatile exudates in mint leaves under control and low light conditions, respectively. The changes in total volatile compounds and monoterpene compositions were investigated. The relative contents of volatile compounds were dramatically reduced compared to the internal standard substance, which was bornanone. The level of the main oxygenated monoterpene, menthol (constituting 32.02% of total volatile compounds), decreased to 12.85% following low light treatment. The proportion of pulegone increased from 0.95 to 8.70%, and that of piperitone increased from 1.23 to 28.31% ( Figure 2). The low light condition caused the accumulation of pulegone and piperitone, which are intermediate compounds in mint leaves, but had no reductive effect on the proportions of total monoterpene hydrocarbons or oxygenated monoterpenes (Table 1).

Illumina HiSeq mRNA Sequencing
To investigate the underlying role of light in menthol synthesis and identify genes with significantly different expression levels under dark conditions, we explored the transcription profiles of mint leaves in the dark for 24 h and light for 24 h, using RNA-seq.

Different Gene Expression Analysis
The expression profiles of Mentha canadensis grown in darkness and light for 24 h each were analyzed separately and were selected based on two criteria: FDR < 0.05 and |log2FC| > 1. Compared to the samples under the control conditions, a total of 10,977 DEGs, including 6296 upregulated and 4681 downregulated genes, were differentially expressed after the 24 h darkness treatment. In total, we counted 8645 DEGs, including 3745 upregulated and 4900 downregulated genes compared to under darkness, after 24 h of light exposure ( Figure 3). All DEGs were subjected to GO, KOG analysis, and KEGG pathway enrichment analysis. Based on the GO enrichment, "metabolic process" and "cellular process", as well as "single-organism process", were the three most highly represented terms in the "biological process" category. In the "cellular component" category, "cell", "cell part" and "organelle" were the most enriched terms. In the "molecular function" category, most DEGs were enriched in the "catalytic activity", "binding", and "transporter activity" subcategories. The enrichment analysis showed that 1919 genes were annotated in 133 KEGG pathways. Among them, 108 genes were involved in the "metabolism of terpenoids and polyketides" pathway ( Figure 4).    The expression profiles of all the DEGs adhered to two response patterns ( Figure  5A  The expression profiles of all the DEGs adhered to two response patterns ( Figure 5A,B). Seven significant expression clusters were observed in all the DEGs under darkness treatment followed by recovery light. The DEGs in cluster 5 were upregulated under darkness and downregulated after recovery (4736 genes), and those in cluster 2 were downregulated under darkness and upregulated under light exposure (3399 genes) ( Figure 5C).
In mint, menthol was the most prevalent chemical compound among the monoterpenoids. In our study, orthologous genes encoding enzymes that catalyzed the biosynthesis of (−)-menthol were identified by multiple sequence alignment, using the reference genes of peppermint and spearmint. The nucleotides of these genes were cloned and sequenced in M. canadensis (Supplementary data). The expression profiles of six identified genes, including limonene synthase (LS), (−)-limonene-3-hydroxylase (L3OH), isopiperitenol dehydrogenase (iPR), isopiperitenone reductase (iPD), pulegone reductase (PR), and neomenthol dehydrogenase (MNR), showed the same patterns: they were significantly downregulated under 24 h of darkness and upregulated after returning to the light for 24 h (Figure 7). Only the MFS gene showed the opposite expression pattern, which was upregulation under darkness and downregulation after removing the covers. Although the GPPS-l, GPPS-s, and MR genes' expression levels were not dramatically altered by light exposure, they were all downregulated as a result of darkness ( Figure 7). In mint, menthol was the most prevalent chemical compound among the monoterpenoids. In our study, orthologous genes encoding enzymes that catalyzed the biosynthesis of (−)-menthol were identified by multiple sequence alignment, using the reference genes of peppermint and spearmint. The nucleotides of these genes were cloned and sequenced in M. canadensis (Supplementary data). The expression profiles of six identified genes, including limonene synthase (LS), (−)-limonene-3-hydroxylase (L3OH), isopiperitenol dehydrogenase (iPR), isopiperitenone reductase (iPD), pulegone reductase (PR), and neomenthol dehydrogenase (MNR), showed the same patterns: they were significantly downregulated under 24 h of darkness and upregulated after returning to the light for 24 h (Figure 7). Only the MFS gene showed the opposite expression pattern, which was upregulation under darkness and downregulation after removing the covers. Although the GPPS-l, GPPS-s, and MR genes' expression levels were not dramatically altered by light exposure, they were all downregulated as a result of darkness ( Figure 7). RT-PCR was used to verify the expression levels of menthol biosynthetic genes (Figure 8). The results indicate that the qRT-PCR expression patterns of 9 genes matched quite well with the RNA-Seq data.  The level of significance, obtained with Student's t-test, is denoted as follows: * p < 0.05; ** p < 0.01. ns, not significant. Each value is the mean of at least three replicates.

Expression of Genes Involved in the Light Signal Transduction Pathway
The gene families of photoreceptors, such as phytochromes, cryptochromes, phototropins, and UVR8, were identified within the DEGs. The expression levels of PHYB (Unigene0004374) and PHOT2 (unigene0003571) were dramatically reduced in darkness and increased slightly when they re-entered light. The CRY gene (Unigene0004413) showed the opposite expression pattern: upregulated under darkness and downregulated after removing the covers. Other photoreceptors' transcription expressions were changed similarly to those of the PHOTs, and these changes were not significant.
COP1, a ring-finger type ubiquitin E3 ligase, functions as a repressor of photomorphogenesis, and was significantly downregulated in darkness and upregulated when returning to light. HY5 and PIFs are the downstream components of photoreceptormediated light signaling. Both of the two HY5 potential homologs, and three potential homologs, were identified. HY5 was significantly downregulated under darkness and upregulated after removing the covers. PIFs are the key regulators in shade-avoidance responses. The transcript levels of PIF1 and PIF4 were dramatically upregulated following the 24 h darkness treatment, while PIF3 was not differentially expressed after the treatment (Figure 9). turning to light. HY5 and PIFs are the downstream components of photoreceptor-mediated light signaling. Both of the two HY5 potential homologs, and three potential homologs, were identified. HY5 was significantly downregulated under darkness and upregulated after removing the covers. PIFs are the key regulators in shade-avoidance responses. The transcript levels of PIF1 and PIF4 were dramatically upregulated following the 24 h darkness treatment, while PIF3 was not differentially expressed after the treatment (Figure 9).

Light-Mediated Regulation of Transcription Factors
Transcription factors play important roles in the response to light, and they regulate multiple genes. Using the PlantTFDB database, 1601 TFs belonging to 57 families were identified. After analyzing the differential expression profiles, totals of 380 and 298 DEGs divided into 44 TF families were differentially expressed. The bHLH, ERF, MYB, and MYB-related, C2H2, bZIP, WRKY, NAC, HD-ZIP, and GRAS families were the top 10 most prevalent TF families. Among them, 235 TFs were upregulated and 145 TFs were downregulated following the 24 h darkness treatment. After removing the covers for 24 h, 166 TFs were upregulated and 132 TFs were downregulated. There were 67 TFs showing similar expression patterns to the terpenoids synthesis pathway genes. The expression levels of these genes were downregulated by darkness and upregulated after returning to light. The MYB family was the largest group, followed by the bHLH family. The opposite expression pattern (upregulated in darkness and downregulated in light) was displayed in 137 TFs. Among them, the ERF family was the largest group (Supplement Figure S4).
Several transcription factors belonging to the bHLH, ERF, MYB, and WRKY families have been determined to control terpenoid biosynthesis. We identified the homologs of the TFs whose functions are known, and which regulate terpenoid biosynthesis, by BLAST

Light-Mediated Regulation of Transcription Factors
Transcription factors play important roles in the response to light, and they regulate multiple genes. Using the PlantTFDB database, 1601 TFs belonging to 57 families were identified. After analyzing the differential expression profiles, totals of 380 and 298 DEGs divided into 44 TF families were differentially expressed. The bHLH, ERF, MYB, and MYB-related, C2H2, bZIP, WRKY, NAC, HD-ZIP, and GRAS families were the top 10 most prevalent TF families. Among them, 235 TFs were upregulated and 145 TFs were downregulated following the 24 h darkness treatment. After removing the covers for 24 h, 166 TFs were upregulated and 132 TFs were downregulated. There were 67 TFs showing similar expression patterns to the terpenoids synthesis pathway genes. The expression levels of these genes were downregulated by darkness and upregulated after returning to light. The MYB family was the largest group, followed by the bHLH family. The opposite expression pattern (upregulated in darkness and downregulated in light) was displayed in 137 TFs. Among them, the ERF family was the largest group (Supplement Figure S4).
Several transcription factors belonging to the bHLH, ERF, MYB, and WRKY families have been determined to control terpenoid biosynthesis. We identified the homologs of the TFs whose functions are known, and which regulate terpenoid biosynthesis, by BLAST searching our transcriptomic data, and identified two MYC1, one MYC2, one WIN1, one ERF1, one WRKY40 and one R2R3-MYB genes. Apart from ERF1, which was significantly upregulated by the darkness treatment, all the genes were downregulated. Interestingly, not all the downregulated genes were upregulated after 24 h light exposure. The transcription levels of the four TFs belonging to the bHLH family and one WIN1 gene were not changed after returning to light for 24 h. The transcription levels of one WIN1, one WRKY40 and one R2R3-MYB were differentially expressed: the WRKY40 was upregulated back to its previous level before the darkness treatment, while the WIN1 and R2R3 MYB remained lower than their initial levels ( Figure 9).

Validation of Expression Patterns of Genes
To verify the reliability of the RNA-Seq data, the expression profiles of eight genes were analyzed using quantitative real-time PCR (qRT-PCR). These genes are listed in Supplemental Table S3. The expression patterns of the selected genes, in terms of the RPKM results, were consistent between the qRT-PCR and RNA-seq ( Figure 10). The following validation results for qRT-PCR indicate that the RNA-Seq data are quite reliable.

Validation of Expression Patterns of Genes
To verify the reliability of the RNA-Seq data, the expression profiles of eight genes were analyzed using quantitative real-time PCR (qRT-PCR). These genes are listed in Supplemental Table S3. The expression patterns of the selected genes, in terms of the RPKM results, were consistent between the qRT-PCR and RNA-seq ( Figure 10). The following validation results for qRT-PCR indicate that the RNA-Seq data are quite reliable.

Discussion
Glandular trichomes, referred to as true cell factories, are able to secrete and store large quantities of specialized metabolites. The monoterpenes of mint essential oil are produced by, and stored in, the peltate glandular trichomes [13]. The number of glandular trichomes per leaf, and their production of chemical compounds, depend on environmental and developmental parameters [14]. The distribution of peppermint gland initiation reflects the basipetal pattern of leaf maturation: relatively immature regions at the leaf's base continue to produce oil glands long after gland production has stopped at the leaf's apex [13]. To accurately determine the densities of mint glandular trichomes under low light conditions, five different areas at the same stage of development, from the basal zone to the top zone of an abaxial leaf, were analyzed. Our results show remarkably reduced densities of peltate glandular trichomes in all five selected leaf zones under low light conditions, which correlates with reduced chemical compound production. The present results suggest that reduced light intensity effects not only peltate glandular trichome density, but also the production and accumulation of essential oil, in Mentha canadensis. This is consistent with the results for Mentra piperita and Mentha arvensis [15,16]. In other plant species, cultivated tomato for instance, increases in photosynthetically active radiation induce the growth of type VI leaf glandular trichomes [17]. In the case of Varronia curassavica

Discussion
Glandular trichomes, referred to as true cell factories, are able to secrete and store large quantities of specialized metabolites. The monoterpenes of mint essential oil are produced by, and stored in, the peltate glandular trichomes [13]. The number of glandular trichomes per leaf, and their production of chemical compounds, depend on environmental and developmental parameters [14]. The distribution of peppermint gland initiation reflects the basipetal pattern of leaf maturation: relatively immature regions at the leaf's base continue to produce oil glands long after gland production has stopped at the leaf's apex [13]. To accurately determine the densities of mint glandular trichomes under low light conditions, five different areas at the same stage of development, from the basal zone to the top zone of an abaxial leaf, were analyzed. Our results show remarkably reduced densities of peltate glandular trichomes in all five selected leaf zones under low light conditions, which correlates with reduced chemical compound production. The present results suggest that reduced light intensity effects not only peltate glandular trichome density, but also the production and accumulation of essential oil, in Mentha canadensis. This is consistent with the results for Mentra piperita and Mentha arvensis [15,16]. In other plant species, cultivated tomato for instance, increases in photosynthetically active radiation induce the growth of type VI leaf glandular trichomes [17]. In the case of Varronia curassavica Jacq., higher irradiance results in an increased essential oil yield and a higher frequency of glandular trichomes [18]. The conclusion is that light density regulates not only the density of mint's glandular trichomes but also the synthesis of essential oil. However, the molecular mechanisms by which light regulates monoterpenes biosynthesis, as well as the reprogramming of the pathway, are still not fully understood. To study the process in more detail, the mRNA expression patterns under changing light conditions were assessed via transcriptome analysis.
Transcriptome sequencing is a powerful tool for screening gene expression patterns and identifying candidate genes. In the current study, high-throughput RNA sequences were generated to analyze the transcriptomes of M. canadensis undergoing a 24 h darkness treatment followed by exposure to light for 24 h. Approximately 72 million high-quality clean reads were generated, and a total of 95,579 unigenes, with an N50 of 1260 bp and average length of 754 bp, were obtained. In a previous study, nine of the orthologous genes involved in the menthol synthesis pathway were identified. Darkness for 24 h down- DEGs involved in the sesquiterpenoid and triterpenoid biosynthesis pathways showed the same transcription patterns, except for unigene0051735 and unigene001095. The significant downregulation of genes involved in the MEP pathway may explain the reduced essential oil metabolite production. Furthermore, changes in the monoterpene profiles show that low light intensity leads to an accumulation of the branchpoint intermediate pulegone and the side product piperitone. In fact, under conditions of salinity and osmotic stresses, pulegone accumulated in Mentha species [19][20][21][22]. This intermediate compound is considered an undesirable component under harsh growth conditions. The downregulation of the gene encoding pulegone reductase (PR), and the consistent expression of the gene encoding (−)-menthol dehydrogenase (MR), may explain the reduced proportions of menthone and menthol. Mahmoud and Croteau [23] reported that high levels of menthofuran reduced pr transcript levels in immature leaves, consequently increasing pulegone content in peppermint. These results suggest that the flux of (+)-pulegone is controlled at the transcriptional level of pr, and they are consistent with the results of our present research. However, our data are presently insufficient to explain the accumulation of piperitone under low light conditions.
The light-dependent transcriptional regulation of MEP pathway genes is conserved in all plant species we investigated. Nevertheless, not much is known about the signal transduction pathways or the transcription factors that convert the environmental information into a transcriptional response [1]. In our transcriptome analysis results, coexpression modules with terpenoid biosynthesis genes from the MEP pathway suggest that the transcription of those genes may be coregulated by the same transcription factor as response to light. In A. annua, HY5, a bZIP transcription factor and a core TF in the light signaling pathway, interacts with the promoter of the monoterpene synthesis gene QH6 [24]. It has also been reported that AaHY5 acts via trichome-specific and light-induced AaGSW1 to control artemisinin formation [8]. Evidence from Arabidopsis has demonstrated that HY5 and PIFs are direct regulators of the light-modulated expression of the DXS and DXR genes encoding the flux-controlling enzymes of the MEP pathway [25]. In the present transcriptome data, homologous genes of HY5 (Unigene0034124), PIF1 (Unigene0086430) and PIF4 (Unigene0065109) were identified in M. Canadensis. The resultant RPKM values indicate that the expression levels of these TF genes were significantly changed under dark conditions. The expression pattern of HY5 was the same as that of the menthol synthesis genes, whereas PIFs showed the opposite expression pattern (upregulated under darkness and downregulated after illumination). These results indicate that a similar regulation mechanism to that targeting the monoterpene synthesis gene (operating via the key-master integrators of light signal transduction HY5 and PIFs) may also operate in M. canadensis. More evidence should be sought in future experiments.
We have established strategies to increase the levels of essential oil in peppermint via the overexpression of unique structural genes in the menthol synthesis pathway, or the overexpression of additional genes that encode enzymes involved in the precursor supply pathway for monoterpene biosynthesis [26]. Transgenic lines were generated, but few of them improved essential oil yields. The overexpression of transcription factors that target and positively control multiple monoterpene biosynthesis, or the suppression of the expression of negatively controlled TFs, could also each represent metabolic engineering strategies for improving essential oil yields. Previously, researchers have elucidated the regulatory roles of several genes encoding the WRKY [27][28][29], MYB [30], bHLH [31], AP2/ERF [32,33], and bZIP [8] classes of TFs in specialized terpenoid biosynthesis. Using RNA-Seq in Solanum lycopersicum, SlMYC1 was identified and shown to transiently transactivate the terpene synthase promoters in Nicotiana benthamiana leaves [34]. Heterologous AtWRKY18, AtWRKY40, and AtMYC2 were overexpressed in S. sclarea hairy roots and proven to regulate the expression of several genes encoding enzymes of the MEP-dependent pathway [35]. The overexpression of SmERF128 increased the expression levels of copalyl diphosphate synthase 1 (SmCPS1), kaurene synthase-like 1 (SmKSL1) and cytochrome P450 monooxygenase [32]. The functional study of TFs showed that the regulated monoterpene genes in mint were inadequate. The TFs MsYABBY5 and MsMYB were reported to negatively regulate the production of monoterpenes in spearmint [11,12]. In our transcriptome data, Unigene0018692 was homologous to MsMYB, and was found to be dramatically differently expressed under dark conditions. Interestingly, given that it was regulated by light, the expression pattern of this gene was similar to that of the negatively regulated genes, which were downregulated by dark and upregulated by light (mean RPKM value of 24.78 in control vs. 0.76 under darkness and 5.69 after illumination). These results suggest that other TFs must be more essential in regulating monoterpene gene expression through light signal transduction pathways. Our present work provides a list of candidate genes for further metabolic engineering research and will be helpful in exploring the molecular mechanism of the transcriptional responses of MEP pathway genes to light signals. However, the present RNA-seq results could not determine whether splice variants could vary between the absence and presence light. Full-length transcriptome analyses and more verification research should be carried out in the future.
In conclusion, the results of this study prove the effects of light on the density of peltate glandular trichomes and essential oil production in cornmint. The low light condition reduced the PGT densities in all five different zones, from the bottom to the top of the mint leaves, by half a percent. We observed dramatically decreased total essential oil and menthol levels, and increased proportions of pulegone (from 0.95 to 8.70%) and piperitone (from 1.23 to 28.31%). Downregulated expression levels of menthol synthesis pathway genes McGpps-L, McGpps-s, McLs, McL3OH, MciPD, MciPR, and McPR and the monoterpenes biosynthesis genes, as well as the sesquiterpenoids biosynthesis pathway genes, were observed following darkness treatment, and may explain the reduced levels of essential oil produced via the regulatory mechanisms at the transcription level. The identified DEGs in the light signal transduction pathway, and the light-mediated regulation of transcription factors, could be used as effective strategies to achieve increased essential oil yields in future metabolic engineering studies.

Plant Materials and Treatment
For phenotype difference experiment, mint (Mentha canadensis) plants were grown on the field in Institute of Botany, Jiangsu Province and Chinese Academy of Sciences, Nanjing, Jiangsu Province, China (28 • 05 N, 113 • 37 E). Low-light treatment group was performed by covering plots with sunshade nets at a height of 2.5m above the ground, which caused a reduction of 50% light. The control group grown at normal Sunshine. Two groups both grown in the same soil with the same 12 h photoperiod. In total, three replications were made for the control and the low light groups, respectively, each replication containing 10 plants. After 30 d low light treatment, leaves were harvest to statistic of peltate glandular trichomes densities and analysis essential oil composition.
For transcriptome analysis, mint (Mentha Canadensis) were propagated from rhizomes and were planted in plastic pots containing a mixture of soil and vermiculite (3:1, v/v). Plants were cultured in growth chamber under control condition at 1000 µmol·m −2 ·s −1 , with the 14 h photoperiod and temperature maintained at 25 • C/22 • C. For darkness treatment, black paper bags was used to cover the whole plants for 24h, then move the bags and exposure the plants to light 24 h. Leaves at the same developing stage were collection at 24h darkness treatment time point and at 24 h illumination time point. Control materials were collected at the same time point as darkness treatment. Three plants were selected as biological replicates. Samples were frozen in liquid nitrogen and stored at −80 • C before RNA extraction.

Statistics of Peltate Glandular Trichomes Densities and Measurement of Essential Oil Composition in Mint Leaves under Darkness Treatment
Leaves at the same development stage were divided into five sampling zones from bottom to top. The number of peltate glandular trichomes in each 1 mm 2 were counted after taking a photograph using fluorescence stereomicroscope MVX 10. The analysis was performed in fifth repetition.
Sampling of volatiles from 0.5 g dry weight leaves were collected by headspace solid-phase microextraction (HS-SPME), followed by gas chromatography/mass spectrometry (GC/MS) determination. The polydimethylsiloxane (PDMS) 65 µm fiber (65 µm Carboxen™-PDMS StableFlex, Supelco) was then exposed to the headspace of the sample for 3 min at 60 • C to adsorption. Then volatiles absorbed were desorbed in the injection port of GC at 250 • C for 3 min. Each sample was spiked with 300 ng bornanone internal standard. The analysis was performed in triplicate.
GC analyses were carried out using a Trace 1300 gas chromatograph equipped with a flame ionization detector (FID) and an electronic pressure control (EPC) injector. GC setting were as follow: inlet 250 • C, the carrier gas was N2 (U) with a split flow 50.0 mL/min and the split ratio was 33:1. Analyses were performed using the following temperature program: oven temperature isotherm at 60 • C for 1 min, from 60 • C to 220 • C at the rate of 3 • C min −1 .
GC-MS analysis was performed on a gas chromatograph Trace 1300HP interfaced with a HP 5972 mass spectrometer with electron impact ionization (70 eV). The column temperature was programmed as follows: initial temperature was 50 • C, increased to 240 • C at 5 • C min-1, and isotherm at 205 • C during 10 min. Helium was used as the carrier gas with a flow rate of 1.0 mL/min −1 . The split ratio was 33:1. Scan time and mass range were 0.2 s and 33-450 m/z, respectively.
Identification of volatile compounds was based on the comparison of their mass spectral fragmentation patterns with the National Institute of Standards and Technology (NIST) mass-spectral library data and other published mass spectra (Adams 2001). The percentage content of compounds was calculated by the area normalization method.

Total RNA Extraction, Library Construction and Transcriptome Sequencing
Total RNA from nine samples including three repetition samples of control, 24 h darkness treated and 24 h exposure to light leave materials were isolated using RNA iso Plus reagent Kit (Takara Bio, Dalian, China). First strand cDNA was synthesized using M-MLV Reverse Transcriptase (Promega, Madison, WI, USA). The cDNA library was constructed and transcriptome sequencing was performed by ORI-GENE Company (Beijing, China). The transcriptome libraries were sequenced using an Illumina HiSeq™ 4000 (Illumina, Inc., San Diego, CA, USA) at the Gene De novo Institute (Guangzhou, China).

De Novo Assembly and Functional Annotation
After removing the adaptor and low-quality sequences, a de novo assembly of the remaining high-quality clean reads was performed using Trinity, as previously described for use without a reference genome [36]. After assembly, unigenes were obtained and functional annotation was carried out using BLASTx program (Available online: http:// www.ncbi.nlm.nih.gov/BLAST/) (accessed on 30 April 2021) with an E-value threshold of 1e−5. Public databases including the NCBI non-redundant protein database (Nr) (Available online: http://www.ncbi.nlm.nih.gov) (accessed on 30 April 2021), the SwissProt database (Available online: http://www.expasy.ch/sprot) (accessed on 30 April 2021), the Kyoto Encyclopedia of Genes and Genomes database (KEGG) (Available online: http://www. genome.jp/kegg) (accessed on 30 April 2021), and the COG/KOG database (Available online: http://www.ncbi.nlm.nih.gov/COG) (accessed on 30 April 2021), were used to annotate the M. canadensis unigenes. TFs were identified by aligning unigenes to the PlantTFDB database [37].

Identification and Analysis of Differentially Expressed Genes
The expression levels of the unigenes were normalized to RPKM [38]. The DEGs were defined as having an absolute fold-change value ≥2 and a false discovery rate (FDR) < 0.05 and |log2FC| > 1. Heatmaps of DEGs were generated using HemI [39].

Real-Time Quantitative PCR Analysis
qRT-PCR was performed to verify the expression of selected genes. First-strand cDNA was synthesized with oligo (dT)18 and M-MLV reverse transcriptase (Promega). qRT-PCR analysis was carried out using the SYBR Universal qPCR Kit (Vazyme) on a qTOWER2.2 Real Time PCR Systems (Analytik, Jena, Germany), according to the manufacturer's instructions, and subjected to qPCR as described in Qi, et al. [40]. Quantification was performed using the 2 −∆∆CT method, and data were normalized to those of the actin gene transcript. Sequences of primers used are listed in additional file Supplementary  Table S3. RT-PCR analysis was conducted with three technical replicates, and the data represent the means ± standard errors (n = 3).
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/plants10050930/s1, Figure S1: GO functional classification of assembled unigenes. Figure  S2: KOG functional classification of assembled unigenes. Figure S3: KEGG pathway annotation of assembled unigenes. Figure S4: Distribution of transcription factor (TF) families. Figure S5: GO terms of DEGs among control light condition vs darkness. Figure S6: GO terms of DEGs among darkness condition vs exposure to light. Table S1: Summary of Illumina HiSeq4000 transcriptome sequencing and de novo assembly results.