Cuili Liang1,2,3†Yitong Yan1,2,3†Yingchao Tan1,2,3Xue Yang1,2,3Jie Cao1,2,3
Chaorong Tang1,2,3,4*
Kaiye Liu1,2,3*- 1National Key Laboratory for Biological Breeding of Tropical Crops, Hainan University, Haikou, China
- 2College of Tropical Crops, Sanya Nanfan Research Institute of Hainan University, Hainan Yazhou Bay Seed Laboratory, Hainan University, Haikou, China
- 3Natural Rubber Cooperative Innovation Center of Hainan Province and Ministry of Education of PR China, Hainan University, Haikou, China
- 4Yunnan Institute of Tropical Crops, Xishuangbanna, China
MicroRNAs (miRNAs) are widely involved in various aspects of plant growth and development. However, how miRNAs and their targets regulate natural rubber metabolism remains unclear in the rubber-producing dandelions, which are being developed as alternative commercial sources of natural rubber. Here, we combined small RNA sequencing, degradome sequencing, target gene prediction, and mRNA sequencing to identify miRNAs and their targets in two dandelion species, the high rubber-yielding Taraxacum kok-saghyz (Tk) and the low rubber-yielding T. spadiceum (Ts). A total of 142 miRNAs, including 108 known and 34 novel ones, were discovered, with 53 identified as differentially expressed (DE) between the latex of Tk and Ts. Degradome sequencing identified 145 targets corresponding to 74 miRNAs. TAPIR and psRNATarget, respectively, predicted 165 and 164 non-redundant targets for the 53 aforementioned DE miRNAs. Gene ontology (GO) enrichment analysis indicated the DE miRNAs and their targets might affect natural rubber production via regulating macromolecular biosynthesis and metabolism in latex. Four critical types of regulatory modules, including miR172-AP2/ERF, miR164-NAC, miR160-ARF, and miRN19-protein kinase, were identified and their interaction networks were constructed, indicating a potential involvement in natural rubber production. The findings and the large miRNA dataset presented here are beneficial to further deciphering the roles of miRNAs in the biosynthesis of natural rubber and medicinal metabolites in dandelion.
1 Introduction
Natural rubber (cis-1,4-polyisoprene, NR) is a strategically important raw material used in the manufacturing of medical, agricultural and industrial products (Cornish, 2001). It is a high molecular polymer with collective elite physical and mechanical properties, such as elasticity, abrasion and impact resistance, heat dispersion and malleability at cold temperature, making it irreplaceable by synthetic alternatives in most applications (Puskas et al., 2014). According to International Rubber Study Group’s (IRSG) report, the global NR demand and production surged by 9.4% and 5.7% in 2021, respectively. Total NR production is forecasted to further grow by 3.7%, reaching 14.80 million tons in 2023. Presently, NR is commercially obtained exclusively from a single tropical tree species: the Para rubber tree (Hevea brasiliensis). However, NR supplies from Hevea tree are not sustainable due to its restriction to specific tropical regions, susceptibility to fungal infections (especially the South American Leaf Blight, SALB), and laborious and skilled harvest work (Mooibroek and Cornish, 2000; van Beilen and Poirier, 2007a; Ahrends et al., 2015; Men et al., 2019). Therefore, it is urgent to develop alternative rubber crops. So far, more than 1800 plant species have been identified as containing rubber in their latex (Metcalfe, 1967). Nevertheless, only a few of these species, including Parthenium argentatum Gray, two dandelion species of Taraxacum kok-saghyz (Tk) and T. brevicorniculatum (Tb), Lactuca sativa, and Ficus bengalensis, produce NR with an average molecular weight over 1000 kg/mol, an essential determinant of high rubber quality (Rousset et al., 2021; Salehi et al., 2021).
The genus Taraxacum, commonly known as the dandelion plants, is a member of the family Asteraceae, subfamily Cichorioideae, and is widely distributed in the temperate zones (Schütz et al., 2006). Taraxacum has long been used as medicinal herbs in traditional Chinese medicine to treat hepatitis and the immune response to upper respiratory infections, bronchitis, and pneumonia (Sweeney et al., 2005; Qadir et al., 2022). In Germany, the application of Taraxacum is recorded in medicating gout, diarrhea, blister, spleen, and liver complaints (Faber, 1958). In North America, the infusions and decoctions of the Taraxacum are applied to treat kidney disease and heartburn (Sweeney et al., 2005). In Mexico and Turkey, the dandelion is used as a laxative and potent anti-diabetic medicine (Schütz et al., 2006). The therapeutic actions of Taraxacum species have been partially ascribed to some sesquiterpenes in their roots and leaves. Recently, the extracts of potential pharmaceutical importance, including a number of sesquiterpenes, triterpenes, phytosterols, and phenolic compounds from dandelion roots, were identified in dandelion plants (Leu et al., 2003; Kisiel and Michalska, 2006). About seventy Taraxacum microspecies are identified in China, of which Tk and Tb are known to produce high-quality NR in their latex, the milky cytoplasm of specialized cells known as laticifers (Post et al., 2012). Tk, also called Russian dandelion or Rubber dandelion, was exploited as the rubber-producing crop in the USSR due to it significant accumulation of high-quality rubber in the root (5% to 24% on a dry weight basis) (van Beilen and Poirier, 2007b). In its roots, the laticifer cells, appearing as long tubular vessels produce and store latex, similar to the laticifers in the rubber tree’s bark (Ramirez-Cadavid et al., 2017). Meanwhile, its known pathway of rubber biosynthesis also corresponds roughly with that of the rubber tree (Tang et al., 2016; Lin et al., 2017). Proteins identified as important for rubber biosynthesis in the rubber tree, including CPT (cis-prenyltransferase), CPTL (cis-prenyltransferase-like), and REF/SRPP (rubber elongation factor/small rubber particle protein), are also recognized as crucial components of the rubber-producing machinery in dandelion species (Salehi et al., 2021). The advantages of Tk as perennial herb, including a relatively simple genome, wide planting area, small plant architecture, ease to transform and harvest, and a relatively short life cycle, make it an ideal model plant for studying natural rubber biosynthesis and production (Lin et al., 2017; Lin et al., 2021). Hitherto, we know little about the regulatory mechanisms of the biosynthesis of secondary metabolites including the high molecular rubber and the component of therapeutic agents in Taraxacum.
miRNAs are endogenous noncoding RNAs with 20-24 nucleotides (nt) that play important roles in regulating plant growth and development as well as biotic and abiotic stress (Song et al., 2019). They are produced through a multistep process including transcription of miRNA genes, precursor processing, and assembly of miRNA-induced silencing complex (miRISC). The mechanisms of miRNA in regulating target gene expression are conserved across different plant species, including cleavage of target genes by base pairing, translational repression, and miRNA-dependent DNA methylation (Jones-Rhoades et al., 2006). The crucial roles of miRNAs have been reported in the regulation of plant growth and development (Chen, 2009). For instance, in Arabidopsis, miR156, miR172, and miR159 are involved in the regulation of vegetative phase change, floral transition, and flowering time (Aukerman and Sakai, 2003; Millar and Gubler, 2005; Schwab et al., 2005; Wu and Poethig, 2006). Additionally, miR160 and miR164 play a role in regulating root development and the emergence of lateral roots (Guo et al., 2005; Mallory et al., 2005). miRNAs are also involved in the responses of plants to environmental stimulation as well as biotic and abiotic stresses (Song et al., 2019). For example, miR396 inhibits leaf growth by targeting GRFs (Growth Regulating Factor) under UV-B radiation damage in Arabidopsis (Casadevall et al., 2013). miR393 targets the auxin receptor genes HvTIR1 and HvAFB, contributing to the inhibition of root elongation under aluminum stress in barley (Bai et al., 2017). Overexpression of miR319 enhances cold tolerance of rice and sugarcane by down-regulating the TCP transcription factors (Thiebaut et al., 2012; Yang et al., 2013). Down-regulation of miR165/166 confers enhanced drought resistance via the elevation of ABA levels in Arabidopsis (Yan et al., 2016). Moreover, miRNAs act as a key regulator in the production of secondary metabolites, which are predicted to protect plants from a series of environmental conditions (Zhang et al., 2022). miR5021 and miR414 regulate the biosynthesis of terpenoid backbone, sesquiterpenoid, and triterpenoid in Xanthium strumarium (Fan et al., 2015). miR2161 regulates the biosynthesis of benzylisoquinoline alkaloids in opium poppy plants by targeting mRNA encoding enzymes such as S-adenosyl-L-methionine, 30-hydroxy-N-methylcoclaurine, and 40-O-methyltransferase (Gupta et al., 2017). Recently, a total of 574 and 396 miRNAs were predicted, respectively, in Tk and T. officinale by homologous retrieval using the publicly available RNA-seq and EST data of these two Taraxacum species (Karimi et al., 2022). Of these miRNAs, miR5021, miR838 and miR1533 were speculated by the KEGG analysis to participate in the terpene biosynthesis (Karimi et al., 2022). So far, the miRNAs and their targets have never been investigated by targeted sequencing analysis in Taraxacum and the roles of miRNA-target modules in the regulation of rubber biosynthesis are not known.
Here, to identify the key miRNA-target modules in the natural rubber biosynthesis pathway of Taraxacum plants, we performed small RNA sequencing, degradome sequencing, target gene prediction, and mRNA sequencing in a pair of Taraxacum species with contrasting ability of rubber production, the elite Tk and the inferior Ts (T. spadiceum). A vast number of miRNAs were identified in the two dandelions species, and most importantly 53 ones revealed to be differentially expressed between the latex of the two dandelion species. The targets of these differentially expressed miRNAs were significantly enriched in the processes of macromolecule biosynthesis and metabolism in latex. The miRNA-target regulatory modules that implicate in the biosynthesis of natural rubber were further investigated. Together, the results we obtained provide the first-hand information on the role of miRNAs in natural rubber biosynthesis and a valuable basis for exploring the detailed functions of miRNAs in rubber-producing dandelions.
2 Materials and methods
2.1 Plant material and growth conditions
Tk and Ts germplasm populations were collected in Xinjiang, China, and the accessions of Ts-01 and Tk-20 from the two microspecies were selected for the study. Tk and Ts plants were propagated by tissue culture following the procedures: the leaves were cut into small pieces of 0.5 cm2 in liquid Murashige & Skoog (MS) medium (4.4 g/L MS, 20 g/L sucrose, pH 5.8) and put on the propagation medium (4.4 g/L MS, 20 g/L sucrose, 0.7 mg/L kinetin, 0.2 mg/L IAA, 9 g/L Agar, pH 5.8) for 3 weeks; the regenerated buds were transferred to the rooting medium (2.2 g/L MS, 30 g/L sucrose, 9 g/L Agar, pH 5.8) for 4 weeks. The resulting plantlets were transplanted to soil and grown in a growth chamber at 16 h light/8 h dark, 24°C, and 60% humidity for 8 weeks. Afterwards, these plants were transferred to a 6°C growth chamber for 4 weeks of vernalization, and then cultivated in the growth chamber at 16 h light/8 h dark, 24°C, and 60% humidity until flowering.
2.2 Histochemical staining and rubber content determination
The roots from five-month-old Tk and Ts plants were fixed in FAA solution for 24 hours, and then sectioned into 6 μm slices. The slices were stained with saturated oil red O solution for 10 minutes, differentiated in 60% isopropanol for 10 seconds, washed with double distilled water and then sealed with glycerin gelatin. The bright field observations were performed using a microscope (SZX16, Olympus, Japan). One gram of dried roots from the Tk and Ts plants was ground into a fine powder and boiled in 20 ml of double-distilled water for an hour. The supernatant was discarded, and 20 ml of 3% KOH was added and boiled for another 2 hours. The coagulated rubber was taken out and dried to a constant weight in a hot air oven at 60°C, and then weighted.
2.3 RNA extraction, library construction and sequencing
During the full flowering stage, total RNA was extracted from the flower, leaf, root and latex of Tk and Ts using the TRNzol Universal Reagent (Catalog Number 4992730) from TIANGEN Biotech Co., Ltd. (Beijing, China) according to the manufacturer’s instructions. The quality and integrity of RNA were evaluated using the INFINTE 200 PRO (TECAN, Switzerland) and Agilent 2100 Bioanalyzer (Agilent Technologies, CA, USA). In total, twelve small RNA libraries were constructed for the both Taraxacum species, including one for each of the tissues of flower, leaf and root, and three for the tissue of latex, and subjected to small-RNA sequencing. Briefly, the libraries were constructed with approximately 2.5 µg total RNA per sample by using the NEBNext® Multiplex Small RNA Library Prep Kit for Illumina (NEB, MA, USA). The small RNA was ligated with 3′ and 5′ end adapters, respectively, followed by reverse transcription to synthesize cDNA and sequenced on the Illumina HiSeq 2500 platform. Six mRNA libraries from latex were constructed by using the NEBNext Ultra RNA Library Prep Kit following the manufacturer’s instructions. mRNA libraries were sequenced by the Illumina NovaSeq™ 6000 platform to generate 150-bp paired-end reads. About 100 µg total RNA from the tissues of flower, leaf, root and latex of each dandelion species were equally pooled for degradome library construction. Degradome cDNA library was also processed on the Illumina HiSeq 2500 platform.
2.4 Identification of miRNA and differentially expressed miRNA
We processed the raw reads of miRNA sequencing using the following steps: the data quality was evaluated by fastp v0.210 (Chen et al., 2018); the adapter sequences, low-quality reads, and inserted reads of less than 18 nt were removed by Cutadapt v2.10 (Martin, 2011); the structural and organelle RNAs, including tRNA, rRNA, snoRNA, and snRNA, were removed by aligning the reads to the Rfam database (https://rfam.org/) and organelle RNA database. To identify known and novel miRNAs, the retained reads were aligned to the miRBase V22 (https://www.mirbase.org) and the reference genome (Unpublished data), respectively. The reads were aligned to the reference genome by using bowtie v1.3.1 with the following parameters: –all -m 20 –best –strata -v 1 (Langmead et al., 2009). The novel miRNAs were further identified by sRNAminer v1.1.1 (https://github.com/kli28/sRNAminer) using the merged data of all the small RNA sequencing libraries (Axtell and Meyers, 2018; Feng et al., 2019). The miRNA expression level was calculated and normalized to transcripts per 10 million (TP10M), and the differentially expressed miRNAs were identified and screened by DESeq2 with log2|FC| > 1, padj < 0.01 and total read counts >100 of six latex libraries (Love et al., 2014).
2.5 Target prediction, degradome, and mRNA sequencing analysis
The targets of differentially expressed miRNA were predicted by TAPIR (score <= 4; free energy ratio >= 0.7) and psRNATarget (expectation <= 2; other parameters were set to default), and only the target genes that were jointly identified by two software programs were used for further analysis (Bonnet et al., 2010; Dai et al., 2018). Degradome sequencing reads were analyzed by the CleaveLand v4.5 pipeline using the parameters: mode 1 and a p-value of less than 0.05 (Addo-Quaye et al., 2009). miRNA targets of different categories were represented in target-plots (T-plots). The mRNA sequencing data were analyzed following the HISAT2 v2.2.1/StringTie v2.1.5 pipeline (Pertea et al., 2016). The gene expression level was normalized by Fragments Per Kilobase of exon model per Million mapped fragments (FPKM) and the differentially expressed genes were determined by DESeq2, with the padj < 0.01.
2.6 GO enrichment analysis
The GO analysis of miRNA targets was performed using the OmicShare tools, a free online platform for data analysis (https://www.omicshare.com/tools).
2.7 miRNA-mRNA regulatory networks construction
The regulatory networks of miRN19, miR160a, miR164a and miR172a, and their targets were constructed and visualized by Cytoscape V3.9.1 (Shannon et al., 2003).
3 Results
3.1 Comparison of rubber production between Tk and Ts
Tk and Ts plants were propagated from leaves by tissue culture, vernalized and grown in pots in the growth chamber. The rubber production phenotypes were measured at the peak flowering time, which was about two months after vernalization (Figure 1A). Compared to Ts, much more latex outflowed in Tk after cutting the main root. Oil-red O staining of transverse root sections showed both Tk and Ts contain ring-shaped and well-dispersed young and mature laticifer cells (Figures 1B, C). Nevertheless, the intensity of oil-red staining revealed that the number of rubber particles in the laticifer cells of Tk was significantly higher than that of Ts (Figures 1D, E). When we broke and pulled apart the dried roots, apparent rubber filaments were observed in Tk but not in Ts, indicating a much higher rubber content of Tk than that of Ts (Figure 1F). Moreover, the rubber content of Tk was quantitatively determined by dry root weight to be 14.07%, whereas that of Ts was barely 0.03%, confirming a striking difference in rubber content between Tk and Ts (Figure 1G).
Figure 1 Phenotypic comparison between Tk and Ts in the flowering period. The representative plants of Tk-20 and Ts-01 in the flowering period (A). In the insets, the outflow of latex in the Tk and Ts was shown after cutting the roots. Oil red O staining was conducted on root cross section of Tk (B, D) and Ts (C, E). Rubber filaments were compared between the dried roots of Tk and Ts (F). Natural rubber contents were measured between Tk and Ts (G). ** indicates p < 0.01.
3.2 Small RNA sequencing in Tk and Ts
In order to achieve more integrated miRNA information of the genus Taraxacum, a total of 12 small RNA (sRNA) libraries were constructed and sequenced from the tissues of flower (FW), leaf (LF), root (RT) and latex (LX) in Tk and Ts. Small RNA-Seq generated raw reads ranging from 9,907,779 to 18,806,129 per library. The raw sequences were computationally analyzed to remove low quality sequences, adaptors, and reads shorter than 18 nt. After filtering, 3,520,292 to 8,209,888 clean reads per library, representing 28%-77% of the initial raw reads, were retained. The sRNA lengths varied extensively from 18 to 30 nt, with the majority ranging from 21 to 24 nt (Figure 2).
Figure 2 Length distribution and abundance of small RNAs in flower (FW), leaf (LF), root (RT) and latex (LX) of Tk and Ts. The sequential length distributions of clean reads ranging from 18 to 30 nt of Tk (A) and Ts (B) are shown. The red and blue bars indicate the lengths of 21 and 24 nt, respectively.
3.3 Identification of miRNAs and differentially expressed miRNAs
The clean reads of these miRNA libraries were subjected to mapping against the miRNABase V22. Furthermore, sRNAmer was used to identify new miRNA by merging all the miRNA sequencing libraries. Totally, we detected 108 known and 34 novel miRNAs in these two dandelion species, with the sequences and genomic locations of miRNA, miRNA star, and miRNA precursor summarized in Supplementary Table 1. To identify the differentially expressed miRNA in the latex between Tk and Ts, we quantified the expression levels of the 142 identified miRNAs (Supplementary Table 2). A total of 53 miRNAs, including 32 known and 21 newly identified ones, were differentially expressed (DE) between the two dandelion accessions (Supplementary Table 2). Among these DE miRNAs, 34 were up-regulated and 19 were down-regulated in the latex of Tk when compared to that of Ts (Figure 3; Supplementary Table 2).
Figure 3 Heatmap of differentially expressed miRNAs. The up-regulated (A) and down-regulated (B) miRNAs are identified in the latex of Tk vs. Ts. The color is based on the log2(TP10M) value.
3.4 miRNA targets analysis
We combined the degradome sequencing and target prediction to identify the miRNA targets. Using the degradome sequencing data, 145 targets corresponding to 74 miRNAs were identified by the CleaveLand4 pipeline. The target plots (T-plots) revealed the comprehensive categories and cleavage information of all these miRNAs and their targets (Supplementary Data 1). The T-plots for six miRNAs and their representative targets are shown in Figure 4A. In order to obtain more complete miRNA-target information, 53 DE miRNAs were also used to predict targets by the software programs TAPIR and psRNATarget. A total of 298 and 306 targets, corresponding to 165 and 164 non-redundant genes, were identified by TAPIR and psRNATarget, respectively (Supplementary Table 3). Among these non-redundant target genes, 127 ones were identified by both types of software (Figure 4B), and subjected further gene ontology (GO) functional annotation. The GO results revealed six statistically significant enrichment for the macromolecule-related pathways, i.e. “regulation of macromolecule biosynthetic process, GO:0010556”, “regulation of macromolecule metabolic process, GO:0060255”, “regulation of cellular macromolecule biosynthetic process, GO:2000112”, “cellular macromolecule metabolic process, GO:0044260”, “macromolecule biosynthetic process, GO:0009059”, and “macromolecule metabolic process, GO:0043170” (Figure 4C; Supplementary Table 4). The root ontology analysis of the enriched GO terms revealed that these pathways related to macromolecules are interconnected and fall within the subcategory of GO:0043170 (Supplementary Figure 1). KEGG enrichment analysis identified the involvement of these target genes in five pathways, particularly the significant enrichment in plant hormone signal transduction (Figure 4D). These observations suggest the DE miRNAs and their targets might participate in the biosynthesis and metabolism of natural rubber, a type of macromolecule.
Figure 4 miRNA target identification and functional enrichment analysis. The six representative miRNA targets are identified by degradome sequencing (A). Venn diagram for targets of differentially expressed miRNAs are predicted by TAPIR and psRNATarget (B). The top 20 GO enrichments are shown for the targets of differentially expressed miRNA (C). The top 5 KEGG enrichments are displayed for the targets of differentially expressed miRNA (D).
3.5 miRNA-target modules implicated in natural rubber metabolism
We compiled the list of miRNA targets involved in macromolecule-related metabolism based on the GO categories and compared their expression levels in the latex of Tk and Ts. Four kinds of target gene families stood out, including AP2/ERF (APETALA2/Ethylene Responsive element binding Factor), ARF (Auxin Response Factor), NAC (NAM, ATAF1/2 and CUC2 domain containing protein), and protein kinase. Of them, six AP2/ERF genes, one NAC gene, and eight protein kinase genes were expressed significantly different between the latex of Tk and Ts (Figure 5A). The DE target genes of AP2/ERF, NAC, ARF, and protein kinase were mainly targeted by miR172, miR164, miR160, and miRN19 families, respectively. To further understand the regulatory networks between these miRNAs and their targets, the miRNA-mRNA interaction network maps were constructed, including the association of miR164 with 15 targets, miRN19a with 28 targets, miR160a with 5 targets, and miR172a with 11 targets. (Figures 5B–E). These DE miRNAs and their target genes might play an important role in the metabolism of natural rubber.
Figure 5 Regulatory networks of miRNA-target involved in the biosynthesis and metabolism of natural rubber macromolecules. Gene families involved in macromolecule biosynthesis and metabolism and their expression levels in Tk and Ts latex (A), the color is based on the log2(FPKM+1) value. ** indicates a significant difference in gene expression between Tk and Ts at the padj-value < 0.01 level. The regulatory networks are shown between miRN19a (B), miR160a (C), miR164a (D) and miR172a (E), and their target genes.
4 Discussion
Plants of the genus Taraxacum, are widespread throughout the world and about 2800 species within 60 sections have been identified (Kirschner et al., 2013). In addition to its long history of use in traditional medicine, the development of rubber-producing dandelions, especially the Russian dandelion T. kok-saghyz (Tk), into alternative rubber crops has regained much attention since the early 2000s (Martinez et al., 2015; Salehi et al., 2021). However, the disadvantages of Tk, such as poor growth vigor and self-incompatibility, make essential the incorporation of other Taraxacum species in Tk breeding and scientific projects (Zeisek et al., 2019; Salehi et al., 2021). In this study, T. spadiceum (Ts), an apomictic dandelion species with strong growth vigor and low root rubber-content, was investigated together with Tk (Figure 1). The findings of miRNA-targets and the vast amount of small RNA datasets obtained here from the two dandelion species will benefit not only the research in rubber production but that in the other important agronomic traits for Tk domestication.
The read distribution of the sRNA sequencing data generally showed that the abundances of 21 and 24 nt sRNAs are the highest in different tissues (Chen, 2009). However, sRNAs of 22 and 23 nt revealed to be highly expressed in our sequencing data (Figure 2). Notably, the abundances of the 22-nt sRNA were the highest in all the three latex samples of Ts when compared to the other length-types of sRNAs (Figure 2B). In Arabidopsis, DCL2 is thought to regulate the biogenesis of viral 22 nt siRNAs, which repress the translation of their cognate mRNAs as well as global translation (Deleris et al., 2006). In addition, the 22 nt sRNA have also been shown play an important role in fine-tuning plant growth and stress responses (Wu et al., 2020). Latex produced by rubber-producing dandelions consists of multiple secondary metabolites and is thought to defend against biotic and abiotic stresses (Ramos et al., 2019). Therefore, the highly expressed 22 nt sRNA in latex might participate in both latex metabolism and plant defense. The mechanisms of generation and physiological functions have never been reported for the type of 23 nt sRNAs in plants. Hence, the exact roles of 22 and 23 nt sRNAs in latex needs to be further investigated, especially their possible involvement in the latex metabolism of the rubber-producing dandelion species. miRNAs and their targets are highly conserved in plants. In this study, 108 conserved miRNAs and 34 dandelion specific miRNAs were identified. Of which, a total of 53 were differentially expressed (DE) between the latex of Tk and Ts, including 32 conserved and 21 newly identified miRNAs (Supplementary Table 2). Among the DE miRNAs, miRN26 and miR319e were, respectively, the most abundant newly identified and conserved miRNAs (Figure 3). These results imply that species-specific miRNAs may play an important role in latex formation and rubber production. We used degradome sequencing to identify the miRNA target genes, and only 145 targets for 74 miRNAs were obtained (Supplementary Data 1). Target genes could be detected for some of the conserved miRNAs, including miR156, miR160, miR166, miR167, miR171, miR319, miR393, and miR396, but not for nearly half of miRNAs identified by degradome sequencing. Four possible reasons might account for this result: the expression abundance of some miRNA targets is too low to be detected; the expression of some target genes exhibits strict spatio-temporal characteristics; some novel and specific miRNAs could have arisen recently and have not necessarily acquired a target; and the sequencing depth of the degradome library needs to be increased. To obtain comprehensive miRNA targets implicated in rubber production, we further used target prediction to identify miRNA target genes for the DE miRNAs between the latex of Tk and Ts. TAPIR and psRNATarget predicted 298 and 307 miRNA targets, respectively, including 165 and 164 non-redundant genes (Supplementary Table 3). Notably, a total of 127 non-redundant targets were detected by both programs. Together, we obtained comprehensive and credible miRNA target information of these DE miRNAs.Terpenoids, including monoterpenes, sesquiterpenes, diterpenes, triterpenes, tetraterpenes and polyterpenes, are types of natural hydrocarbons that are widely distributed in plants. The C5 phosphates, isopentenyl diphosphate (IPP) and dimethylallyl diphosphate (DMAPP), are the common initial molecules for synthesizing of both terpenoid and natural rubber in plants (Hossain et al., 2022). Currently, miRNAs have been widely reported in the regulation of terpenoids biosynthesis. For instance, miR4995 plays a promoter function in the biosynthesis of terpenoids that ultimately affects the production of picroside-I in Picrorhiza kurroa Royle (Vashisht et al., 2015). The miR5021 was predicted to target the enzymes that affect the synthesis of terpenoid indole alkaloids in Catharanthus roseus (L.) G. Don (Pani and Mahapatra, 2013). Overexpression of miRStv_11 and anti-miR319 enhances the biosynthesis of steviol glycosides by suppressing the expression levels of KO, KS, UGT86C2 and KAH in Stevia rebaudiana (Saifi et al., 2019). In this work, GO analysis on the targets of the differentially expressed miRNAs between the latex of Tk and Ts revealed significant enrichment for six pathways involved in the regulation of macromolecule biosynthesis and metabolism (Figure 4; Supplementary Table 4). These target genes might play an important role in the regulation of natural rubber biosynthesis in these two dandelion species with distinct rubber content in their roots (Figure 1). KEGG analysis revealed that the target genes were enriched in five pathways, with plant hormone signal transduction being the most significant pathway. Production practices have demonstrated that extrinsically applying ethylene on rubber trees can enhance latex metabolism and rubber flow time, ultimately leading to a significant increase in rubber yield (Domiciano et al., 2018). JA (jasmonic acid) signaling involves the nature rubber production by regulating the laticifer differentiation in rubber tree (Chao et al., 2023). These results indicate that plant hormones play an important role in the biosynthesis of natural rubber. Four types of miRNA-target modules are predicted to be widely involved in macromolecule-related processes (Figure 5). Notably, the miRN19-protein kinase module is a module specific in the rubber-producing dandelions (Supplementary Table 3). The expression levels of many protein kinase genes in latex were significantly different between Tk and Ts, reinforcing the regulatory roles of the DE miRNAs identified to target protein kinases. Whether the species-specific miRNA-target regulatory modules in Taraxacum participate directly or indirectly in natural rubber production warrants further investigation.
5 Conclusion
In summary, 142 miRNAs were discovered by miRNA sequencing in the two dandelion species, Tk and Ts, with striking differences in rubber production, and 53 are differentially expressed in their latex. The targets of these differentially expressed miRNAs were significantly enriched in the processes of macromolecular biosynthesis and metabolism. Four types of miRNA-target modules, including miR172-AP2/ERF, miR164-ARF, miR160-NAC and miRN19-protein kinase, might play an important role in rubber production for their differential expressions in the latex of the two dandelions. Taken together, the findings and extensive datasets presented here provide a basis for the deeper understanding of the regulatory roles of miRNA in dandelion rubber production.
Data availability statement
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.
Author contributions
CL: Investigation, Resources, Software, Validation, Writing – original draft. YY: Data curation, Formal Analysis, Investigation, Methodology, Writing – original draft. YT: Data curation, Formal Analysis, Methodology, Writing – original draft. XY: Resources, Visualization, Writing – original draft. JC: Methodology, Visualization, Writing – original draft. CT: Funding acquisition, Supervision, Writing – review & editing. KL: Conceptualization, Funding acquisition, Supervision, Writing – original draft.
Funding
The author(s) declare financial support was received for the research, authorship, and/or publication of this article. This work was supported by the Hainan Provincial Natural Science Foundation of China (Grant No. 322RC577), the National Natural Science Foundation of China (Grant No. 32160383, 31825007), Yunnan Academician and Expert Workstation Project (202105AF150017) and the Starting Research Fund of Hainan University.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2023.1287318/full#supplementary-material
Supplementary Figure 1 | Root ontology analysis of the enriched GO terms. The GO terms involving macromolecule-related pathways are labeled in red and displayed in an enlarged font.
References
Addo-Quaye, C., Miller, W., Axtell, M. J. (2009). Cleaveland: A pipeline for using degradome data to find cleaved small RNA targets. Bioinformatics 25, 130–131. doi: 10.1093/bioinformatics/btn604
Ahrends, A., Hollingsworth, P. M., Ziegler, A. D., Fox, J. M., Chen, H., Su, Y., et al. (2015). Current trends of rubber plantation expansion may threaten biodiversity and livelihoods. Global Environ. Change 34, 48–58. doi: 10.1016/j.gloenvcha.2015.06.002
Aukerman, M. J., Sakai, H. (2003). Regulation of flowering time and floral organ identity by a microRNA and its APETALA2-like target genes. Plant Cell 15, 2730–2741. doi: 10.1105/tpc.016238
Axtell, M. J., Meyers, B. C. (2018). Revisiting criteria for plant microRNA annotation in the era of big data. Plant Cell 30, 272–284. doi: 10.1105/tpc.17.00851
Bai, B., Bian, H., Zeng, Z., Hou, N., Shi, B., Wang, J., et al. (2017). MiR393-mediated auxin signaling regulation is involved in root elongation inhibition in response to toxic aluminum stress in barley. Plant Cell Physiol. 58, 426–439. doi: 10.1093/pcp/pcw211
Bonnet, E., He, Y., Billiau, K., Van de Peer, Y. (2010). TAPIR, a web server for the prediction of plant microRNA targets, including target mimics. Bioinformatics 26, 1566–1568. doi: 10.1093/bioinformatics/btq233
Casadevall, R., Rodriguez, R. E., Debernardi, J. M., Palatnik, J. F., Casati, P. (2013). Repression of growth regulating factors by the microRNA396 inhibits cell proliferation by UV-B radiation in arabidopsis leaves. Plant Cell 25, 3570–3583. doi: 10.1105/tpc.113.117473
Chao, J., Wu, S., Shi, M., Xu, X., Gao, Q., Du, H., et al. (2023). Genomic insight into domestication of rubber tree. Nat. Commun. 14, 4651. doi: 10.1038/s41467-023-40304-y
Chen, S., Zhou, Y., Chen, Y., Gu, J. (2018). Fastp: An ultra-fast all-in-one FASTQ preprocessor. Bioinformatics 34, i884–i890. doi: 10.1093/bioinformatics/bty560
Chen, X. (2009). Small RNAs and their roles in plant development. Annu. Rev. Cell Dev. Biol. 25, 21–44. doi: 10.1146/annurev.cellbio.042308.113417
Cornish, K. (2001). Similarities and differences in rubber biochemistry among plant species. Phytochemistry 57, 1123–1134. doi: 10.1016/S0031-9422(01)00097-8
Dai, X., Zhuang, Z., Zhao, P. X. (2018). psRNATarget: A plant small RNA target analysis server, (2017 release). Nucleic Acids Res. 46, W49–W54. doi: 10.1093/nar/gky316
Deleris, A., Gallego-Bartolome, J., Bao, J., Kasschau, K. D., Carrington, J. C., Voinnet, O. (2006). Hierarchical action and inhibition of plant Dicer-like proteins in antiviral defense. Science 313, 68–71. doi: 10.1126/science.1128214
Domiciano, D., Carvalho, P., Oliveira, L., Cárdenaz, H., Prudente, D. (2018). Ethrel stimulant increases the activity of soluble invertase isoforms in rubber tree (hevea brasiliensis) bark tissues. Aust. J. Crop Sci. 12, 99–107. doi: 10.21475/ajcs.18.12.01.pne770
Fan, R., Li, Y., Li, C., Zhang, Y. (2015). Differential microRNA analysis of glandular trichomes and young leaves in xanthium strumarium L. reveals their putative roles in regulating terpenoid biosynthesis. PLoS One 10, e0139002. doi: 10.1371/journal.pone.0139002
Feng, L., Xia, R., Liu, Y. (2019). Comprehensive characterization of miRNA and PHAS loci in the diploid strawberry (Fragaria vesca) genome. Hortic. Plant J. 5, 255–267. doi: 10.1016/j.hpj.2019.11.004
Guo, H.-S., Xie, Q., Fei, J.-F., Chua, N.-H. (2005). MicroRNA directs mRNA cleavage of the transcription factor NAC1 to downregulate auxin signals for arabidopsis lateral root development. Plant Cell 17, 1376–1386. doi: 10.1105/tpc.105.030841
Gupta, O. P., Karkute, S. G., Banerjee, S., Meena, N. L., Dahuja, A. (2017). Contemporary understanding of miRNA-based regulation of secondary metabolites biosynthesis in plants. Front. Plant Sci. 8, 374. doi: 10.3389/fpls.2017.00374
Hossain, R., Quispe, C., Saikat, A. S. M., Jain, D., Habib, A., Janmeda, P., et al. (2022). Biosynthesis of secondary metabolites based on the regulation of microRNAs. BioMed. Res. Int. 2022, 9349897. doi: 10.1155/2022/9349897
Jones-Rhoades, M. W., Bartel, D. P., Bartel, B. (2006). MicroRNAs and their regulatory roles in plants. Annu. Rev. Plant Biol. 57, 19–53. doi: 10.1146/annurev.arplant.57.032905.105218
Karimi, A. A., Naghavi, M. R., Peyghambari, S. A., Sobhani, A., Rasoulnia, A. (2022). Identification of miRNAs and their target genes in Taraxacum spp. J. Agr. Sci. Tech. 24, 1457–1471. doi: 10.52547/jast.24.6.1457
Kirschner, J., Štěpánek, J., Černý, T., De Heer, P., van Dijk, P. J. (2013). Available exsitu germplasm of the potential rubber crop Taraxacum koksaghyz belongs to a poor rubber producer, T. brevicorniculatum (compositae–crepidinae). Genet. Resour. Crop Evol. 60, 455–471. doi: 10.1007/s10722-012-9848-0
Kisiel, W., Michalska, K. (2006). Matricarin-type guaianolides from Taraxacum bessarabicum and their chemotaxonomic significance. Biochem. Syst. Ecol. 34, 356–359. doi: 10.1016/j.bse.2005.09.006
Langmead, B., Trapnell, C., Pop, M., Salzberg, S. L. (2009). Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 10, R25. doi: 10.1186/gb-2009-10-3-r25
Leu, Y.-L., Shi, L.-S., Damu, A. G. (2003). Chemical constituents of Taraxacum formosanum. Chem. Pharm. Bull. 51, 599–601. doi: 10.1248/cpb.51.599
Lin, T., Xu, X., Du, H., Fan, X., Chen, Q., Hai, C., et al. (2021). Extensive sequence divergence between the reference genomes of Taraxacum kok-saghyz and Taraxacum mongolicum. Sci. China Life Sci. 65, 515–528. doi: 10.1007/s11427-021-2033-2
Lin, T., Xu, X., Ruan, J., Liu, S., Wu, S., Shao, X., et al. (2017). Genome analysis of Taraxacum kok-saghyz rodin provides new insights into rubber biosynthesis. Natl. Sci. Rev. 5, 78–87. doi: 10.1093/nsr/nwx101
Love, M. I., Huber, W., Anders, S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with Deseq2. Genome Biol. 15, 550. doi: 10.1186/s13059-014-0550-8
Mallory, A. C., Bartel, D. P., Bartel, B. (2005). MicroRNA-directed regulation of arabidopsis AUXIN RESPONSE FACTOR17 is essential for proper development and modulates expression of early auxin response genes. Plant Cell 17, 1360–1375. doi: 10.1105/tpc.105.031716
Martin, M. (2011). Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet J. 17, 10–12. doi: 10.14806/ej.17.1.200
Martinez, M., Poirrier, P., Chamy, R., Prüfer, D., Schulze-Gronover, C., Jorquera, L., et al. (2015). Taraxacum officinale and related species—an ethnopharmacological review and its potential as a commercial medicinal plant. J. Ethnopharmacol. 169, 244–262. doi: 10.1016/j.jep.2015.03.067
Men, X., Wang, F., Chen, G.-Q., Zhang, H.-B., Xian, M. (2019). Biosynthesis of natural rubber: Current state and perspectives. Int. J. Mol. Sci. 20, 50. doi: 10.3390/ijms20010050
Metcalfe, C. R. (1967). Distribution of latex in the plant kingdom. Econ. Bot. 21, 115–127. doi: 10.1007/BF02897859
Millar, A. A., Gubler, F. (2005). The arabidopsis GAMYB-like genes, MYB33 and MYB65, are microRNA-regulated genes that redundantly facilitate anther development. Plant Cell 17, 705–721. doi: 10.1105/tpc.104.027920
Mooibroek, H., Cornish, K. (2000). Alternative sources of natural rubber. Appl. Microbiol. Biotechnol. 53, 355–365. doi: 10.1007/s002530051627
Pani, A., Mahapatra, R. K. (2013). Computational identification of microRNAs and their targets in Catharanthus roseus expressed sequence tags. Genom. Data 1, 2–6. doi: 10.1016/j.gdata.2013.06.001
Pertea, M., Kim, D., Pertea, G. M., Leek, J. T., Salzberg, S. L. (2016). Transcript-level expression analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nat. Protoc. 11, 1650–1667. doi: 10.1038/nprot.2016.095
Post, J., van Deenen, N., Fricke, J., Kowalski, N., Wurbs, D., Schaller, H., et al. (2012). Laticifer-specific cis-prenyltransferase silencing affects the rubber, triterpene, and inulin content of Taraxacum brevicorniculatum. Plant Physiol. 158, 1406–1417. doi: 10.1104/pp.111.187880
Puskas, J. E., Chiang, K., Barkakaty, B. (2014). “Natural rubber (NR) biosynthesis: Perspectives from polymer chemistry,” in Chemistry, Manufacture and Applications of Natural Rubber. Eds. Kohjiya, S., Ikeda, Y. (Cambridge, England: Woodhead Publishing). doi: 10.1533/9780857096913.1.30
Qadir, I., Nazir, S., Sheikh, M. A., Naaz, F., Bashir, S., Ovais, S., et al. (2022). “Taraxacum officinale: the esculent dandelion as herbal medicine,” in Edible Plants in Health and Diseases : Volume ii : Phytochemical and Pharmacological Properties. Eds. Masoodi, M. H., Rehman, M. U. (Singapore: Springer).
Ramirez-Cadavid, D. A., Cornish, K., Michel, F. C. (2017). Taraxacum kok-saghyz (TK): compositional analysis of a feedstock for natural rubber and other bioproducts. Ind. Crops Prod. 107, 624–640. doi: 10.1016/j.indcrop.2017.05.043
Ramos, M. V., Demarco, D., da Costa Souza, I. C., de Freitas, C. D. T. (2019). Laticifers, latex, and their role in plant defense. Trends Plant Sci. 24, 553–567. doi: 10.1016/j.tplants.2019.03.006
Rousset, A., Amor, A., Punvichai, T., Perino, S., Palu, S., Dorget, M., et al. (2021). Guayule (Parthenium argentatum A. Gray), a renewable resource for natural polyisoprene and resin: composition, processes and applications. Molecules 26, 664. doi: 10.3390/molecules26030664
Saifi, M., Yogindran, S., Nasrullah, N., Nissar, U., Gul, I., Abdin, M. Z. (2019). Co-expression of anti-miR319g and miRstv-11 lead to enhanced steviol glycosides content in Stevia rebaudiana. BMC Plant Biol. 19, 274. doi: 10.1186/s12870-019-1871-2
Salehi, M., Cornish, K., Bahmankar, M., Naghavi, M. R. (2021). Natural rubber-producing sources, systems, and perspectives for breeding and biotechnology studies of Taraxacum kok-saghyz. Ind. Crops Prod. 170, 113667. doi: 10.1016/j.indcrop.2021.113667
Schütz, K., Carle, R., Schieber, A. (2006). Taraxacum—a review on its phytochemical and pharmacological profile. J. Ethnopharmacol. 107, 313–323. doi: 10.1016/j.jep.2006.07.021
Schwab, R., Palatnik, J. F., Riester, M., Schommer, C., Schmid, M., Weigel, D. (2005). Specific effects of microRNAs on the plant transcriptome. Dev. Cell 8, 517–527. doi: 10.1016/j.devcel.2005.01.018
Shannon, P., Markiel, A., Ozier, O., Baliga, N. S., Wang, J. T., Ramage, D., et al. (2003). Cytoscape: A software environment for integrated models of biomolecular interaction networks. Genome Res. 13, 2498–2504. doi: 10.1101/gr.1239303
Song, X., Li, Y., Cao, X., Qi, Y. (2019). MicroRNAs and their regulatory roles in plant–environment interactions. Annu. Rev. Plant Biol. 70, 489–525. doi: 10.1146/annurev-arplant-050718-100334
Sweeney, B., Vora, M., Ulbricht, C., Basch, E. (2005). Evidence-based systematic review of dandelion (Taraxacum officinale) by natural standard research collaboration. J. Herb. Pharmacother. 5, 79–93. doi: 10.1080/J157v05n01_09
Tang, C., Yang, M., Fang, Y., Luo, Y., Gao, S., Xiao, X., et al. (2016). The rubber tree genome reveals new insights into rubber production and species adaptation. Nat. Plants 2, 16073. doi: 10.1038/nplants.2016.73
Thiebaut, F., Rojas, C. A., Almeida, K. L., Grativol, C., Domiciano, G. C., Lamb, C. R. C., et al. (2012). Regulation of miR319 during cold stress in sugarcane. Plant Cell Environ. 35, 502–512. doi: 10.1111/j.1365-3040.2011.02430.x
van Beilen, J. B., Poirier, Y. (2007a). Establishment of new crops for the production of natural rubber. Trends Biotechnol. 25, 522–529. doi: 10.1016/j.tibtech.2007.08.009
van Beilen, J. B., Poirier, Y. (2007b). Guayule and Russian dandelion as alternative sources of natural rubber. Crit. Rev. Biotechnol. 27, 217–231. doi: 10.1080/07388550701775927
Vashisht, I., Mishra, P., Pal, T., Chanumolu, S., Singh, T. R., Chauhan, R. S. (2015). Mining NGS transcriptomes for miRNAs and dissecting their role in regulating growth, development, and secondary metabolites production in different organs of a medicinal herb, Picrorhiza kurroa. Planta 241, 1255–1268. doi: 10.1007/s00425-015-2255-y
Wu, H., Li, B., Iwakawa, H.-o., Pan, Y., Tang, X., Ling-hu, Q., et al. (2020). Plant 22-nt siRNAs mediate translational repression and stress adaptation. Nature 581, 89–93. doi: 10.1038/s41586-020-2231-y
Wu, G., Poethig, R. S. (2006). Temporal regulation of shoot development in Arabidopsis thaliana by miR156 and its target SPL3. Development 133, 3539–3547. doi: 10.1242/dev.02521
Yan, J., Zhao, C., Zhou, J., Yang, Y., Wang, P., Zhu, X., et al. (2016). The miR165/166 mediated regulatory module plays critical roles in ABA homeostasis and response in Arabidopsis thaliana. PloS Genet. 12, e1006416. doi: 10.1371/journal.pgen.1006416
Yang, C., Li, D., Mao, D., Liu, X. U. E., Ji, C., Li, X., et al. (2013). Overexpression of microRNA319 impacts leaf morphogenesis and leads to enhanced cold tolerance in rice (Oryza sativa L.). Plant Cell Environ. 36, 2207–2218. doi: 10.1111/pce.12130
Zeisek, V., Kirschner, J., van Dijk, P. J., Štěpánek, J., Černý, T., Kotek, J., et al. (2019). Analysis of wild Taraxacum bicorne dahlst. (compositae–crepidinae) as a potential alternative natural rubber crop. Genet. Resour. Crop Evol. 66, 1341–1361. doi: 10.1007/s10722-019-00788-4
Keywords: miRNA, miRNA targets, regulatory modules, natural rubber production, Taraxacum
Citation: Liang C, Yan Y, Tan Y, Yang X, Cao J, Tang C and Liu K (2023) Identification of miRNAs and their targets in two Taraxacum species with contrasting rubber-producing ability. Front. Plant Sci. 14:1287318. doi: 10.3389/fpls.2023.1287318
Received: 05 September 2023; Accepted: 25 October 2023;
Published: 08 November 2023.
Edited by:
Peng Wang, Jiangsu Province and Chinese Academy of Sciences, ChinaReviewed by:
Sylvain Legrand, Université de Lille, FranceGulshan Kumar, University of Kentucky, United States
Copyright © 2023 Liang, Yan, Tan, Yang, Cao, Tang and Liu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Kaiye Liu, kaiyeliu@hainanu.edu.cn; Chaorong Tang, chaorongtang@hainanu.edu.cn
†These authors have contributed equally to this work