Skip to main content

Integrated analysis of mRNA-seq and miRNA-seq reveals the advantage of polyploid Solidago canadensis in sexual reproduction

Abstract

Background

The invasion of Solidago canadensis probably related to polyploidy, which may promotes its potential of sexual reproductive. S. canadensis as an invasive species which rapidly widespread through yield huge numbers of seed, but the mechanism remains unknown. To better understand the advantages of sexual reproduction in hexaploid S. canadensis, transcriptome and small RNA sequencing of diploid and hexaploid cytotypes in flower bud and fruit development stages were performed in this study.

Results

The transcriptome analysis showed that in the flower bud stage, 29 DEGs were MADS-box related genes with 14 up-regulated and 15 down-regulated in hexaploid S. canadensis; 12 SPL genes were detected differentially expressed with 5 up-regulated and 7 down-regulated. In the fruit development stage, 26 MADS-box related genes with 20 up-regulated and 6 down-regulated in hexaploid S. canadensis; 5 SPL genes were all up-regulated; 28 seed storage protein related genes with 18 were up-regulated and 10 down-regulated. The weighted gene co-expression network analysis (WGCNA) identified 19 modules which consisted of co-expressed DEGs with functions such as sexual reproduction, secondary metabolism and transcription factors. Furthermore, we discovered 326 miRNAs with 67 known miRNAs and 259 novel miRNAs. Some of miRNAs, such as miR156, miR156a and miR156f, which target the sexual reproduction related genes.

Conclusion

Our study provides a global view of the advantages of sexual reproduction in hexaploid S. canadensis based on the molecular mechanisms, which may promote hexaploid S. canadensis owing higher yield and fruit quality in the process of sexual reproduction and higher germination rate of seeds, and finally conductive to diffusion, faster propagation process and enhanced invasiveness.

Peer Review reports

Background

Biological invasion acted as a global threat to biodiversity [1]. Many invasive species pose a serious threat to local biodiversity, the invasion of alien plants will affect and simplify the new state of ecosystem [2]. The invasive plant show significantly negative impacts on the native ecosystems. Thus, plant invasion will ultimately affect human and society include reduction of agricultural production, obstacle of sustainable development and influence on human health [3]. Many studies have been dedicated to explore the mechanism of plant invasion [4,5,6]. At present, no common mechanism has been found to explain and predict that why some alien species rapidly expand beyond their local range and become dominant in non-native habitats [7, 8]. Thus, understanding the factors that support the rapid spread and successful growth in invasive species is necessary to design effective control strategies.

Reproduction as an essential role for invaders when they are introduced into a new habitat [9]. Thus, reproductive traits especially for the sexual reproduction generally have developed into important determinants of invasion [10, 11]. The advantages of sexual reproduction include the increased offspring diversity, reduced intraspecific competition and wider distribution potential for rapidly initiating new populations far from the parental plants [12]. Invasive plants have greater advantages in seed production capacity or seed size in the introduced range to speed up naturalization and spread [13]. Some invasive species spread rapidly via seed reproduction in the invasion process [11, 14]. Therefore, knowing the priority of sexual reproduction mechanism in the invasion success of exotic species is essential to understanding invasion strategies for invasive species.

Polyploid plants have been developed a series of drastic competitive traits. Polyploidy usually are morphological different to diploid [15]. This significant morphological variation will enhance plant vigor that may confer polyploids a pre-adaptive advantage, making them more likely to become invasive weeds than their diploid ancestors [16, 17]. Polyploid plant possessed a higher competitive advantage in the early stages of invasion due to the higher germination rate, faster growth and more robust seedlings. A greater sexual and vegetative reproduction ability may be helpful to establish and expand polyploid population [18]. Therefore, polyploidy might be one of the potential determinants in invasion success through a greater competitive ability than diploids [17, 19]. The comparative study of invasive species and native congeners is a common means to elucidate the evolution of competitiveness [20, 21], and many studies have revealed differences in the competitiveness of cytotypes in some alien invasive species [22,23,24]. However, a comprehensive overview of ploidy levels underlying molecular mechanism across the most invasive plant species is still lacking. Therefore, the studies of the difference of alien plant with different cytotypes from molecular basis levels will provide a better understanding of its molecular mechanisms about invasiveness for polyploids. Furthermore, such information may be essential for establishing the importance of ploidy in determining invasion success.

The whole genome duplication will lead to large-scale changes in the gene expression level [25,26,27]. These large scale changes are mostly including chromatin remodeling, methylation, and the production and function of small RNA [27,28,29]. MicroRNAs (miRNAs) are a class of small noncoding RNAs, which consisting of approximately 20-24 nucleotides in eukaryotes. In general, miRNA negatively regulate gene expression through complementarity to target mRNAs at the post-translational level [30, 31].

In recent years, small RNAs, especially miRNAs, have become important regulators of plant growth, development and stress response [32,33,34]. Some studies have shown that miRNAs play an important role in the sexual reproductive process of plants, which including the regulation of flowering time and flower development [34,35,36,37], and fruit ripening, seed yield and development [37,38,39,40]. In addition, the regulation function and expression pattern of miRNAs were changed and elucidated important molecular mechanisms in polyploid plant [41,42,43].

Solidago canadensis was a perennial herb with three cytotypes, diploid (2n = 2x = 18), tetraploid (2n = 4x = 36) and hexaploid (2n = 6x = 54), which has been introduced into many areas in the world from North American and has developed into invasive species including China [44]. Based on the different cytotypes, S. canadensis was also showed the different invasiveness ability. Polyploidy facilitate S. canadensis becoming competitive ability and promotes its successful formation of invasive species in China [22]. For example, polyploidization can promote S. canadensis to adapt in new environment [45]. The allelopathy potential of S. canadensis were enhanced by polyploidization, therefore, making the introduced polyploid S. canadensis more competitive than diploid, which helps its successful invasion [46]. Besides, reproductive traits of invasive plants are important determinant of their invasion. S. canadensis as an invasive plant can spread by sexual reproduction and asexual reproduction by rhizome. Studies suggested that sexual reproduction facilitates new S. canadensis populations establishment [44]. The advantage of invasiveness of S. canadensis has owing to the abilities to yield a large number of seeds, which is closely related to the number of inflorescences [47]. The molecular basis for the allelopathic metabolite synthesis has suggested that the altered related gene expression may enhanced invasive potential of polyploid S. canadensis [48]. While few studies fully clarified the molecular mechanism for the sexual reproduction in S. canadensis.

In this study, we explored the expression of reproduction related genes and miRNAs between hexaploid and diploid S. canadensis. Based on the gene and miRNA expression characteristics in the inflorescence of flower bud and fruit development stage between two cytotypes of S. canadensis, we elucidate the biological function of these genes and miRNAs. Furthermore, we combined the expression profile of gene with miRNA to construct miRNA-mRNA interaction network to explore the regulatory action of miRNAs. Our work will finally provided the theoretical basis for further exploration of the invasion mechanism of polyploid S. canadensis.

Result

Gene expression patterns in diploid and hexaploid cytotypes of S. canadensis

Total RNA was extracted from flower bud stage and fruit development in hexaploid (HFa and HFb) and diploid (DFa and DFb) samples with both of each group including 3 biological replicates. To explore the gene expression patterns in two cytotypes of S. canadensis, the raw data which obtained by RNA sequencing of the flower bud and fruit development stage in the hexaploid and diploid S. canadensis (a total of 12 samples) were filtered. The clean data for each sample was not less than 6 G with more than 99% Q20 bases and more than 96% Q30 bases were acquired (Additional file 1: Table S1). Based on the assembly of transcripts and the removal of low-abundance expressed genes, a total of 121,278 expressed genes were detected, of which 67,095, 73,360, 65,976 and 71,551 were detected in HFa, HFb, DFa and DFb respectively. Among all the detected genes, 30,602 genes were expressed in all four groups, 8035 genes were specifically expressed in HFa, 14,098 genes were specifically expressed in HFb, 5839 genes were specifically expressed in DFa, and 12,478 genes were specifically expressed in DFb (Fig. 1a). Furthermore, based on the alignment of differentially expressed genes (DEGs) between two cytotypes with the Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases, we found that 7892 genes were differentially expressed in the flower bud stage and 8990 genes in the fruit development stage. Among these differentially expressed genes, 4197 genes were detected in both flower bud and fruit development stage, 3695 genes were only detected in flower bud stage, while 4793 genes were only detected in fruit development stage (Fig. 1b).

Fig. 1
figure1

The distribution of expressed genes in each sample and differentially expressed genes between two cytotypes. a The distribution of expressed genes in the flower bud and fruit development stages of hexaploid (HFa and HFb) and diploid S. canadensis (DFa and DFb). b The distribution of differentially expressed genes at flower bud stage (DFa-vs-HFa) and fruit development stage (DFb-vs-HFb)

GO functional analysis of DEGs

To explore the insight into functional categories of the DEGs in the flower bud and fruit development stages, a GO classification was performed. In the flower bud stage, DEGs were distributed in 46 GO terms, which were mainly included in three categories: biology process, cellular component and molecular function. In the biology process category, most of DEGs were distributed in the subcategories of ‘cellular process’, ‘metabolic process’ and ‘biological regulation’ (Additional file 2: Fig. S1). In the cellular component, most of DEGs were distributed in the subcategories of ‘cell’, ‘cell part’ and ‘membrane’ (Additional file 2: Fig. S1). In the molecular function, most of DEGs were distributed in the subcategories of ‘catalytic activity’, ‘binding’ and ‘transporter activity’ (Additional file 2: Fig. S1). Furthermore, GO functional enrichment analysis was used to collect genes that play important roles in the flower bud stage. Based on the functionally de-redundant analysis of significantly enriched GO terms by REVIGO, we screened out many GO terms that are involved in biological functions such as ‘nucleosome’, ‘response to stress’, ‘mature ribosome assembly’, ‘protein catabolism’, ‘double fertilization forming a zygote and endosperm’ and ‘ethylene biosynthesis’ (Fig. 2a). The expression of genes which involved in these GO terms might be altered and played a vital role in the flower bud stage of hexaploid S. canadensis. In the fruit development stage, DEGs were distributed in 50 GO terms, which were also mainly included in three categories. In the biology process and cellular component category, most of DEGs were distributed similarly with flower bud stage (Additional file 3: Fig. S2). In the molecular function, most of DEGs were distributed in the subcategories of ‘catalytic activity’, ‘binding’ and ‘structural molecule activity’ (Additional file 3: Fig. S2). GO functional enrichment and functionally de-redundant analysis in the fruit development stage were also performed. Compared with the flower bud stage, there were no significant changes in the functional categories of DEGs in fruit development stage, while some changes were appeared in the enriched GO terms (Fig. 2b).

Fig. 2
figure2

GO analysis of DEGs. a The DFa-vs-HFa comparision; b The DFb-vs-HFb comparision. Each rectangle indicates the enriched GO terms

Pathway functional analysis of DEGs

The KEGG pathway classification and functional enrichment analysis were carried out to further identify the function of DEGs in flower bud and fruit development stages. The DEGs in the flower bud stage were mapped into 131 KEGG pathways, which most of genes were contained into the pathway of ‘Carbon metabolism’, ‘Ribosome’, ‘RNA transport’ and ‘Ubiquitin mediated proteolysis’ (Additional file 4 Table S2). The top 20 enriched pathways ranked by Q-value were mainly about ‘plant hormone signal transduction’, ‘circadian rhythm-plant’ ‘biosynthesis of amino acids’ and some of secondary metabolism related pathways such as ‘flavonoid biosynthesis’, ‘flavone and flavonol biosynthesis’ and ‘anthocyanin biosynthesis’ (Fig. 3a). In the fruit development stage, DEGs were also mapped into 131 KEGG pathways and most of genes were contained into the same pathways with flower bud stage (Additional file 5 Table S3). Some of the top 20 enriched pathways which ranked by Q-value were similar with the flower bud stage, such as ‘carbon metabolism’, ‘flavonoid biosynthesis’ and ‘biosynthesis of amino acids’. Some pathways, such as ‘ribosome’, ‘phenylpropanoid biosynthesis’ and ‘linoleic acid metabolism’ were specially collected (Fig. 3b). The DEGs which enriched in these pathways may play an important role in the reproductive development of hexaploid S. canadensis.

Fig. 3
figure3

KEGG pathway enrichment analysis of DEGs. a The DFa-vs-HFa comparision; b The DFb-vs-HFb comparision. The Y-axis indicates the KEGG pathway, the X-axis indicates the rich factor. The bubble size indicates the number of DEGs of the pathway, and the colour bar indicates the q-value

Secondary metabolism related DEGs

According to the MapMan analysis, many DEGs which involved in secondary metabolite synthesis between two cytotypes of S. canadensis were collected in flower bud and fruit development stage (Fig. 4). In the flower bud stage, most of DEGs which mapped into the ‘phenlypropanoids’, ‘flavonoids’ and ‘terpenoids’ pathways were down-regulated in hexaploid S. canadensis (Fig. 4a). For example, the PAL (phenylalanine ammonia lyase), C4H (cinnamate 4-hydroxylase) and 4CL (4-coumaroyl: CoA ligase) gene, which play a key role in ‘phenlypropanoids’ pathway, were down regulated (Additional file 6: Fig. S3). In addition, many DEGs related with terpenoid biosynthesis and metabolism were up-regulated. For example, some key enzyme genes which involve in the upstream of MVA pathway, such as ACAA, HMGCR, and PMVK genes were up regulated in hexaploid S. canadensis (Additional file 7: Fig. S4). In the fruit development stage, many DEGs mapped into the ‘phenlypropanoids’, ‘flavonoids’ and ‘terpenoids’ were down-regulated in hexaploid S. canadensis (Fig. 4b). While, in the phenylpropanoid biosynthesis pathway, PAL, C4H, CHS and CHI genes were up-regulated in hexaploid S. canadensis (Additional file 8: Fig. S5). In terpenoid synthesis related pathways, most of terpenoid biosynthesis related DEGs were down-regulated in hexaploid S. canadensis (Additional file 9: Fig. S6).

Fig. 4
figure4

Mapman analysis of DEGs. a The DFa-vs-HFa comparision; b The DFb-vs-HFb comparision. The color bar indicates the normalized log2 transformed fold change value

Expression profiling of transcription factor (TF) genes

According to the differential expression analysis of TF genes between two cytotypes of S. canadensis in flower bud stage, 433 differentially expressed TF genes were obtained to be classified into 45 TF gene families. The most number of TF genes were MYB genes, with 27 were up-regulated and 25 were down regulated in hexaploid S. canadensis. As follows, 43 AP2-EREBP TF genes with 21 were up-regulated and 22 were down-regulated; 32 NAC TF genes with 12 were up-regulated and 20 were down-regulated; 31 bHLH TF genes with 16 were up-regulated and 15 were down-regulated; 29 MADS TF genes with 14 were up-regulated and 15 were down-regulated (Fig. 5a). In the fruit development stage, 448 differentially expressed TF genes were classified into 48 TF gene families. The most number of TF genes were AP2-EREBP genes with 23 were up-regulated and 37 were down-regulated in hexaploid S. canadensis. As follows were 40 MYB TF genes, of which 21 were up-regulated and 19 were down-regulated. In addition, other TF genes were also have different expression patterns compared with flower bud stage, such as NAC, WRKY and C3H TF genes (Fig. 5b). All of these differentially expressed TFs may play a major role in the growth and development in S. canadensis. Compared with diploid, the regulation mechanism in hexaploid S. canadensis may be altered.

Fig. 5
figure5

The number of differentially expressed TFs. a The DFa-vs-HFa comparision; b The DFb-vs-HFb comparision. The red colour indicates up-regulated TFs, the blue colour indicates down-regulated TFs

Sexual reproductive development related DEGs

Based on the pathway enrichment analysis, 79 DEGs were collected in plant-circadian rhythm pathway in the flower bud stage (Table 1). The expression pattern of these genes were various between two cytotypes S. canadensis. For example, there were 7 PHYB genes with 5 were up-regulated in hexaploid S. canadensis, and 2 down-regulated; 2 GI genes were up-regulated, a CCA1 gene was down-regulated. Furthermore, 6 LHY genes with 3 up-regulated and 3 down-regulated; 7 CO genes with 5 up-regulated and 2 down-regulated. These DEGs may be involved in regulating the flowering process and flowering time in S. canadensis.

Table 1 Differentially expressed genes in the plant-circadian during flower bud stage

We have identified 29 differentially expressed MADS-box related genes in the flower bud stage, of which 14 were up-regulated in hexaploid S. canadensis and 15 were down-regulated. Among them were mainly about AGL (AGL11, AGL8 and AGL9 etc.), AP1 and SEP1 genes (Table 2). In the fruit development stage, 26 differentially expressed MADS-box related genes were detected, of which 20 were up-regulated and 6 were down-regulated in hexaploid S. canadensis (Table 3). There were 16 genes co-differentially expressed between two stages, 13 genes were differentially expressed specially in the flower bud stage, and 10 genes were differentially expressed specially in fruit development stage. According to these differentially expressed genes, which belong to MADS TF gene families, we inferred that hexaploid S. canadensis may be different with diploid in flower morphogenesis, fruit ripening and development.

Table 2 Differentially expressed MADS-box genes in the flower bud stage
Table 3 Differentially expressed MADS-box genes in the fruit development stage

Besides, 12 differentially expressed SPL genes were detected in flower bud stage, 5 of which were up-regulated and 7 down-regulated in hexaploid S. canadensis. In the fruit development stage, the detected DEGs with 5 SPL genes were up-regulated. In addition, 28 seed storage protein related genes were also detected in fruit development stage, 18 of which were up-regulated and 10 down-regulated in hexaploid. These genes played an important role in improving seed yield and quality. Therefore, the increased expression of these genes may play an important role of sexual reproduction in hexaploid S. canadensis.

Gene expression network analysis

A total of 12,685 DEGs from two comparison groups in flower bud and fruit development stage were analyzed by the weighted gene co-expression network analysis (WGCNA). We identified 19 gene expression modules, which was classified the number of gene from 30 to 2484 (Fig. 6a). Based on the correlation analysis of the module and trait, many co-regulatory modules were collected (correlation coefficient > 0.6, p < 0.05) (Fig. 6b).

Fig. 6
figure6

WGCNA analysis. a Clustering dendrogram of genes, with dissimilarity based on topological overlap, together with assigned merged module colors and theory final module colors. b Module–trait associations. The colour bar indicates the correlation of each module with each sample. Each cell contains the corresponding correlation and p-value

In the HFa group, four modules with cyan, grey60, green and red were collected. To examine the transcript profiles of these modules, we depicted the heatmaps of genes across all samples. The cyan module with 61 genes were collected, most of which were up-regulated in DFa group and HFa group and down-regulated in DFb and HFb group (Fig. 7a). The network analysis showed that many of functional genes were included in this module (Fig. 8a), such as starch synthase gene (CL15636.Contig7_All), histone H3 gene (CL18639.Contig3_All, Unigene44875_All), histone H2A gene (CL19683.Contig1_All, CL19683.Contig5_All), histone H2B gene (CL8238.Contig2_All), auxin-responsive protein gene (CL7139.Contig3_All), gibberellin receptor GID1 gene (Unigene56138_All), MADS-box TF gene (CL9917.Contig5_All), MYB TF gene (Unigene52378_All) and FLOWERING LOCUS T gene (Unigene832_All). In the green module, 1337 genes were collected, most of these genes were up-regulated in HFa and HFb group and down-regulated in DFa and DFb group (Fig. 7b). The network representation of 150 genes with WGCNA edge weight > 0.46 for further analyses (Fig. 8b), which included ubiquitin related genes (CL2636.Contig2_All, Unigene17130_All, CL2880.Contig5_All, CL4786.Contig3_All, CL1625.Contig11_All, CL1677.Contig2_All and CL13448.Contig2_All), histone related genes (CL1658.Contig4_All and Unigene37108_All), translation initiation factor genes (Unigene75377_All, CL17267.Contig3_All, CL17350.Contig2_All, Unigene28891_All and CL16220.Contig2_All) and some TF genes such as MYB gene (CL17002.Contig1_All), SBP gene (CL3869.Contig1_All) and TCP gene (CL2760.Contig4_All). In the red module, 1077 genes were collected, most of these genes were up-regulated in HFa group and down-regulated in DFa, DFb and HFb group (Fig. 7c). The network representation of 154 genes with WGCNA edge weight > 0.47 were used for further analyses (Fig. 8c), which included flavone and flavonol biosynthesis genes (CL11.Contig1_All and CL11.Contig2_All), ubiquitin related genes (CL19520.Contig1_All, Unigene17221_All, CL6885.Contig6_All, CL1548.Contig2_All, CL1203.Contig2_All, CL12120.Contig5_All, CL7593.Contig1_All and CL13587.Contig1_All), MADS-box TF gene (CL8783.Contig1_All, Unigene24590_All), ARF gene (Unigene25133_All), CHS gene (CL14895.Contig7_All) and SBP gene (CL2806.Contig1_All).

Fig. 7
figure7

The heatmap of gene expression level in all samples of selected WGCNA modules. The colour bar indicates the relative expression of module genes, red denotes up-regulation, green denotes down-regulation

Fig. 8
figure8

Interaction analyses of the selected modules. The bubble size indicates connect degree of each gene

In the HFb group, three modules with salmon, lightcyan, and turquoise were collected. In the lightcyan module, 41 genes were collected, most of these genes were down-regulated in DFa and DFb group, and up-regulated in HFa and HFb group (Fig. 7d). The network analysis showed that some TF genes and secondary metabolite synthesis related genes were included in this module (Fig. 8d), such as cinnamyl-alcohol dehydrogenase (CAD) gene (Unigene18113_All and Unigene5751_All), flavonol synthase (FLS) gene (CL6377.Contig4_All), carotenoid cleavage dioxygenase (CCD) gene (CL11186.Contig2_All), Tify (CL2464.Contig12_All) and MYB (CL12049.Contig1_All) gene. In the turquoise module, 2484 genes were collected, most of these genes were up-regulated in HFa and HFb group, and down-regulated in DFa and DFb group (Fig. 7e). The network representation of 102 genes with WGCNA edge weight > 0.45 were used for further analyses (Fig. 8e), which included ubiquitin related genes (CL15905.Contig2_All, CL5915.Contig4_All and CL779.Contig7_All), 4CL gene (Unigene28994_All), calcium-dependent protein kinase gene (Unigene26697_All), FLS gene (CL11344.Contig2_All) and germacrene D synthase gene (CL16768.Contig3_All). These collected genes may be vital for the regulation of sexual reproduction, metabolism and hormone signal transduction related processes in hexaploid S. canadensis.

Small RNA profiling in diploid and hexaploid cytotypes of S. canadensis

Despite the numerous mRNAs involved in sexual reproduction, miRNAs, which, regulate mRNA expressions at the post-transcriptional level, also play vital roles in floral development. To investigate the small RNAs component and the dynamic changes of miRNAs between two cytotypes of S. canadensis, twelve small RNA libraries were generated from the flower bud and fruit development stage samples. After the deep sequencing, and removing low quality reads, adaptor, insert and poly (A) contaminations, average about 25.2, 26.7, 27.2 and 27.1 million clean reads were obtained in three replicate samples of DFa, DFb, HFa and HFb respectively (Additional file 10: Table S4). The length of distributions of miRNAs were range from 18 to 30 nucleotides, and the lengths of miRNAs in all samples were concentrated between 21 and 24 nt (Additional file 11: Fig. S7).

A total 67 known miRNAs were finally identified in the flower bud and fruit development stage of S. canadensis, of which 61 were identified in DFa, 64 were identified in DFb and HFa, 63 were identified in HFb. There were 60 miRNAs in all libraries, 2 specifically expressed in HFa, and 1 specifically expressed in DFb (Fig. 9a). In addition, based on the structure of miRNAs, 259 novel miRNAs were predicted, of which, 241 were predicted in DFa, 235 in DFb, 238 in HFa and 239 in HFb. There were 212 novel miRNAs expressed in all libraries, 2 specifically expressed in DFb, HFa and HFb respectively, and 3 specifically expressed in DFa (Fig. 9b).

Fig. 9
figure9

The distribution of expressed miRNAs in each sample. a The distribution of known miRNAs. b The distribution of novel miRNAs

Differential expression analysis of miRNAs and target prediction of miRNAs

Based on differential expression analysis of miRNAs, 114 differentially expressed miRNAs (DEMs) were detected in the flower bud stage, of which 52 were up-regulated and 62 down-regulated in hexaploid S. canadensis. There were 3365 genes putatively targeted by multiple miRNAs in the flower bud stage. The GO functional classification revealed that these target genes were distributed in 43 GO terms. In the biology process category, most of DEGs were distributed in the subcategories of ‘cellular process’, ‘metabolic process’ and ‘biological regulation’ (Additional file 12: Fig. S8). Based on the GO enrichment analysis of these target genes and functional de-redundant analysis of the enriched GO terms, we found that many target genes were involved in some biological functions such as ‘response to biotic stimulus’, ‘thiol−dependent ubiquitinyl hydrolase activity’, ‘stamen development’ and ‘auxin influx transmembrane transporter activity’ (Fig. 10a).

Fig. 10
figure10

GO analysis of target genes of differentially expressed miRNAs. a The DFa-vs-HFa comparision; b The DFb-vs-HFb comparision

Furthermore, 163 DEMs were detected in the fruit development stage, of which 109 were up-regulated and 54 down-regulated. There were 4003 genes putatively targeted by multiple miRNAs in the fruit development stage. These target genes were also distributed in 43 GO terms and most of genes distributed into the same GO terms with flower stage (Additional file 13: Fig. S9). The GO enrichment analysis and functional de-redundant analysis revealed that many target genes were involved in ‘response to biotic stimulus’, ‘histamine N-methyltransferase activity’ and ‘coenzyme binding’ (Fig. 10b).

Visualization of miRNA-mRNA interaction network

Based on the negative regulatory role of miRNA in the process of mRNA expression, in this study, the interaction network between miRNA and target genes in the flower bud and fruit development stage of S. canadensis were constructed. In the flower bud stage, 36 DEMs were detected to negatively regulate 383 target DEGs, of which 17 down-regulated miRNAs negatively regulate 36 target genes (Fig. 11a). In the fruit development stage, 56 DEMs were detected to negatively regulate 490 target DEGs, of which 17 down-regulated miRNAs negatively regulate 40 target genes (Fig. 11b). These target genes include SPL genes which involved in promoting growth and development, regulating flowering time and seed yield, TF genes and some enzyme-related genes. In the flower bud stage, some miRNAs were putative to target functional genes, such as miR156b_2 targeted with SPL12 (CL2806.Contig2_All) and SPL13A (CL6003.Contig2_All) gene, miR172a_2 targeted with ERF (Unigene44257_All) TF genes, novel_mir50 targeted with MYB (CL824.Contig1_All) and zinc-finger protein (Unigene77225_All) TF gene. In the fruit development stage, some miRNAs were putative to target the same functional genes with that of flower bud stage, furthermore, these miRNAs were also have specific regulation effects. For example, except for targeting SPL12 and SPL13 genes, miR156 also targeted with SPL10 (Unigene40563_All) gene. In addition, the specific differentially expressed miR156 and miR156a during this period targeted with SPL12 (CL2806.Contig4_All) gene, miR156f targeted with SPL10 (Unigene40563_All), SPL13B (CL6003.Contig2_All) and SPL18 (CL17383.Contig4_All) gene. These DEMs may play a vital role in the flower bud and fruit development stage in hexaploid S. canadensis, and these miRNAs have differences regulation in the different stages.

Fig. 11
figure11

The interaction network of miRNAs with target genes. a The DFa-vs-HFa comparision; b The DFb-vs-HFb comparision

Confirmation the expression profile data of mRNAs and miRNAs by the quantitative real-time PCR (qRT-PCR)

To verify the accuracy of the expression of mRNA and miRNA obtained from mRNA-seq and miRNA-seq, we randomly selected 12 genes (Additional file 14: Fig. S10) and 7 miRNAs (Additional file 15: Fig. S11) for the qRT-PCR analyses. The relative expression levels of selected genes and miRNAs in each groups, which determined by qRT-PCR, were consistent with RNA-seq. The primer sequences used in this study were listed in the Additional file 16: Table S5.

Discussion

S. canadensis spread rapidly in the invasion area owing to its great capacity for sexual reproduction [47]. Hexaploid S. canadensis have developed into an invasive plant may related to polyploidization [17], which induced the gene expression level alteration [48]. This study depicted mRNA and miRNA integrated analyses in the flower bud and fruit development stage of the S. canadensis. The results showed that hexaploid S. canadensis has significant differences in the regulation of gene expression related to hormone signal transduction, reproductive development and the synthesis of secondary metabolites. These differentially expressed genes in hexaploid S. canadensis may provide a certain advantage in sexual reproduction, and lay a foundation for its rapid spread and formation of invasive plants.

Plant circadian rhythm related genes may affect the growth and development in hexaploid S. canadensis

Plants have capacities to adapt its survival environment by its internal circadian rhythm and regulate a series of physiological processes such as photoperiod, flowering time, hormone signal transduction, plant growth and metabolism, and response to biological and abiotic stresses. Studies have shown that CIRCADIAN CLOCK ASSOCIATED 1 (CCA1) protein can inhibit the expression of phytochrome B (PHYB) activation-tagged suppressor 1 (BAS1) gene to regulate the synthesis of brassinosteroids [49]. Cryptochrome 1 (CRY1) and phyB genes can directly interacted with auxin response factor (ARF) gene to regulate auxin signaling pathway [50]. The circadian rhythm system will affect flowering time in plants by sensing altered surrounding environment conditions (such as light and temperature) [51, 52]. GI and CO act as the key genes in the photoperiodic pathway, besides, the interaction of phyA, CRY1 and CRY2 genes will affected the expression of GI gene, which promote the expression of CO gene. The CO gene coding proteins which activate the FT gene, and then affected the expression of LFY and AP1 genes. These genes constructed a series of complex networks that will eventually promote the floral organs development. Polyploidy can induced a series of related gene expression regulation changes, and then significantly affected the plant morphological vigor [53]. In this study, most of related genes were differentially expressed in flower bud stage. The key genes with a large number of differentially expressed genes, such as phyB and PIF3, were up-regulated in hexaploid. They may play an important role in the photoperiod signaling and downstream regulation of hexaploid. In addition, PRR7, PRR5 and CO genes involved in flower development were also up-regulated in hexaploid. These genes may play a certain role in the growth and development of hexaploid S. canadensis, and the advantages in reproductive level.

Flower and fruit development related genes may promote the reproductive process of hexaploid S. canadensis

MADS-box genes are widely distributed in plants and play an important role in the regulation of growth and development, particularly in the development of floral organs, flowering time regulation and fruit development and ripeness in angiosperms [54,55,56]. In Arabidopsis, AGL23 gene participate in the regulation of the development of female gametophytes and the formation of organelles during embryonic development [57]; AGL61 cooperated with AGL80 gene to participate in the differentiation of female gametophytes [58]; AGL62 gene can be involved in the formation of endosperm cells; SEP gene involved in the formation of floral organs [59]; PHE1 gene plays an important role in the seed development and nutrient storage [60]; AGL61 and AGL80 genes were also involved in the development of central cells and endosperm [61, 62]. In addition, AP1 gene will interferes with the specification of floral organs from common primordia into floral organs, and the true conversion of flowers into inflorescences in legumes [63]. MADS1 and MADS7 genes in orchids play a role in the development of stamens and ovary [64]. In this study, many MADS genes were differentially expressed between two cytotypes of S. canadensis in the flower bud and fruit development stage. In the flower bud stage, most of the up-regulated genes were AGL8 and AP1 genes, while some AGL genes, such as AGL12 and AGL16 were down-regulated in hexaploid S. canadensis. These different expression patterns in flower bud stage indicate that there were different regulatory modes between two cytotypes, and may be a certain bias in functional selection, which may causing in the formation and development of flowers. In fruit stage, most of MADS-box-related genes were up-regulated in hexaploid S. canadensis, which indicates that hexaploid S. canadensis may have a significant effect compared with diploid during fruit formation and development. The number of differentially expressed AGL8 and PMADS2 genes was decreased when compared with flower bud stage, while the number of genes such as AP1, AGL9 and AGL65 was not altered. This indicated that these related genes may have a temporal and spatial bias among the flower bud and fruit development stage in S. canadensis. For example, AGL8 and PMADS2 genes may play a more important role in the flower bud stage compared with the fruit development stage. While, some genes such as AG2, MADS15, AGL11, PHE1, and SEP1 were detected up-regulated in hexaploid fruit development stage. Some of these genes were specifically and differentially expressed in the hexaploid S. canadensis, such as AG2 and MADS15, which may be vital to the development of the fruit. In addition, plant seeds accumulate a large amount of storage protein to improve seed quality during maturation. In this study, a large number of seed storage protein genes were differentially expressed between two cytotypes and most of them were up-regulated in hexaploid S. canadensis, which was likely to improve seed quality and provide material basis for higher germination efficiency after dissemination.

miRNA may play a vital role in the reproductive development of hexaploid S. canadensis

MiRNAs in plants play an important role in various regulatory processes such as growth, development, and stress resistance. Altered the regulation of miRNAs with target genes may lead to phenotypic changes in plants. OsmiR156 directed regulate the expression of SPL14 gene and activates the expression of downstream genes, which affected the panicle branching and grain yield in rice [65, 66]. miR172 affected the flowering time and flower organ development in plants by regulating the expression of AP2 genes [67]. Blue light alters the expression of miR167 and the target genes of auxin response factor genes in Arabidopsis, which finally affected the process of plant growth and generative development [68]. Research on the miRNA of S. canadensis shown that many DEMs were detected in flower bud and fruit development stages. Some DEMs have the same trend of expression in flower bud and fruit development stages, for example, the expression level of miR156b_2, and miR167d_1 were down-regulated in both stages of hexaploid S. canadensis. While, the expression level of miR156k_2, miR167d-5p and miR168a were up-regulated. Some miRNAs owing opposite trends of expression level between flower bud and fruit development stage. For example, miR156a-3p was up-regulated in flower bud stage, however, it down-regulated in fruit development stage; miR166m_2 was down-regulated in flower bud stage, however, it up-regulated in fruit development stage. There were also some miRNAs that were specifically differentially expressed in flower bud stage, such as miR162-3p, miR167f-5p and miR396b, and fruit development stage, such as miR156f, miR156f-5p, miR166-3p and miR172c-3p. The number of DEMs was different between flower bud and fruit development stages, which indicated that miRNAs may play special roles in different periods of hexaploid S. canadensis. Some of DEMs were predicted to target genes which involved in sexual reproduction, such as, miR156, miR156a and miR156b_2. These results suggest that DEMs may be involved in the growth, development and reproduction of S. canadensis.

Conclusions

We combined transcriptome and miRNA sequence technologies to study the molecular mechanism of invasiveness in hexaploid S. canadensis under the sexual reproduction. In this study, based on the screening of differentially expressed genes and miRNAs in flower bud and fruit development stages of hexaploid S. canadensis, we found that many genes involved in the reproductive development process were up-regulated in hexaploid. In addition, some miRNAs involved in the expression regulation process may play an important role in the reproductive development of hexaploid S. canadensis. The result may lead to hexaploid S. canadensis owing higher yield and fruit quality in the process of sexual reproduction and higher germination rate of seeds which conductive to diffusion, faster propagation process and enhanced invasiveness.

Materials and methods

Plant materials, cDNA and small RNA library construction and sequencing

The inflorescence of hexaploid cytotype (2n = 2x = 54), which contained in flower bud (HFa) and fruit development (HFb) stages, were collected from Wuhan (30°32′N, 114°25′E), Hubei Province, China. The rhizomes of diploid cytotype (2n = 2x = 18) were transplant into the Wuhan University open-air garden from Kunming (24°55′N, 102°47′E), Yunnan Province, China [69]. The diploid were used for flower bud (DFa) and fruit development (DFb) preparation. Three replicates for each sample of hexaploid and diploid were harvested from three independent individuals. The fresh tissues between two cytotypes were collected and immediately frozen in liquid nitrogen and stored at − 80 °C for further investigation. S. canadensis as an invasive plant in China and the specimen has been deposited in many publicly available herbarium, such as Herbarium of Chen shan Botanical Garden (0007639), Wuhan Botanical Garden Herbarium, Chinese Academy of Sciences (0005503). The plant material used in the study was consistent with Xu et al. 2019, and the formal identification was undertaken by Xu.

Total RNAs were extracted using TRIzol Reagent and then treated with DNase I. The RNA quality was verified by Agilent 2100 Bioanalyzer (Agilent RNA 6000 Nano Kit) with RIN number > 7.0. The cDNA and small RNA library were constructed by the methods provided by Beijing Genomics Institute (BGI, Shenzhen, China).

De novo assembly and unigene annotation

All of the constructed libraries were sequenced on an Illumina Hiseq X ten (Illumina Inc., MI, USA) platform. The raw reads of sequencing which include low-quality, contained adaptors, high content of unknown bases (more than 5%) and low-quality bases (more than 20% of the bases with a quality score less than 15) were removed. The trinity (version: v2.0.6) software with parameters as follows: --min_contig_length 150 --CPU 8 --min_kmer_cov 5 --min_glue 5 --bfly_opts’-V 5 --edge-thr = 0.1 --stderr’ [70] were used to assemble the obtained clean reads. Then the TGICL (version: v2.0.6) software with parameters as follows: -l 40 -c 10 -v 25 -O′-repeat_stringency 0.95 -minmatch 35 -minscore 35′ [71] were used to further cluster the transcripts to remove the redundant Trinity-generated transcripts. Finally, the “All-Unigene” was obtained for subsequent analysis. The “All-Unigene” sequences were aligned with the Kyoto Encyclopedia of Genes and Genomes (KEGG) public databases [72] by Blast (version: v2.2.23) and Gene Ontology (GO) by Blast2GO (version: v2.5.0) software [73] with the default parameter respectively.

Transcription factor (TF)-encoding gene prediction

To identify candidate genes for TF, getorf (version: EMBOSS: 6.5.7.0) software with parameters as follows: -minsize 150 was used to detect the ORF of unigene [74], and then used hmmsearch (version: v3.0) with default parameter [75] to identify the ability of the TF gene family according to the characteristics described by the PlantTFDB database.

Quantification of gene expression level and analysis of differently expressed genes

The Bowtie2 (version: v2.2.5) software [76] was used to align clean reads to assembled “All-Unigene” with the parameters as follows: -q --phred64 --sensitive --dpad 0 --gbar 99,999,999 --mp 1,1 --np 1 --score-min L,0,-0.1 -p 16 -k 200, and then based on the fragments per kilobase of transcript per million mapped reads (FPKM) analysis, RSEM (version: v1.2.8) software [77] was used to normalize the gene expression levels of each sample with the default parameter. The DEGseq [78] was used to identify differentially expressed genes (DEGs) by the value of |log2 Ratio| > 1.00 and adjusted p-value < 0.001 which was corrected by FDR for comparisons between diploid and hexaploid cytotypes with three biological replicates. Based on the result of annotation, GO and KEGG enrichment analyses were performed by phyper function. The hypergeometric test with the threshold of FDR (Q-value < 0.05) to find the significantly enriched terms (or pathways) in DEGs compared with the whole background.

Weighted gene co-expression network analysis (WGCNA)

WGCNA was used in gene co-expression network identification, which has been widely applied in various biological contexts for gene expression studies, and can be used for finding highly correlated genes from clusters or modules [79]. All the DEGs of the flower bud and fruit development stage between two cytotypes were used to construct gene network by R package WGCNA. The clustered modules were collected numbers of genes, which were assembled by unique color. For each module which satisfied with correlation coefficient > 0.6 and p < 0.05 as the significant correlation. The co-expressed genes in these modules were extracted and visualized by Cytoscape 3.7.1 for network construction [80].

Identification of known and novel miRNAs in S. canadensis

The clean reads were obtained by removing low-quality contaminated reads and adaptors. The length of clean reads which range from 18 to 30 nt was chosen for further analysis. The unique reads were immediately used to search against the miRBase 22.0 database by using the BLASTn program to annotate the conserved miRNAs. In addition, the miRNA precursor can be characterized by its hairpin structure, which was used to predict novel miRNA. Here, we used the software miRA (V1.2.0) [81] to predicted novel miRNAs through exploring the secondary structures and subsequent precursor.

Differentially expressed miRNAs (DEMs) analysis

The expression profiles for small RNAs were calculated by using Transcripts Per Kilobase Million (TPM), which based on the formula as follows: Normalized expression = Actual miRNA count × 106/Total count of clean reads. After normalization, the calculated small RNA expression level data was used directly to perform the differential expression analysis between the samples by the DEGseq R package. The threshold which satisfied with |log2 Ratio| > 1.00 and adjusted p-value < 0.001 were acted as significant DEMs.

Target gene prediction

To find more accurate targeted genes of miRNAs, multiple types of software were used. We used psRobot (version: V1.2) software [82] with parameter as follows: -gl 17 -p 8 -gn 1 and TargerFinder (version: V1.0) software [83] with parameter as follows: -c 4 to predict miRNA targets. The function of these target genes were also detected by using similar method based on GO and KEGG database.

MiRNA-mRNA interaction network analysis

The annotation and GO enrichment analysis of target genes was performed to explore the biological and critical functions of DEMs. Based on the integrated analysis of DEMs and the target genes, the Cytoscape 3.2.0 was used to construct the miRNA-mRNA regulatory network.

Quantitative real-time PCR (qRT-PCR) validation

To validate the sequencing result, qRT-PCR was used to detect the expression patterns of DEGs and DEMs of each stage of diploid and hexaploid S. canadensis. The primer sequences of qRT-PCR were designed with Primer 5 software. The GAPDH gene was used as an internal control. The selected genes and miRNAs were calculated by the 2 Ct method for relative expression with the PCR reaction condition was set as first denaturation at 95 °C in 10 min, then followed by 40 cycles of denaturation at 95 °C for 10 s, annealing and extension at 60 °C for 30 S.

Availability of data and materials

All RNA-seq data associated with this study have been submitted to the NCBI Sequence Read Archive, (Accession number: PRJNA726048, https://dataview.ncbi.nlm.nih.gov/object/PRJNA726048?reviewer=j2krcdp2ntv9h5e7hhsr7e5qg5).

Abbreviations

DFa:

Flower bud stage in diploid

DFb:

Fruit development in diploid

HFa:

Flower bud stage in hexaploid

HFb:

Fruit development in hexaploid

DEG:

Differentially Expressed Gene

TF:

Transcription Factor

KEGG:

Kyoto Encyclopedia of Genes and Genomes

GO:

Gene Ontology

TPM:

Transcripts Per Kilobase Million

DEM:

Differentially Expressed miRNA

qRT-PCR:

Quantitative real-time PCR

References

  1. 1.

    Bellard C, Rysman JF, Leroy B, Claud C, Mace GM. A global picture of biological invasion threat on islands. Nat Ecol Evol. 2017;1(12):1862–9.

    PubMed  Article  PubMed Central  Google Scholar 

  2. 2.

    McCary MA, Mores R, Farfan MA, Wise DH. Invasive plants have different effects on trophic structure of green and brown food webs in terrestrial ecosystems: a meta-analysis. Ecol Lett. 2016;19(3):328–35.

    PubMed  Article  PubMed Central  Google Scholar 

  3. 3.

    Mazza G, Tricarico E, Genovesi P, Gherardi F. Biological invaders are threats to human health: an overview. Ethol Ecol Evol. 2014;26(2-3):112–29.

    Article  Google Scholar 

  4. 4.

    Kettenring KM, McCormick MK, Baron HM, Whigham DF. Mechanisms of Phragmites australis invasion: feedbacks among genetic diversity, nutrients, and sexual reproduction. J Appl Ecol. 2011;48(5):1305–13.

    Article  Google Scholar 

  5. 5.

    Moron D, Skorka P, Lenda M. Disappearing edge: the flowering period changes the distribution of insect pollinators in invasive goldenrod patches. Insect Conserv Diver. 2019;12(2):98–108.

    Article  Google Scholar 

  6. 6.

    Sokol NW, Kuebbing SE, Bradford MA. Impacts of an invasive plant are fundamentally altered by a co-occurring forest disturbance. Ecology. 2017;98(8):2133–44.

    PubMed  Article  PubMed Central  Google Scholar 

  7. 7.

    Callaway RM, Maron JL. What have exotic plant invasions taught us over the past 20 years? Trends Ecol Evol. 2006;21(7):369–74.

    PubMed  Article  PubMed Central  Google Scholar 

  8. 8.

    Liebhold AM, Brockerhoff EG, Kalisz S, Nunez MA, Wardle DA, Wingfield MJ. Biological invasions in forest ecosystems. Biol Invasions. 2017;19:3437–58.

    Article  Google Scholar 

  9. 9.

    Cao YS, Xiao YA, Zhang SS, Hu WH. Simulated warming enhances biological invasion of Solidago canadensis and Bidens frondosa by increasing reproductive investment and altering flowering phenology pattern. Sci Rep-Uk. 2018;8(1):16073.

    Article  CAS  Google Scholar 

  10. 10.

    Ayres DR, Zaremba K, Sloop CM, Strong DR. Sexual reproduction of cordgrass hybrids (Spartina foliosa x alterniflora) invading tidal marshes in San Francisco Bay. Divers Distrib. 2008;14(2):187–95.

    Article  Google Scholar 

  11. 11.

    Liu H, Lin Z, Zhang M, Qi X. Relative importance of sexual and asexual reproduction for range expansion of Spartina alterniflora in different tidal zones on Chinese coast. Estuar Coast Shelf S. 2017;185(5):22–30.

    CAS  Article  Google Scholar 

  12. 12.

    Chen XS, Li YF, Xie YH, Deng ZM, Li X, Li F, et al. Trade-off between allocation to reproductive ramets and rhizome buds in Carex brevicuspis populations along a small-scale elevational gradient. Sci Rep-Uk. 2015;5:12688.

    CAS  Article  Google Scholar 

  13. 13.

    Correia M, Montesinos D, French K, Rodriguez-Echeverria S. Evidence for enemy release and increased seed production and size for two invasive Australian acacias. J Ecol. 2016;104(5):1391–9.

    Article  Google Scholar 

  14. 14.

    Jesse LC, Nason JD, Obrycki JJ, Moloney KA. Quantifying the levels of sexual reproduction and clonal spread in the invasive plant, Rosa multiflora. Biol Invasions. 2010;12(6):1847–54.

    Article  Google Scholar 

  15. 15.

    Fialova M, Jandova M, Ohryzek J, Duchoslav M. Biology of the polyploid geophyte Allium oleraceum (Amaryllidaceae): variation in size, sexual and asexual reproduction and germination within and between tetra-, penta- and hexaploid cytotypes. Flora. 2014;209(7):312–24.

    Article  Google Scholar 

  16. 16.

    Nagy UD, Stranczinger S, Godi A, Weisz A, Rosche C, Suda J, et al. Does higher ploidy level increase the risk of invasion? A case study with two geo-cytotypes of Solidago gigantea Aiton (Asteraceae). J Plant Ecol. 2017;11(2):317–27.

    Article  Google Scholar 

  17. 17.

    te Beest M, Le Roux JJ, Richardson DM, Brysting AK, Suda J, Kubesova M, et al. The more the better? The role of polyploidy in facilitating plant invasions. Ann Bot. 2012;109(1):19–45.

    Article  Google Scholar 

  18. 18.

    Schlaepfer DR, Edwards PJ, Billeter R. Why only tetraploid Solidago gigantea (Asteraceae) became invasive: a common garden comparison of ploidy levels. Oecologia. 2010;163(3):661–73.

    PubMed  Article  PubMed Central  Google Scholar 

  19. 19.

    Wani GA, Shah MA, Reshi ZA, Dar MA. Polyploidy determines the stage of invasion: clues from Kashmir Himalayan aquatic flora. Acta Physiol Plant. 2018;40(3):1–12.

    CAS  Article  Google Scholar 

  20. 20.

    Čertner M, Sudová R, Weiser M, Suda J, Kolář FA-O. Ploidy-altered phenotype interacts with local environment and may enhance polyploid establishment in Knautia serpentinicola (Caprifoliaceae). New Phytol. 2019;221(2):1117–27.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  21. 21.

    Anderson LB. A study of some seedling characters and the effects of competition on seedlings in diploid and tetraploid red clover (Trifolium pratense L.). New Zeal J Agr Res. 1971;14(3):563–71.

    Article  Google Scholar 

  22. 22.

    Cheng JL, Yang XH, Xue LF, Yao BB, Lu H, Tian ZS, et al. Polyploidization contributes to evolution of competitive ability: a long term common garden study on the invasive Solidago canadensis in China. Oikos. 2020;129(5):700–13.

    Article  Google Scholar 

  23. 23.

    Sandra PM, Martin L, Sylvain G. Competitive ability of Capsella species with different mating systems and ploidy levels. Ann Bot. 2018;121:1257–64.

    Article  Google Scholar 

  24. 24.

    Sun Y, Müller-Schärer H, Schaffner U. Neighbour origin and ploidy level drive impact of an alien invasive plant species in a competitive environment. PLoS One. 2016;11(5):e0155712.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  25. 25.

    Voshall A, Moriyama EN. Next-generation transcriptome assembly and analysis: impact of ploidy. Methods. 2020;176:14–24.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  26. 26.

    Song QX, Ando A, Jiang N, Ikeda Y, Chen ZJ. Single-cell RNA-seq analysis reveals ploidy-dependent and cell-specific transcriptome changes in Arabidopsis female gametophytes. Genome Biol. 2020;21:178.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  27. 27.

    Coate JE, Doyle JJ. Quantifying whole transcriptome size, a prerequisite for understanding transcriptome evolution across species: an example from a plant allopolyploid. Genome Biolo Evol. 2010;2:534–46.

    Article  CAS  Google Scholar 

  28. 28.

    Wendel JF, Lisch D, Hu G, Mason AS. The long and short of doubling down: polyploidy, epigenetics, and the temporal dynamics of genome fractionation. Curr Opin Genet Dev. 2018;49:1–7.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  29. 29.

    Bottani S, Zabet NR, Wendel JF, Veitia RA. Gene expression dominance in allopolyploids: hypotheses and models. Trends Plant Sci. 2018;23(5):393–402.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  30. 30.

    Pablo AM, Seong WY, Javier P. Keep calm and carry on: miRNA biogenesis under stress. Plant J. 2019;99(5):1360–75.

    Google Scholar 

  31. 31.

    Song X, Yan X, Cao Y. MicroRNAs and their regulatory roles in plant-environment interactions. Annu Rev Plant Biol. 2019;70:489–525.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  32. 32.

    Curaba J, Talbot M, Li Z, Helliwell C. Over-expression of microRNA171 affects phase transitions and floral meristem determinancy in barley. BMC Plant Biol. 2013;13:6.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  33. 33.

    Ahmed W, Li R, Xia Y, Bai G, HM, Siddique K, et al. Comparative analysis of miRNA expression profiles between heat-tolerant and heat-sensitive genotypes of flowering chinese cabbage under heat stress using high-throughput sequencing. Genes. 2020;11:264.

    CAS  PubMed Central  Article  Google Scholar 

  34. 34.

    Glazinska P, Kulasek M, Glinkowski W, Wojciechowski W, Kosinski J. Integrated analysis of small RNA, transcriptome and degradome sequencing provides new insights into floral development and abscission in yellow lupine (Lupines luteus L.). Int J Mol Sci. 2019;20:5122.

    CAS  PubMed Central  Article  Google Scholar 

  35. 35.

    Silva GFF, Silva EM, Correa JPO, Vicente MH, Jiang N, Notini MM, et al. Tomato floral induction and flower development are orchestrated by the interplay between gibberellin and two unrelated microRNA-controlled modules. New Phytol. 2019;221:1328–44.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  36. 36.

    Waheed S, Zeng L. The critical role of miRNAs in regulation of flowering time and flower development. Genes. 2020;11:319.

    CAS  PubMed Central  Article  Google Scholar 

  37. 37.

    Wang H, Li Y, Chern M, Zhu Y, Zhang L-L, Lu J-H, et al. Suppression of rice miR168 improves yield, flowering time and immunity. Nat Plants. 2021;7(2):129–36.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  38. 38.

    Das A, Nigam D, Junaid A, Tribhuvan KU, Gaikwad K. Expressivity of the key genes associated with seed and pod development is highly regulated via lncRNAs and miRNAs in Pigeonpea. Sci Rep. 2019;9(1):18191.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  39. 39.

    Carbone F, Bruno L, Perrotta G, Bitonti MB, Chiappetta A. Identification of miRNAs involved in fruit ripening by deep sequencing of Olea europaea L. transcriptome. PLoS One. 2019;14(8):e0221460.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  40. 40.

    Liu XX, Luo XF, Luo KX, Liu YL, Pan T. Small RNA sequencing reveals dynamic microRNA expression of important nutrient metabolism during development of Camellia oleifera fruit. Int J Biol Sci. 2019;15(2):416–29.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  41. 41.

    Suo Y, Min Y, Dong C, Wang Y, Cheng S, Kang X. MicroRNA expression changes following synthesis of three full-sib Populus triploid populations with different heterozygosities. Plant Mol Biol. 2017;95:215–25.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  42. 42.

    Fu Y, Xiao M, Yu H, Mason AS, Yin J, Li J, et al. Small RNA changes in synthetic Brassica napus. Planta. 2016;244:607–22.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  43. 43.

    Shen Y, Zhao Q, Zou J, Wang W, Gao Y, Meng J, et al. Characterization and expression patterns of small RNAs in synthesized Brassica hexaploids. Plant Mol Biol. 2014;85:287–99.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  44. 44.

    Dong M, Lu BR, Zhang HB, Chen JK, Li B. Role of sexual reproduction in the spread of an invasive clonal plant Solidago canadensis revealed using intersimple sequence repeat markers. Plant Spec Biol. 2006;21(1):13–8.

    Article  Google Scholar 

  45. 45.

    Lu H, Xue LF, Cheng JL, Yang XH, Xie HJ, Song XL, et al. Polyploidization-driven differentiation of freezing tolerance in Solidago canadensis. Plant Cell Environ. 2020;43(6):1394–403.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  46. 46.

    Yang XH, Cheng JL, Yao BB, Lu H, Zhang Y, Xu JX, et al. Polyploidy-promoted phenolic metabolism confers the increased competitive ability of Solidago canadensis. Oikos. 2021;130(6):1014–25.

    CAS  Article  Google Scholar 

  47. 47.

    Hua H, Chen G, Guo S. Reproductive biology in an invasive plant Solidago canadensis. Front Biol China. 2007;2(2):196–204.

    Article  Google Scholar 

  48. 48.

    Wu M, Ge Y, Xu C, Wang J. Metabolome and transcriptome analysis of hexaploid Solidago canadensis roots reveals its invasive capacity related to polyploidy. Genes. 2020;11:187.

    CAS  PubMed Central  Article  Google Scholar 

  49. 49.

    Peng H, Neff MM. CIRCADIAN CLOCK ASSOCIATED 1 and ATAF2 differentially suppress cytochrome P450-mediated brassinosteroid inactivation. J Exp Bot. 2020;71(3):970–85.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  50. 50.

    Mao Z, He S, Xu F, Wei X, Jiang L, Liu Y, et al. Photoexcited CRY1 and PHYB interact directly with ARF6 and ARF8 to regulate their DNA-binding activity and auxin-induced hypocotyl elongation in Arabidopsis. New Phytol. 2020;225(2):848–65.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  51. 51.

    Blümel M, Dally N, Jung C. Flowering time regulation in crops—what did we learn from Arabidopsis? Curr Opin Biotech. 2015;32:121–9.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  52. 52.

    Frédéric B, Guillaume L, Pierre T, Claire P. FLOR-ID: an interactive database of flowering-time gene networks in Arabidopsis thaliana. Nucleic Acids Res. 2015;44:D1167–71.

    Google Scholar 

  53. 53.

    Ni Z, Kim ED, Ha M, Lackey E, Liu J, Zhang Y, et al. Altered circadian rhythms regulate growth vigour in hybrids and allopolyploids. Nature. 2009;457(7227):327–31.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  54. 54.

    Alvarez-Buylla ER, Pelaz S, Liljegren SJ, Gold SE, Burgeff CD, Gary S, et al. An ancestral MADS-box gene duplication occurred before the divergence of plants and animals. Proc Natl Acad Sci U S A. 2000;97:5328–33.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  55. 55.

    Saedler H, Becker A, Winter KU, Kirchner C, Theissen G. MADS-box genes are involved in floral development and evolution. Acta Biochim Pol. 2001;48(2):351–8.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  56. 56.

    Sheng XG, Zhao ZQ, Wang JS, Yu HF, Shen YS, Zeng XY, et al. Genome wide analysis of MADS-box gene family in Brassica oleracea reveals conservation and variation in flower development. BMC Plant Biol. 2019;19:106.

    PubMed  PubMed Central  Article  Google Scholar 

  57. 57.

    Colombo M, Masiero S, Vanzulli S, Lardelli P, Kater MM, Colombo L. AGL23, a type I MADS-box gene that controls female gametophyte and embryo development in Arabidopsis. Plant J. 2008;54(6):1037–48.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  58. 58.

    Bemer M, Wolters-Arts M, Grossniklaus U, Angenent GC. The MADS domain protein DIANA acts together with AGAMOUS-LIKE80 to specify the central cell in Arabidopsis ovules. Plant Cell. 2008;20(8):2088–101.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  59. 59.

    Pelaz S, Ditta GS, Baumann E, Wisman E, Yanofsky MF. B and C floral organ identity functions require SEPALLATA MADS-box genes. Nature. 2000;405(6783):200–3.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  60. 60.

    Fatihi A, Zbierzak AM, Dormann P. Alterations in seed development gene expression affect size and oil content of Arabidopsis seeds. Plant Physiol. 2013;163(2):973–85.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  61. 61.

    Portereiko MF, Lloyd A, Steffen JG, Punwani JA, Otsuga D, Drews GN. AGL80 is required for central cell and endosperm development in Arabidopsis. Plant Cell. 2006;18(8):1862–72.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  62. 62.

    Steffen JG, Kang I-H, Portereiko MF, Lloyd A, Drews GN. AGL61 interacts with AGL80 and is required for central cell development in Arabidopsis. Plant Physiol. 2008;148(1):259–68.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  63. 63.

    Benlloch R, d'Erfurth I, Ferrándiz C, Cosson V, Beltrán J, Cañas L, et al. Isolation of mtpim proves Tnt1 a useful reverse genetics tool in Medicago truncatula and uncovers new aspects of AP1-like functions in legumes. Plant Physiol. 2006;142:972–83.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  64. 64.

    Chen YY, Lee PF, Hsiao YY, Wu WL, Pan ZJ, Lee YI, et al. C- and D-class MADS-box genes from Phalaenopsis equestris (Orchidaceae) display functions in gynostemium and ovule development. Plant Cell Physiol. 2012;53(6):1053–67.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  65. 65.

    Jiao Y, Wang Y, Xue D, Wang J, Yan M, Liu G, et al. Regulation of OsSPL14 by OsmiR156 defines ideal plant architecture in rice. Nat Genet. 2010;42(6):541–4.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  66. 66.

    Miura K, Ikeda M, Matsubara A, Song XJ, Ito M, Asano K, et al. OsSPL14 promotes panicle branching and higher grain productivity in rice. Nat Genet. 2010;42(6):545–9.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  67. 67.

    Qian-Hao Z, Helliwell CA. Regulation of flowering time and floral patterning by miR172. J Exp Bot. 2010;62(2):487–95.

    Google Scholar 

  68. 68.

    Pashkovskiy PP, Kartashov AV, Zlobin IE, Pogosyan SI, Kuznetsov VV. Blue light alters miR167 expression and microRNA-targeted auxin response factor genes in Arabidopsis thaliana plants. Plant Physiol Bioch. 2016;104:146–54.

    CAS  Article  Google Scholar 

  69. 69.

    Xu CC, Ge YM, Wang JB. Molecular basis underlying the successful invasion of hexaploid cytotypes of Solidago canadensis L.: insights from integrated gene and miRNA expression profiling. Ecol Evol. 2019;9(8):4820–52.

    PubMed  PubMed Central  Article  Google Scholar 

  70. 70.

    Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, et al. Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat Biotechnol. 2011;29(7):644–52.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  71. 71.

    Pertea G, Huang X, Liang F, Antonescu V, Sultana R. TIGR gene indices clustering tools (TGICL): a software system for fast clustering of large EST datasets. Bioinformatics. 2003;19(5):651–2.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  72. 72.

    Kanehisa M, Goto S. KEGG: Kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 2000;28(1):27–30.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  73. 73.

    Conesa A, Gotz S, Garcia-Gomez JM, Terol J, Talon M, Robles M. Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics. 2005;21(18):3674–6.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  74. 74.

    Rice P, Longden I, Bleasby A. EMBOSS: the European molecular biology open software suite. Trends Genet. 2000;16(6):276–7.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

  75. 75.

    Mistry J, Finn RD, Eddy SR, Bateman A, Punta M. Challenges in homology search: HMMER3 and convergent evolution of coiled-coil regions. Nucleic Acids Res. 2013;41(12):e121.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  76. 76.

    Langmead B, Salzberg SL. Fast gapped-read alignment with bowtie 2. Nat Methods. 2012;9(4):357–9.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  77. 77.

    Li B, Dewey CN. RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinformatics. 2011;12(1):93–9.

    CAS  Article  Google Scholar 

  78. 78.

    Wang L, Feng Z, Wang X, Wang X, Zhang X. DEGseq: an R package for identifying differentially expressed genes from RNA-seq data. Bioinformatics. 2010;26(1):136–8.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  79. 79.

    Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008;9:559.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  80. 80.

    Saito R, Smoot M, Ono K, Ruscheinski J, Wang PL, Lotia R, et al. A travel guide to Cytoscape plugins. Nat Methods. 2012;9:1069–76.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  81. 81.

    Evers M, Huttner M, Dueck A, Meister G, Engelmann JC. miRA: adaptable novel miRNA identification in plants using small RNA sequencing data. BMC Bioinformatics. 2015;16:370.

    PubMed  PubMed Central  Article  Google Scholar 

  82. 82.

    Wu HJ, Ma YK, Chen T, Wang M, Wang XJ. PsRobot: a web-based plant small RNA meta-analysis toolbox. Nucleic Acids Res. 2012;40:W22–8.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  83. 83.

    Fahlgren N, Carrington JC. miRNA target prediction in plants. Methods Mol Biol. 2010;592:51–7.

    CAS  PubMed  Article  PubMed Central  Google Scholar 

Download references

Acknowledgements

Not applicable.

Funding

This work was supported by the Science and technology research project of Henan Province (agricultural field): Research on preparation, function and health food of bioactive peptide from Dioscorea opposita Thunb (182102110093). The funding providers was not involved in designing the study; collecting, analyzing or interpreting the data; or deciding to submit the paper for publication, but provide basic expenses during the experiment.

Author information

Affiliations

Authors

Contributions

MW conceived and designed the study. MW and HL performed the experimental analyses. MW wrote the manuscript and responsible for planting and collecting materials. BL and TZ provided some constructive suggestions. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Miao Wu.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Wu, M., Liu, H., Li, B. et al. Integrated analysis of mRNA-seq and miRNA-seq reveals the advantage of polyploid Solidago canadensis in sexual reproduction. BMC Plant Biol 21, 462 (2021). https://0-doi-org.brum.beds.ac.uk/10.1186/s12870-021-03240-x

Download citation

Keywords

  • Polyploidy
  • Solidago canadensis
  • Plant invasion
  • Sexual reproduction
  • Transcriptome
  • miRNA