# MicroRNA and Putative Target Discoveries in Chrysanthemum Polyploidy Breeding

Authors: Fengjiao Zhang; Jingya Zhao; Sujuan Xu; Weimin Fang; Fadi Chen; Nianjun Teng
Journal: International Journal of Genomics
 (2017)

---

## Abstract

MicroRNAs (miRNAs), around 22 nucleotides (nt) in length, are a class of endogenous and noncoding RNA molecule that play an essential role in plant development, either by suppressing the transcription of target genes at a transcriptional level or inhibiting translation at a posttranscriptional level. To understand the roles of miRNAs and their target genes in chrysanthemum polyploidy breeding, three sRNA libraries of normal and abnormal embryos after hybridization were performed by RNA-Seq. As a result, a total of 170 miRNAs were identified and there are 41 special miRNAs in cross of paternal chromosome doubling, such as miR169b, miR440, and miR528-5p. miR164c and miR159a were highly expressed in a normal embryo at 18 days after pollination, suggesting the regulatory role at the late stage of embryonic development. miR172c was only detected in the normal embryo at 18 days after pollination, which means that miR172c mainly mediates gene expression in postembryonic development and these genes may promote embryo maturation. Other miRNAs, including miR414, miR2661, and miR5021, may regulate the genes participated in pathways of auxin response and energy metabolism; then they regulate the complex embryonic development together.


---

## Body

## 1. Introduction

In plants, microRNAs (miRNAs) are a major class of small noncoding RNAs (sRNAs) with 20–22 nt. sRNAs have been identified to control the developmental processes in plants by regulating gene expression [1]. They have the potential to regulate a gene by two ways: (1) posttranscriptional gene silencing (PTGS) by binding to 3′ untranslated region (UTR) of messenger RNAs (mRNAs) and repressing the translation of target mRNAs; (2) transcriptional gene silencing (TGS) by epigenetic modifications [2, 3]. PTGS is the main strategy used by miRNAs to regulate gene expression in plant development. Plant mature miRNAs are generated from miRNA precursors that are processed by a ribonuclease DICER-LIKE1, which negatively regulates specific target mRNAs [2].A number of studies have investigated the key regulatory role of miRNAs in a wide range of growth and development processes in plants, including the regulation of embryonic development [4, 5]. Plant embryo development includes two stages, embryo morphogenesis and seed maturation. In the model plant Arabidopsis thaliana, embryonic pattern formation has been mainly concerned and studied [5]. Based on the embryo shape, the cell division goes through several stages from preglobular to mature embryo, tightly regulated by multiple genes [6]. As a class of small regulatory RNA, the function of miRNAs during plant embryo development has been considerably reported recently [7]. 28,645 mature miRNAs have been discovered and deposited in the public miRNA database miRBase (Release 21, http://www.mirbase.org/) and hundreds of miRNA target in plant embryo. A large number of Arabidopsis mutants of miRNA biogenesis genes have revealed the crucial roles of miRNA during seed morphogenesis and maturation. Willmann et al. reported the earlier timing of embryo maturation in Arabidopsis mutant for strong alleles of DCL1 (DICER-LIKE1) that are required for miRNA biogenesis and demonstrated the negative regulatory role of specific miRNAs during early embryogenesis and later in embryonic development and that miRNAs are key regulators during seed maturation program [8].Embryonic miRNAs mediate plant embryo development by regulating transcription factors located in special spatial organization and other key developmental regulators [5]. In plants, miR165/166 is one of the best characterized miRNA families, which regulates five class III homeodomain leucine zipper genes (HD-ZIP IIIs) (PHB, PHV, REV, CNA, and ATHB8) [9, 10]. The HD-Zip III gene family regulates apical embryo patterning and organ polarity as well as controls shoot and root apical meristem (SAM and RAM) formation in embryogenesis [11–13]. The family of miR160/miR167 regulates the target mRNAs of auxin response factors (ARFs) associated with auxin homeostasis [14, 15]. Auxin response is an important signaling pathway during embryonic pattern formation, embryo development, and seed maturation [16].Chrysanthemum (Chrysanthemum morifolium) is an economically important flower around the world, with the increase of chrysanthemum consumption, breeders are driven to improve the cultivars’ traits, such as color, size, shape, and tolerance [17]. Artificial distant hybridization is one of the most effective methods to improve and create new cultivars. However, embryo abortion commonly happens during hybridization [18, 19]. Although previous studies analyzed many different reasons leading to plant embryo abortion, including maternal genotypes [20], parent ploidy [21, 22], and gene regulation [6, 23], the molecular mechanisms regulating embryo development are poorly understood. Zhang et al. revealed the gene, protein, and miRNA change in the stage of chrysanthemum embryo development [19, 24] and provided numerous information of embryo abortion in chrysanthemum hybridization breeding. Evidence supporting a major role for chromosome doubling in overcoming chrysanthemum embryo abortion has been obtained from some studies [22, 25], but the functions of miRNA remain largely uncharacterized.Next generation sequencing (NGS) approaches have been employed to identify individual miRNAs in various samples, and bioinformatics analyses have offered the technical support to predict the miRNA targets [26]. In the present study, we sequenced three small RNA libraries from chrysanthemum embryo in cross C. morifolium “Yuhualuoying” × tetraploid C. nankingense and identified the key miRNAs and targets that may facilitate embryo development in cross of paternal chromosome doubling.


## 2. Materials and Methods

### 2.1. Plant Materials and Artificial Hybridization

Artificial hybridization was performed in crossC. morifolium “Yuhualuoying” (♀, 2n = 6X = 54) × tetraploid C. nankingense (♂, 2n = 4X = 36). Here, the male is an autopolyploid generated by colchicine doubling of the diploid C. nankingense (2n = 2X = 18) [19]. After pollination, three samples were collected, corresponding to the normal embryo at 12 days after pollination (DAP) (NE12), normal embryo at 18 DAP (NE18), and abnormal embryo at 18 DAP (AE18). For each sample, we collected 0.2 g independent biological replicates for verification and mixed the rest of triplicate samples (~0.5 g) for RNA-Seq. All of the samples were immediately frozen in liquid nitrogen and stored at −80°C.


### 2.2. RNA Extraction, Small RNA Library Preparation, and Illumina Sequencing

Total RNA was extracted with TRIzol reagent (Takara Bio Inc., Otsu, Japan) according to the manufacturer’s protocol. Small RNAs with 18–30 nt fragments were enriched by 15% denaturing polyacrylamide gel electrophoresis. After purification, they were ligated to 5′ and 3′ adaptors and reversed transcribed into cDNA by reverse transcription-PCR (RT-PCR). Three small RNA libraries were constructed and sequenced using the Illumina HiSeq™ 2000 by the Beijing Genomics Institute (BGI) (Shenzhen, Guangdong Province, China).


### 2.3. Conserved and Novel miRNA Prediction

Firstly, the data cleaning analysis was performed by getting rid of low-quality reads, reads with 5′ primer contaminants and poly A, reads without 3′ primer and the insert tag, and reads shorter than 18 nt. The small RNA tags with miRNA, rRNA, snRNA, snoRNA, and tRNA were annotated by aligning to GenBank (http://www.ncbi.nlm.nih.gov/genbank/) and Rfam database (http://rfam.xfam.org/) using all clean reads of 18~30 nt. The number and proportion of each type of sRNAs were calculated in three libraries. Then, to identify the miRNAs in chrysanthemum embryo, miRBase 19.0 (http://www.mirbase.org/) was used to search the conserved miRNAs by BlantN. Only 90% matched sequences were considered to be conserved miRNAs. To allow the unambiguous mapping of small RNAs to annotations, the priority rule was followed: rRNA (in which GenBank > Rfam) > known miRNA > repeat > exon > intron. Finally, the novel miRNAs were predicted using the Mireap software (https://sourceforge.net/projects/mireap/); here, the chrysanthemum embryo transcriptome library obtained from the same sample with the present study (the NCBI accession number PRJNA315793) was used as the reference database. In the chrysanthemum embryo, a total of 99,119 unigenes were assembled with a mean length of 550–580 nt, which were used for prediction in the present study.


### 2.4. Target Prediction and Functional Annotation for miRNAs

The potential target genes of the known miRNAs were predicted by the web tool psRNATarget (http://plantgrn.noble.org/psRNATarget/) with parameters suggested by Allen et al. [27]. The target genes were identified in chrysanthemum embryo transcriptome dataset [22], and the function of these potential target genes was annotated using the two protein databases, Gene Ontology (GO) (http://geneontology.org/) and Kyoto Encyclopedia of Genes and Genomes (KEGG) (http://www.genome.jp/kegg/).


### 2.5. Differential Expression of Known miRNA and Their Targets

To find out the differentially expressed miRNAs in three samples, we normalized the expression of miRNA in three samples and obtained the expression of transcript per million (TPM). Normalized expression = actual miRNA count/total count of clean reads∗1000000. If the miRNA count was zero, it was revised to 0.01 for analysis of differential expression. The fold-change of log2 (sample 1/sample 2) and p value from the normalized counts were calculated to determine significant expression changes. Finally, those miRNAs with fold-change > 1 and p value <0.05 were considered to be differentially expressed in the two samples. Heatmap is a visual tool reflecting the expression of differential miRNAs. In the present study, the OmicShare tool, a free online platform for data analysis (http://www.omicshare.com/tools), was plotted. Based on the transcriptome library, we searched the expression pattern of these target genes regulated by the differentially expressed miRNAs.


### 2.6. Real-Time Quantitative PCR (qRT-PCR) Validation of miRNAs and Target Genes

We randomly selected 12 miRNAs and 9 target genes to validate the reliability of sRNA sequencing by qRT-PCR. RNA samples used for sequencing were reverse transcribed using PrimeScript miRNA qPCR starter kit ver 2.0 (Takara Bio, Dalian, China). qRT-PCR was performed using the SYBR Premix EX Taq Kit (Takara, Dalian, China), and the PCR amplification was done as described by Zhang et al. [24]. Three biological replicates were performed for each sample, and relative expression levels were calculated by the 2−△△CT method. The chrysanthemum gene EF1α (elongation factor 1a) (GenBank accession number KF305681) was used as a reference, which is stably expressed in chrysanthemum [19]. Special primers were designed using PRIMER3 RELEASE 2.3.4. All of the primers were shown in Table S1.


## 2.1. Plant Materials and Artificial Hybridization

Artificial hybridization was performed in crossC. morifolium “Yuhualuoying” (♀, 2n = 6X = 54) × tetraploid C. nankingense (♂, 2n = 4X = 36). Here, the male is an autopolyploid generated by colchicine doubling of the diploid C. nankingense (2n = 2X = 18) [19]. After pollination, three samples were collected, corresponding to the normal embryo at 12 days after pollination (DAP) (NE12), normal embryo at 18 DAP (NE18), and abnormal embryo at 18 DAP (AE18). For each sample, we collected 0.2 g independent biological replicates for verification and mixed the rest of triplicate samples (~0.5 g) for RNA-Seq. All of the samples were immediately frozen in liquid nitrogen and stored at −80°C.


## 2.2. RNA Extraction, Small RNA Library Preparation, and Illumina Sequencing

Total RNA was extracted with TRIzol reagent (Takara Bio Inc., Otsu, Japan) according to the manufacturer’s protocol. Small RNAs with 18–30 nt fragments were enriched by 15% denaturing polyacrylamide gel electrophoresis. After purification, they were ligated to 5′ and 3′ adaptors and reversed transcribed into cDNA by reverse transcription-PCR (RT-PCR). Three small RNA libraries were constructed and sequenced using the Illumina HiSeq™ 2000 by the Beijing Genomics Institute (BGI) (Shenzhen, Guangdong Province, China).


## 2.3. Conserved and Novel miRNA Prediction

Firstly, the data cleaning analysis was performed by getting rid of low-quality reads, reads with 5′ primer contaminants and poly A, reads without 3′ primer and the insert tag, and reads shorter than 18 nt. The small RNA tags with miRNA, rRNA, snRNA, snoRNA, and tRNA were annotated by aligning to GenBank (http://www.ncbi.nlm.nih.gov/genbank/) and Rfam database (http://rfam.xfam.org/) using all clean reads of 18~30 nt. The number and proportion of each type of sRNAs were calculated in three libraries. Then, to identify the miRNAs in chrysanthemum embryo, miRBase 19.0 (http://www.mirbase.org/) was used to search the conserved miRNAs by BlantN. Only 90% matched sequences were considered to be conserved miRNAs. To allow the unambiguous mapping of small RNAs to annotations, the priority rule was followed: rRNA (in which GenBank > Rfam) > known miRNA > repeat > exon > intron. Finally, the novel miRNAs were predicted using the Mireap software (https://sourceforge.net/projects/mireap/); here, the chrysanthemum embryo transcriptome library obtained from the same sample with the present study (the NCBI accession number PRJNA315793) was used as the reference database. In the chrysanthemum embryo, a total of 99,119 unigenes were assembled with a mean length of 550–580 nt, which were used for prediction in the present study.


## 2.4. Target Prediction and Functional Annotation for miRNAs

The potential target genes of the known miRNAs were predicted by the web tool psRNATarget (http://plantgrn.noble.org/psRNATarget/) with parameters suggested by Allen et al. [27]. The target genes were identified in chrysanthemum embryo transcriptome dataset [22], and the function of these potential target genes was annotated using the two protein databases, Gene Ontology (GO) (http://geneontology.org/) and Kyoto Encyclopedia of Genes and Genomes (KEGG) (http://www.genome.jp/kegg/).


## 2.5. Differential Expression of Known miRNA and Their Targets

To find out the differentially expressed miRNAs in three samples, we normalized the expression of miRNA in three samples and obtained the expression of transcript per million (TPM). Normalized expression = actual miRNA count/total count of clean reads∗1000000. If the miRNA count was zero, it was revised to 0.01 for analysis of differential expression. The fold-change of log2 (sample 1/sample 2) and p value from the normalized counts were calculated to determine significant expression changes. Finally, those miRNAs with fold-change > 1 and p value <0.05 were considered to be differentially expressed in the two samples. Heatmap is a visual tool reflecting the expression of differential miRNAs. In the present study, the OmicShare tool, a free online platform for data analysis (http://www.omicshare.com/tools), was plotted. Based on the transcriptome library, we searched the expression pattern of these target genes regulated by the differentially expressed miRNAs.


## 2.6. Real-Time Quantitative PCR (qRT-PCR) Validation of miRNAs and Target Genes

We randomly selected 12 miRNAs and 9 target genes to validate the reliability of sRNA sequencing by qRT-PCR. RNA samples used for sequencing were reverse transcribed using PrimeScript miRNA qPCR starter kit ver 2.0 (Takara Bio, Dalian, China). qRT-PCR was performed using the SYBR Premix EX Taq Kit (Takara, Dalian, China), and the PCR amplification was done as described by Zhang et al. [24]. Three biological replicates were performed for each sample, and relative expression levels were calculated by the 2−△△CT method. The chrysanthemum gene EF1α (elongation factor 1a) (GenBank accession number KF305681) was used as a reference, which is stably expressed in chrysanthemum [19]. Special primers were designed using PRIMER3 RELEASE 2.3.4. All of the primers were shown in Table S1.


## 3. Results

### 3.1. Sequencing Analysis of sRNA

To explore the regulation of miRNAs in chrysanthemum embryo development when paternal chromosomes were doubled, three cDNA libraries of small RNAs were constructed, which were named NE12, NE18, and AE18. All clean reads were obtained by filtering the low-quality sequences, adapter sequences, and poly-A sequences shorter than 18 nt, which altogether resulted in more than 99.79% of raw reads. When these small RNA tags were mapped to genome, about 20% of reads in each sample were matched. The unique sRNAs matched to genome in the three libraries were 3,864,037, 3,780,667, and 3,734,164, accounting for 9.62%, 10.04%, and 10.29% of all unique sRNAs (Table1).Table 1
Quantity and type distribution of small RNAs in chrysanthemum embryo.
NE12
NE18
AE18
Category
Unique sRNAs
Percent (%)
Total sRNAs
Percent (%)
Unique sRNAs
Percent (%)
Total sRNAs
Percent (%)
Unique sRNAs
Percent (%)
Total sRNAs
Percent (%)
Total
3,864,037
100
11,535,216
100
3,780,667
100
10,805,231
100
3,734,164
100
11,220,674
100
Mapping to genome
371,560
9.62
2,251,582
19.52
379,579
10.04
2,194,212
20.31
384,290
10.29
2,383,150
21.24
miRNA
9900
0.26
410,884
3.56
9192
0.24
343,417
3.18
9871
0.26
360,354
3.21
rRNA
57,663
1.49
547,936
4.75
68,080
1.80
662,520
6.13
81,807
2.19
852,436
7.60
snRNA
1423
0.04
3088
0.03
1475
0.04
3034
0.03
1564
0.04
3278
0.03
snoRNA
421
0.01
757
0.01
461
0.01
789
0.01
613
0.02
1205
0.01
tRNA
7221
0.19
72,063
0.62
7332
0.19
97,766
0.90
9506
0.25
118,373
1.05
Unannotated
3,787,409
98.02
10,500,488
91.03
3,694,127
97.71
9,697,705
89.75
3,630,803
97.23
9,885,028
88.10When these unique sRNAs were aligned to the GenBank and miRBase, many types of sRNAs were identified, including miRNA, rRNA, snRNA, snoRNA, and tRNA, but the vast majority of sequences were unannotated. The length of different types of sRNAs was discrepant; since in chrysanthemum embryo, the most abundant classes of sRNA showed the length of 24 nt (dominant siRNA), then 21 nt (mainly miRNAs), and 22 nt (Figure1). In the present study, the unique miRNAs in each library were taken into consideration for subsequent analysis, and the proportion of miRNAs was 0.26%, 0.24%, and 0.26% in libraries of NE12, NE18, and AE18 (Table 1).Figure 1
Length distribution of sRNAs in chrysanthemum embryo.x-axis is the length of sRNA distribution, and y-axis is the proportion of the sRNAs of different lengths. NE12, NE18, and AE18 mean the normal embryo at 12 DAP (days after pollination), normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.


### 3.2. miRNAs and Target Genes Identified in Three Libraries

In all samples, a total of 170 conversed miRNAs were identified in miRBase, and 130, 131, and 132 miRNAs were expressed in NE12, NE18, and AE18, respectively. 100 miRNAs (accounting for 58.8%) were detected in three samples; however, some were in two different samples, such as 10 miRNAs in NE12 and NE18, 7 miRNAs in NE18 and AE18, and 6 miRNAs in NE12 and AE18. Venn diagram (Figure2) presented the quantity distribution of conserved miRNAs in chrysanthemum embryo.Figure 2
Venn diagram presented the quantity distribution of conserved miRNAs in chrysanthemum embryos. There were 130, 131, and 132 miRNAs expressed in NE12, NE18, and AE18, respectively. In the middle, 100 means that the number of miRNAs was detected in three samples. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.Expression level is an important feature to explain the regulation function of miRNAs, and they commonly varied greatly in different samples. Here, the expression of 170 miRNAs with sequences was normalized and analyzed. The most abundant miRNAs in three samples were miR156a, miR157a, and miR166a, with the expression more than a thousand (TableS2). However, some miRNAs highly expressed in a particular sample, such as the expression of miR398b-5p, are 3524 in AE18 but neither in NE12 nor NE18, miR5721 only expressed in NE18, and miR5662 only in NE12.To study the biological function of miRNA in chrysanthemum embryo development, the sequences of target mRNA were paid attention. A total of 770 target genes were identified and regulated by 88 miRNAs (51.76% of all miRNAs), and miRNA414 had the most targets (347 unigenes), followed were miR5293 (61 targets) and miR5021 (44 targets) (TableS3). 34 target genes were regulated by two miRNAs, such as CL1999. Contig2 was targeted by miR165a and miR166a; unigene10412 was the target of miR156a and miR157a (Table S3).


### 3.3. Functional Annotation of Target Genes

Gene Ontology (GO) is a standardized classification system of gene used to describe the characteristics of genes and gene products in organisms. The result showed that 517 target genes were classified into three gene ontology categories: biological processes, cellular components, and molecular functions (Figure3). Of the 40 functional categories, the number of genes in each sample is different, but in most of categories, NE12 had the largest number but least in AE18. In the third ontology of molecular function, some target genes were expressed in a specific sample, such as the genes related to “electron carrier activity” only expressed in NE12, two genes in the category of “enzyme regulator activity” expressed in NE18, but there was no gene regulated “molecular transducer activity” in NE18 (Figure 3).Figure 3
GO classification of target genes in three samples.x-axis means three categories of GO: biological process, cellular component, and molecular function. y-axis is the number of genes in three samples associated with each subcategory. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.To understand the biological function of target genes, KEGG analysis, as with GO, was also used to analyze candidate targets. In NE12, a total of 337 target genes had the biological function on 207 KEGG pathways, but there were less target genes and pathways in NE18 and AE18. 211 target genes with 118 pathway annotation in NE18 and 179 annotated target genes with 113 pathways in AE18 (TableS4). Some target genes were annotated in dozens of pathways that are only in library of NE12, including some energy metabolism pathways, such as ko00280 (valine, leucine, and isoleucine degradation), ko00310 (lysine degradation), and ko00410 (beta-alanine metabolism).


### 3.4. Differentially Expressed miRNAs in Chrysanthemum Embryo

The aim of this study is to analyze the different miRNAs possibly involved in chrysanthemum embryo development, so we identified 112 differentially expressed miRNAs with at least 1.5-fold after standardized expression. The heatmap showed the expression pattern of these miRNAs (Figure4). They were assembled in three groups depending on the expression trend. In group NE18/AE18, the number of upregulated miRNAs was almost the same as downregulated, and similar quantity distribution occurred in NE12/NE18. However, the difference is in NE12/AE18, in which only 1/3 of the miRNAs were upregulated in NE12 and means more negative genes in AE18.Figure 4
Heatmap showed the expression pattern of differentially expressed miRNAs in chrysanthemum embryos. The bar represents the ratio of expression levels of miRNAs between two samples. Three ratios of NE18/AE18, NE12/NE18, and NE12/AE18 indicated the expression pattern of each miRNA in three samples. Red rectangle means upregulation, and green means downregulation. All information for each miRNA list can be found in TableS2. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.


### 3.5. Characteristics of Target Genes in Chrysanthemum Embryo Development

miRNAs regulate the plant development by mediating the expression of target genes. Apart from those redundancy regulated by several miRNAs, a total of 770 target genes were identified and annotated by chrysanthemum transcriptome dataset. The expression level and annotation were presented in TableS5. Most of them were regulated by one miRNA, and some were negatively regulated by miRNAs, such as the transcription factor MYB11 (Unigene2183) highly expressed in NE12, but the targeted miR858b was the lowest expression level in it. Another transcription factor WRKY48 (Unigene25838) was expressed to be the highest in AE18 and lowest in NE12, showing the negative regulation by miR414 (Figure 5, Tables S2 and S5). Several targets were regulated by two miRNAs, and there was no difference in the level of expression among samples. For example, unigene21756 was the target of miR5215 and miR5373, and its FPKM is near the three samples, without negative regulation by these two miRNAs.Figure 5
Validation of target genes in chrysanthemum embryos by qRT-PCR. The orange histogram represented the qRT-PCR validation, and the value of relative expression was shown on the left ofy-axis. The blue slope lines represented the FPKM of each unigene obtained using RNA-Seq [19], and the value of FPKM was shown on the right of y-axis. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.In order to study the target genes that regulate the embryo development, some differentially expressed target genes (the fold of FPKM > 1.5) were selected according to their Nr/Nt annotation, which may be involve in chrysanthemum embryo development (Table2). They contained some transcription factor, genes related to energy metabolism and protein synthesis, and some uncharacterized protein. The expression of transcription factor was various; WRKY and NAC were highest in AE18, but MYB was lowest. Some of the genes associated with auxin and ATP synthesis were downregulated in AE18 (Table 2 and Figure 5).Table 2
Target genes identified in chrysanthemum embryo.
miRNA
ID of target genes
NE12_FPKM
NE18_FPKM
AE18_FPKM
Nr/Nt annotation
miR159a
Unigene13056_All
22.8776
14.0501
14.2128
50S ribosomal protein
miR160a
CL8743.Contig2_All
5.2658
7.4008
5.3906
Auxin response factor
miR160a
Unigene9514_All
17.7131
8.5644
12.0923
Auxin response factor 10
miR164c
Unigene24220_All
1.0397
6.6612
11.3245
NAC transcription factor
miR166a
Unigene33904_All
2.1379
3.8597
1.4223
Homeobox-leucine zipper protein
miR167c-3p
CL15007.Contig2_All
16.0571
29.3484
35.2898
AMP-activated protein kinase
miR169b
CL3042.Contig2_All
2.9511
8.0941
7.2676
Nuclear transcription factor Y
miR169b
CL3042.Contig3_All
2.5107
7.0631
5.6611
Nuclear transcription factor Y
miR169b
Unigene26767_All
11.8519
30.0804
26.5044
Nuclear transcription factor Y
miR172c
Unigene23005_All
1.827
3.0981
4.1985
Floral homeotic protein APETALA2
miR2653a
Unigene7445_All
0.4387
4.6181
26.8109
Solanum lycopersicum RAP2.6-like protein gene
miR2661
Unigene35846_All
1.8822
4.5798
0.9556
ATP binding protein
miR2868
Unigene46_All
10.1108
17.9793
28.7056
Protein NLP6-like
miR396a-3p
CL936.Contig4_All
10.6134
9.8829
7.7245
MAK16 protein-like protein
miR397a
Unigene30076_All
29.7323
13.9125
13.4372
Glucan endo-1,3-beta-glucosidase 7-like
miR397a
Unigene4252_All
1.9732
5.4082
2.1149
Glycine max laccase-4-like
miR403
CL516.Contig5_All
26.1176
43.3925
45.8309
Argonaute protein group
miR414
CL4550.Contig3_All
12.7195
19.889
8.8036
ADP-glucose pyrophosphorylase small subunit 2
miR414
CL2586.Contig3_All
2.6852
5.6627
7.2585
ATP binding protein
miR414
CL12204.Contig1_All
35.1412
15.8001
13.8476
CAPIP2
miR414
Unigene19498_All
3.0041
2.9203
1.0519
F-box protein At5g07610
miR414
Unigene12916_All
31.2274
55.6863
66.4644
F-box/LRR-repeat protein
miR414
Unigene10707_All
0.158
2.7039
4.4112
Glutamate receptor 2.3
miR414
CL7369.Contig4_All
15.1331
11.8677
23.1584
Glycine max surfeit locus protein 2-like
miR414
CL4274.Contig1_All
0.1141
5.574
1.6227
Histone H2B.2
miR414
Unigene32715_All
2.8017
9.3813
2.9311
Minichromosome maintenance factor
miR414
CL7491.Contig4_All
56.838
23.8057
28.7017
NA
miR414
CL7491.Contig3_All
40.3425
15.2757
19.8302
NA
miR414
CL7491.Contig1_All
40.1769
18.4149
17.7948
Proline-rich protein
miR414
CL11069.Contig1_All
1.594
2.2067
4.0463
Protein NLP6-like
miR414
Unigene30228_All
1.7963
4.9734
9.0057
Uncharacterized protein
miR414
Unigene31458_All
33.3706
106.7705
151.1735
Uncharacterized protein
miR414
Unigene28799_All
0.6622
1.2
2.6895
Uncharacterized protein
miR414
CL15100.Contig1_All
2.7599
11.7908
10.9824
Uncharacterized protein
miR414
CL5368.Contig4_All
23.3918
71.719
84.0485
Uncharacterized protein
miR414
CL5528.Contig1_All
3.1358
4.9787
1.9595
Uncharacterized protein
miR414
Unigene35998_All
2.9444
0.2695
0.8154
Uncharacterized protein
miR414
Unigene19911_All
3.7645
2.6847
0.7964
Uncharacterized protein
miR414
Unigene25838_All
2.1786
8.7075
24.971
Transcription factor WRKY
miR414
Unigene27592_All
109.2811
89.6706
72.8937
ATP synthase subunit O
miR5021
CL13854.Contig1_All
20.2409
11.1804
9.5598
ATP binding protein
miR5021
Unigene49221_All
0.4791
0.8039
3.8915
NA
miR5021
Unigene8421_All
8.1632
6.9848
3.84
PAPS-reductase-like protein
miR5021 miR5293
CL11776.Contig1_All
39.399
34.2386
34.2098
Uncharacterized protein
miR5049c
CL559.Contig1_All
33.8479
33.5268
45.5856
Chrysanthemum x morifolium DREBa
miR5049c miR5052a
Unigene1703_All
1.4503
2.1902
4.4178
NA
miR5077
CL14985.Contig3_All
2.6774
4.7455
5.6638
Global transcription factor group
miR5227
CL4864.Contig1_All
58.2913
17.2209
11.1415
ATP-citrate synthase beta chain protein
miR5270a
Unigene3852_All
2.5276
5.0261
13.2176
Uncharacterized protein
miR5293
CL7596.Contig2_All
999.1282
2027.7828
2112.3254
Dehydrin family protein
miR5293
Unigene16818_All
9.169
7.1012
3.1716
Glycine max reticulon-like protein B21-like
miR5293
Unigene31007_All
34.275
32.3294
106.0944
Uncharacterized protein
miR5298b
Unigene4147_All
12.7734
20.3587
42.7602
Uncharacterized protein
miR5298b miR5380a
CL5088.Contig2_All
22.2504
20.7529
20.0829
Nucleosome assembly protein 1-like
miR5380a
Unigene7407_All
6.0977
88.9232
149.3867
Cytochrome P450
miR5380a
Unigene19991_All
0.6069
3.055
1.0955
Ferric reduction oxidase
miR5215
CL965.Contig5_All
0
0.5661
1.5225
Disease-resistance protein NRSA1
miR5215 miR5373
Unigene21756_All
5.9645
5.6689
6.2018
ATP binding protein
miR5373
CL561.Contig3_All
39.0762
27.122
25.8912
Glycyl-tRNA synthetase 1
miR5640
Unigene11631_All
1.3529
38.198
4.7781
Glutathione S-transferase 1
miR5641
Unigene3971_All
21.8502
14.3017
11.807
Glycosyl transferase
miR5641
Unigene7799_All
1.2821
1.6384
5.4077
RPL5A-related protein
miR6111-5p
Unigene29644_All
5.3413
11.3752
17.9804
Cytokinin-regulated kinase 1
miR6220-3p miR818a
CL1844.Contig1_All
9.9143
3.9489
3.9549
Benzoyl CoA benzoic acid benzoyltransferase
miR845b-5p
Unigene2353_All
10.1133
40.3692
61.9776
NA
miR858b
Unigene2183_All
19.3548
9.2995
5.5141
Transcription factor MYB11


### 3.6. Validation of qRT-PCR

qRT-PCR is an efficient and accurate way to examine the result of RNA-Seq. In the present study, a total of 12 miRNAs were double tested by qRT-PCR, suggesting the high compatibility of the expression pattern between the two methods (Figure6). Seven of them had the highest expression in AE18, such as miR169b, miR159a, and miR858b. miR167c-3p was verified with the highest expression in NE18 by two methods. miR414 that regulated the most target genes with various expression levels was expressed dominantly in NE12 (Figure 6).Figure 6
Validation of miRNAs in chrysanthemum embryos by qRT-PCR. The orange histogram represented the qRT-PCR validation, and the left ofy-axis was the value of relative expression in three samples. The blue slope lines represented the expression level (on the right of y-axis) detected by RNA-Seq. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.


## 3.1. Sequencing Analysis of sRNA

To explore the regulation of miRNAs in chrysanthemum embryo development when paternal chromosomes were doubled, three cDNA libraries of small RNAs were constructed, which were named NE12, NE18, and AE18. All clean reads were obtained by filtering the low-quality sequences, adapter sequences, and poly-A sequences shorter than 18 nt, which altogether resulted in more than 99.79% of raw reads. When these small RNA tags were mapped to genome, about 20% of reads in each sample were matched. The unique sRNAs matched to genome in the three libraries were 3,864,037, 3,780,667, and 3,734,164, accounting for 9.62%, 10.04%, and 10.29% of all unique sRNAs (Table1).Table 1
Quantity and type distribution of small RNAs in chrysanthemum embryo.
NE12
NE18
AE18
Category
Unique sRNAs
Percent (%)
Total sRNAs
Percent (%)
Unique sRNAs
Percent (%)
Total sRNAs
Percent (%)
Unique sRNAs
Percent (%)
Total sRNAs
Percent (%)
Total
3,864,037
100
11,535,216
100
3,780,667
100
10,805,231
100
3,734,164
100
11,220,674
100
Mapping to genome
371,560
9.62
2,251,582
19.52
379,579
10.04
2,194,212
20.31
384,290
10.29
2,383,150
21.24
miRNA
9900
0.26
410,884
3.56
9192
0.24
343,417
3.18
9871
0.26
360,354
3.21
rRNA
57,663
1.49
547,936
4.75
68,080
1.80
662,520
6.13
81,807
2.19
852,436
7.60
snRNA
1423
0.04
3088
0.03
1475
0.04
3034
0.03
1564
0.04
3278
0.03
snoRNA
421
0.01
757
0.01
461
0.01
789
0.01
613
0.02
1205
0.01
tRNA
7221
0.19
72,063
0.62
7332
0.19
97,766
0.90
9506
0.25
118,373
1.05
Unannotated
3,787,409
98.02
10,500,488
91.03
3,694,127
97.71
9,697,705
89.75
3,630,803
97.23
9,885,028
88.10When these unique sRNAs were aligned to the GenBank and miRBase, many types of sRNAs were identified, including miRNA, rRNA, snRNA, snoRNA, and tRNA, but the vast majority of sequences were unannotated. The length of different types of sRNAs was discrepant; since in chrysanthemum embryo, the most abundant classes of sRNA showed the length of 24 nt (dominant siRNA), then 21 nt (mainly miRNAs), and 22 nt (Figure1). In the present study, the unique miRNAs in each library were taken into consideration for subsequent analysis, and the proportion of miRNAs was 0.26%, 0.24%, and 0.26% in libraries of NE12, NE18, and AE18 (Table 1).Figure 1
Length distribution of sRNAs in chrysanthemum embryo.x-axis is the length of sRNA distribution, and y-axis is the proportion of the sRNAs of different lengths. NE12, NE18, and AE18 mean the normal embryo at 12 DAP (days after pollination), normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.


## 3.2. miRNAs and Target Genes Identified in Three Libraries

In all samples, a total of 170 conversed miRNAs were identified in miRBase, and 130, 131, and 132 miRNAs were expressed in NE12, NE18, and AE18, respectively. 100 miRNAs (accounting for 58.8%) were detected in three samples; however, some were in two different samples, such as 10 miRNAs in NE12 and NE18, 7 miRNAs in NE18 and AE18, and 6 miRNAs in NE12 and AE18. Venn diagram (Figure2) presented the quantity distribution of conserved miRNAs in chrysanthemum embryo.Figure 2
Venn diagram presented the quantity distribution of conserved miRNAs in chrysanthemum embryos. There were 130, 131, and 132 miRNAs expressed in NE12, NE18, and AE18, respectively. In the middle, 100 means that the number of miRNAs was detected in three samples. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.Expression level is an important feature to explain the regulation function of miRNAs, and they commonly varied greatly in different samples. Here, the expression of 170 miRNAs with sequences was normalized and analyzed. The most abundant miRNAs in three samples were miR156a, miR157a, and miR166a, with the expression more than a thousand (TableS2). However, some miRNAs highly expressed in a particular sample, such as the expression of miR398b-5p, are 3524 in AE18 but neither in NE12 nor NE18, miR5721 only expressed in NE18, and miR5662 only in NE12.To study the biological function of miRNA in chrysanthemum embryo development, the sequences of target mRNA were paid attention. A total of 770 target genes were identified and regulated by 88 miRNAs (51.76% of all miRNAs), and miRNA414 had the most targets (347 unigenes), followed were miR5293 (61 targets) and miR5021 (44 targets) (TableS3). 34 target genes were regulated by two miRNAs, such as CL1999. Contig2 was targeted by miR165a and miR166a; unigene10412 was the target of miR156a and miR157a (Table S3).


## 3.3. Functional Annotation of Target Genes

Gene Ontology (GO) is a standardized classification system of gene used to describe the characteristics of genes and gene products in organisms. The result showed that 517 target genes were classified into three gene ontology categories: biological processes, cellular components, and molecular functions (Figure3). Of the 40 functional categories, the number of genes in each sample is different, but in most of categories, NE12 had the largest number but least in AE18. In the third ontology of molecular function, some target genes were expressed in a specific sample, such as the genes related to “electron carrier activity” only expressed in NE12, two genes in the category of “enzyme regulator activity” expressed in NE18, but there was no gene regulated “molecular transducer activity” in NE18 (Figure 3).Figure 3
GO classification of target genes in three samples.x-axis means three categories of GO: biological process, cellular component, and molecular function. y-axis is the number of genes in three samples associated with each subcategory. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.To understand the biological function of target genes, KEGG analysis, as with GO, was also used to analyze candidate targets. In NE12, a total of 337 target genes had the biological function on 207 KEGG pathways, but there were less target genes and pathways in NE18 and AE18. 211 target genes with 118 pathway annotation in NE18 and 179 annotated target genes with 113 pathways in AE18 (TableS4). Some target genes were annotated in dozens of pathways that are only in library of NE12, including some energy metabolism pathways, such as ko00280 (valine, leucine, and isoleucine degradation), ko00310 (lysine degradation), and ko00410 (beta-alanine metabolism).


## 3.4. Differentially Expressed miRNAs in Chrysanthemum Embryo

The aim of this study is to analyze the different miRNAs possibly involved in chrysanthemum embryo development, so we identified 112 differentially expressed miRNAs with at least 1.5-fold after standardized expression. The heatmap showed the expression pattern of these miRNAs (Figure4). They were assembled in three groups depending on the expression trend. In group NE18/AE18, the number of upregulated miRNAs was almost the same as downregulated, and similar quantity distribution occurred in NE12/NE18. However, the difference is in NE12/AE18, in which only 1/3 of the miRNAs were upregulated in NE12 and means more negative genes in AE18.Figure 4
Heatmap showed the expression pattern of differentially expressed miRNAs in chrysanthemum embryos. The bar represents the ratio of expression levels of miRNAs between two samples. Three ratios of NE18/AE18, NE12/NE18, and NE12/AE18 indicated the expression pattern of each miRNA in three samples. Red rectangle means upregulation, and green means downregulation. All information for each miRNA list can be found in TableS2. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.


## 3.5. Characteristics of Target Genes in Chrysanthemum Embryo Development

miRNAs regulate the plant development by mediating the expression of target genes. Apart from those redundancy regulated by several miRNAs, a total of 770 target genes were identified and annotated by chrysanthemum transcriptome dataset. The expression level and annotation were presented in TableS5. Most of them were regulated by one miRNA, and some were negatively regulated by miRNAs, such as the transcription factor MYB11 (Unigene2183) highly expressed in NE12, but the targeted miR858b was the lowest expression level in it. Another transcription factor WRKY48 (Unigene25838) was expressed to be the highest in AE18 and lowest in NE12, showing the negative regulation by miR414 (Figure 5, Tables S2 and S5). Several targets were regulated by two miRNAs, and there was no difference in the level of expression among samples. For example, unigene21756 was the target of miR5215 and miR5373, and its FPKM is near the three samples, without negative regulation by these two miRNAs.Figure 5
Validation of target genes in chrysanthemum embryos by qRT-PCR. The orange histogram represented the qRT-PCR validation, and the value of relative expression was shown on the left ofy-axis. The blue slope lines represented the FPKM of each unigene obtained using RNA-Seq [19], and the value of FPKM was shown on the right of y-axis. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.In order to study the target genes that regulate the embryo development, some differentially expressed target genes (the fold of FPKM > 1.5) were selected according to their Nr/Nt annotation, which may be involve in chrysanthemum embryo development (Table2). They contained some transcription factor, genes related to energy metabolism and protein synthesis, and some uncharacterized protein. The expression of transcription factor was various; WRKY and NAC were highest in AE18, but MYB was lowest. Some of the genes associated with auxin and ATP synthesis were downregulated in AE18 (Table 2 and Figure 5).Table 2
Target genes identified in chrysanthemum embryo.
miRNA
ID of target genes
NE12_FPKM
NE18_FPKM
AE18_FPKM
Nr/Nt annotation
miR159a
Unigene13056_All
22.8776
14.0501
14.2128
50S ribosomal protein
miR160a
CL8743.Contig2_All
5.2658
7.4008
5.3906
Auxin response factor
miR160a
Unigene9514_All
17.7131
8.5644
12.0923
Auxin response factor 10
miR164c
Unigene24220_All
1.0397
6.6612
11.3245
NAC transcription factor
miR166a
Unigene33904_All
2.1379
3.8597
1.4223
Homeobox-leucine zipper protein
miR167c-3p
CL15007.Contig2_All
16.0571
29.3484
35.2898
AMP-activated protein kinase
miR169b
CL3042.Contig2_All
2.9511
8.0941
7.2676
Nuclear transcription factor Y
miR169b
CL3042.Contig3_All
2.5107
7.0631
5.6611
Nuclear transcription factor Y
miR169b
Unigene26767_All
11.8519
30.0804
26.5044
Nuclear transcription factor Y
miR172c
Unigene23005_All
1.827
3.0981
4.1985
Floral homeotic protein APETALA2
miR2653a
Unigene7445_All
0.4387
4.6181
26.8109
Solanum lycopersicum RAP2.6-like protein gene
miR2661
Unigene35846_All
1.8822
4.5798
0.9556
ATP binding protein
miR2868
Unigene46_All
10.1108
17.9793
28.7056
Protein NLP6-like
miR396a-3p
CL936.Contig4_All
10.6134
9.8829
7.7245
MAK16 protein-like protein
miR397a
Unigene30076_All
29.7323
13.9125
13.4372
Glucan endo-1,3-beta-glucosidase 7-like
miR397a
Unigene4252_All
1.9732
5.4082
2.1149
Glycine max laccase-4-like
miR403
CL516.Contig5_All
26.1176
43.3925
45.8309
Argonaute protein group
miR414
CL4550.Contig3_All
12.7195
19.889
8.8036
ADP-glucose pyrophosphorylase small subunit 2
miR414
CL2586.Contig3_All
2.6852
5.6627
7.2585
ATP binding protein
miR414
CL12204.Contig1_All
35.1412
15.8001
13.8476
CAPIP2
miR414
Unigene19498_All
3.0041
2.9203
1.0519
F-box protein At5g07610
miR414
Unigene12916_All
31.2274
55.6863
66.4644
F-box/LRR-repeat protein
miR414
Unigene10707_All
0.158
2.7039
4.4112
Glutamate receptor 2.3
miR414
CL7369.Contig4_All
15.1331
11.8677
23.1584
Glycine max surfeit locus protein 2-like
miR414
CL4274.Contig1_All
0.1141
5.574
1.6227
Histone H2B.2
miR414
Unigene32715_All
2.8017
9.3813
2.9311
Minichromosome maintenance factor
miR414
CL7491.Contig4_All
56.838
23.8057
28.7017
NA
miR414
CL7491.Contig3_All
40.3425
15.2757
19.8302
NA
miR414
CL7491.Contig1_All
40.1769
18.4149
17.7948
Proline-rich protein
miR414
CL11069.Contig1_All
1.594
2.2067
4.0463
Protein NLP6-like
miR414
Unigene30228_All
1.7963
4.9734
9.0057
Uncharacterized protein
miR414
Unigene31458_All
33.3706
106.7705
151.1735
Uncharacterized protein
miR414
Unigene28799_All
0.6622
1.2
2.6895
Uncharacterized protein
miR414
CL15100.Contig1_All
2.7599
11.7908
10.9824
Uncharacterized protein
miR414
CL5368.Contig4_All
23.3918
71.719
84.0485
Uncharacterized protein
miR414
CL5528.Contig1_All
3.1358
4.9787
1.9595
Uncharacterized protein
miR414
Unigene35998_All
2.9444
0.2695
0.8154
Uncharacterized protein
miR414
Unigene19911_All
3.7645
2.6847
0.7964
Uncharacterized protein
miR414
Unigene25838_All
2.1786
8.7075
24.971
Transcription factor WRKY
miR414
Unigene27592_All
109.2811
89.6706
72.8937
ATP synthase subunit O
miR5021
CL13854.Contig1_All
20.2409
11.1804
9.5598
ATP binding protein
miR5021
Unigene49221_All
0.4791
0.8039
3.8915
NA
miR5021
Unigene8421_All
8.1632
6.9848
3.84
PAPS-reductase-like protein
miR5021 miR5293
CL11776.Contig1_All
39.399
34.2386
34.2098
Uncharacterized protein
miR5049c
CL559.Contig1_All
33.8479
33.5268
45.5856
Chrysanthemum x morifolium DREBa
miR5049c miR5052a
Unigene1703_All
1.4503
2.1902
4.4178
NA
miR5077
CL14985.Contig3_All
2.6774
4.7455
5.6638
Global transcription factor group
miR5227
CL4864.Contig1_All
58.2913
17.2209
11.1415
ATP-citrate synthase beta chain protein
miR5270a
Unigene3852_All
2.5276
5.0261
13.2176
Uncharacterized protein
miR5293
CL7596.Contig2_All
999.1282
2027.7828
2112.3254
Dehydrin family protein
miR5293
Unigene16818_All
9.169
7.1012
3.1716
Glycine max reticulon-like protein B21-like
miR5293
Unigene31007_All
34.275
32.3294
106.0944
Uncharacterized protein
miR5298b
Unigene4147_All
12.7734
20.3587
42.7602
Uncharacterized protein
miR5298b miR5380a
CL5088.Contig2_All
22.2504
20.7529
20.0829
Nucleosome assembly protein 1-like
miR5380a
Unigene7407_All
6.0977
88.9232
149.3867
Cytochrome P450
miR5380a
Unigene19991_All
0.6069
3.055
1.0955
Ferric reduction oxidase
miR5215
CL965.Contig5_All
0
0.5661
1.5225
Disease-resistance protein NRSA1
miR5215 miR5373
Unigene21756_All
5.9645
5.6689
6.2018
ATP binding protein
miR5373
CL561.Contig3_All
39.0762
27.122
25.8912
Glycyl-tRNA synthetase 1
miR5640
Unigene11631_All
1.3529
38.198
4.7781
Glutathione S-transferase 1
miR5641
Unigene3971_All
21.8502
14.3017
11.807
Glycosyl transferase
miR5641
Unigene7799_All
1.2821
1.6384
5.4077
RPL5A-related protein
miR6111-5p
Unigene29644_All
5.3413
11.3752
17.9804
Cytokinin-regulated kinase 1
miR6220-3p miR818a
CL1844.Contig1_All
9.9143
3.9489
3.9549
Benzoyl CoA benzoic acid benzoyltransferase
miR845b-5p
Unigene2353_All
10.1133
40.3692
61.9776
NA
miR858b
Unigene2183_All
19.3548
9.2995
5.5141
Transcription factor MYB11


## 3.6. Validation of qRT-PCR

qRT-PCR is an efficient and accurate way to examine the result of RNA-Seq. In the present study, a total of 12 miRNAs were double tested by qRT-PCR, suggesting the high compatibility of the expression pattern between the two methods (Figure6). Seven of them had the highest expression in AE18, such as miR169b, miR159a, and miR858b. miR167c-3p was verified with the highest expression in NE18 by two methods. miR414 that regulated the most target genes with various expression levels was expressed dominantly in NE12 (Figure 6).Figure 6
Validation of miRNAs in chrysanthemum embryos by qRT-PCR. The orange histogram represented the qRT-PCR validation, and the left ofy-axis was the value of relative expression in three samples. The blue slope lines represented the expression level (on the right of y-axis) detected by RNA-Seq. NE12, NE18, and AE18 mean the normal embryo at 12 DAP, normal embryo at 18 DAP, and abnormal embryo at 18 DAP, respectively.


## 4. Discussion

### 4.1. Identified miRNAs in Chrysanthemum Embryo

miRNAs can regulate the plant embryonic development, which important and diverse roles have been studied in various species; however, the specific function of individual miRNA is still uncharacterized for stage-specific embryo [5]. Next-generation sequencing makes it easier to identify individual miRNA families from different organs or treated plants. In Pinus taeda [15] and Brassica napus [28], miRNAs were identified from zygotic embryos at late developmental stages. More and more researchers have provided abundant evidences that miRNAs are required for the majority of embryonic cell differentiation and development in Arabidopsis [29, 30]. In chrysanthemum cross breeding, embryonic development is a crucial stage for seed formation, but embryo abortion is prevalent in chrysanthemum distant hybridization [31, 32]. Previous studies explored this barrier at the level of cell structure, gene expression, and miRNA regulation [19, 24], in which 227 miRNAs were identified in hybrid embryos from cross C. morifolium and diploid C. nankingense. Because of the importance of chromosome doubling in cross breeding, in the present study, we performed the hybridization using C. morifolium and tetraploid C. nankingense; 179 miRNAs were identified in three embryo samples (Table S2). Compared with two crosses, less miRNAs were identified after paternal chromosome doubling; as a result, 135 miRNAs were identified simultaneously in two crosses, 44 new miRNAs in the present cross, such as miR169b, miR440, and miR528-5p, but 92 miRNAs expressed in the previous cross were not detected here, such as miR172a, miR172b, and miR391 [24]. These similarities and differences suggest that most of the identified miRNAs in two crosses are the main factor regulating embryonic development in chrysanthemum and those expressed in specific cross may regulate the development depending on whether the male chromosome doubled. To some extent, it implicated that chromosome doubling has an effect on embryo development in distant hybridization by regulation of miRNAs.


### 4.2. The Role of miRNA Chrysanthemum Embryo

In plant embryo, miRNAs can mediate their downstream targets and regulate the expression of some transcription factors or other key developmental regulators [5]. The overexpression of miRNAs and their target genes have allowed assignment of developmental roles in embryonic, vegetative, and floral organ boundary formation [33, 34], such as miR159, miR164, miR165/166, miR172, and miR319 families. It has demonstrated that miR164 directly regulates NAC domain genes for making the function of normal plant morphogenesis and normal embryonic [34]. In the present study, there was no miR164 detected, but miR164c was expressed in three samples and higher expression level in abnormal embryo (Table S2). Previous study showed the highest expression of miR164c in normal embryo at 18 DAP when the male was diploid C. nankingense [24], suggesting that miR164c may have the important regulatory role at the late stage of embryonic development, and the lower expression level of miR164c in normal embryo at 18 DAP may be beneficial to chrysanthemum embryo development normally. Similar situation happened in miR159a; there was no miR159 detected in chrysanthemum embryo, but the miR159a was differentially expressed in normal and abnormal embryo at 18 DAP, suggesting the connected function with miR164a during chrysanthemum embryonic development.miR172 showed the regulatory function for early flowering and floral organ identity defects by regulatingAP2 and TOE [35]. miR172c was detected in this study, and the expression was zero in samples of NE12 and AE18, which means that miR172c mainly mediates gene expression in postembryonic development and these genes may promote embryo maturation. DCL1 was a key gene for embryo development by embryo lethality, which was regulated by miR163 [36, 37]. However, there was no miR163 identified in chrysanthemum embryo whether the cross has paternal chromosome doubling or not. This result exhibits that DCL1 was not regulated by miR163 and participated in embryonic lethality in chrysanthemum.


### 4.3. miRNA-Mediated Target Genes in Chrysanthemum Embryo

Research in the last decade has demonstrated that miRNAs make the crucial roles during plant embryogenesis by regulating various genes and pathways [38]. It contains the process of spatial control of differentiation, regulation of auxin responses, and temporal control of differentiation [5]. A study on miR160-resistant ARF17 transgenes showed the defected cotyledons, suggesting that miR160 negatively regulated genes involved in auxin signaling that is critical for proper development of the embryo and cotyledons [39]. Auxin response is a critical biological pathway for embryonic development, and it has also been reported in chrysanthemum [22]. In the present study, miR160 was identified in NE12 and NE18, highly expressed in NE18 (Table S2). The result of target gene prediction showed that the targets of miR160 were auxin response factor with differential expression among samples (Table 2). In the late embryonic development (18 DAP), the target genes (unigene9514) of auxin response factor downregulated compared with NE12. Besides, there are a portion of target genes involved in energy metabolism according to Nr/Nt annotation, and their expression was downregulated in an abnormal embryo (Table 2). Transcriptome provided the evidence of importance of energy synthesis for normal embryo development [22]; here, these identified target genes were regulated by miRNAs, such as miR414, miR2661, and miR5021, also support the significance of energy metabolism for chrysanthemum embryo development.


## 4.1. Identified miRNAs in Chrysanthemum Embryo

miRNAs can regulate the plant embryonic development, which important and diverse roles have been studied in various species; however, the specific function of individual miRNA is still uncharacterized for stage-specific embryo [5]. Next-generation sequencing makes it easier to identify individual miRNA families from different organs or treated plants. In Pinus taeda [15] and Brassica napus [28], miRNAs were identified from zygotic embryos at late developmental stages. More and more researchers have provided abundant evidences that miRNAs are required for the majority of embryonic cell differentiation and development in Arabidopsis [29, 30]. In chrysanthemum cross breeding, embryonic development is a crucial stage for seed formation, but embryo abortion is prevalent in chrysanthemum distant hybridization [31, 32]. Previous studies explored this barrier at the level of cell structure, gene expression, and miRNA regulation [19, 24], in which 227 miRNAs were identified in hybrid embryos from cross C. morifolium and diploid C. nankingense. Because of the importance of chromosome doubling in cross breeding, in the present study, we performed the hybridization using C. morifolium and tetraploid C. nankingense; 179 miRNAs were identified in three embryo samples (Table S2). Compared with two crosses, less miRNAs were identified after paternal chromosome doubling; as a result, 135 miRNAs were identified simultaneously in two crosses, 44 new miRNAs in the present cross, such as miR169b, miR440, and miR528-5p, but 92 miRNAs expressed in the previous cross were not detected here, such as miR172a, miR172b, and miR391 [24]. These similarities and differences suggest that most of the identified miRNAs in two crosses are the main factor regulating embryonic development in chrysanthemum and those expressed in specific cross may regulate the development depending on whether the male chromosome doubled. To some extent, it implicated that chromosome doubling has an effect on embryo development in distant hybridization by regulation of miRNAs.


## 4.2. The Role of miRNA Chrysanthemum Embryo

In plant embryo, miRNAs can mediate their downstream targets and regulate the expression of some transcription factors or other key developmental regulators [5]. The overexpression of miRNAs and their target genes have allowed assignment of developmental roles in embryonic, vegetative, and floral organ boundary formation [33, 34], such as miR159, miR164, miR165/166, miR172, and miR319 families. It has demonstrated that miR164 directly regulates NAC domain genes for making the function of normal plant morphogenesis and normal embryonic [34]. In the present study, there was no miR164 detected, but miR164c was expressed in three samples and higher expression level in abnormal embryo (Table S2). Previous study showed the highest expression of miR164c in normal embryo at 18 DAP when the male was diploid C. nankingense [24], suggesting that miR164c may have the important regulatory role at the late stage of embryonic development, and the lower expression level of miR164c in normal embryo at 18 DAP may be beneficial to chrysanthemum embryo development normally. Similar situation happened in miR159a; there was no miR159 detected in chrysanthemum embryo, but the miR159a was differentially expressed in normal and abnormal embryo at 18 DAP, suggesting the connected function with miR164a during chrysanthemum embryonic development.miR172 showed the regulatory function for early flowering and floral organ identity defects by regulatingAP2 and TOE [35]. miR172c was detected in this study, and the expression was zero in samples of NE12 and AE18, which means that miR172c mainly mediates gene expression in postembryonic development and these genes may promote embryo maturation. DCL1 was a key gene for embryo development by embryo lethality, which was regulated by miR163 [36, 37]. However, there was no miR163 identified in chrysanthemum embryo whether the cross has paternal chromosome doubling or not. This result exhibits that DCL1 was not regulated by miR163 and participated in embryonic lethality in chrysanthemum.


## 4.3. miRNA-Mediated Target Genes in Chrysanthemum Embryo

Research in the last decade has demonstrated that miRNAs make the crucial roles during plant embryogenesis by regulating various genes and pathways [38]. It contains the process of spatial control of differentiation, regulation of auxin responses, and temporal control of differentiation [5]. A study on miR160-resistant ARF17 transgenes showed the defected cotyledons, suggesting that miR160 negatively regulated genes involved in auxin signaling that is critical for proper development of the embryo and cotyledons [39]. Auxin response is a critical biological pathway for embryonic development, and it has also been reported in chrysanthemum [22]. In the present study, miR160 was identified in NE12 and NE18, highly expressed in NE18 (Table S2). The result of target gene prediction showed that the targets of miR160 were auxin response factor with differential expression among samples (Table 2). In the late embryonic development (18 DAP), the target genes (unigene9514) of auxin response factor downregulated compared with NE12. Besides, there are a portion of target genes involved in energy metabolism according to Nr/Nt annotation, and their expression was downregulated in an abnormal embryo (Table 2). Transcriptome provided the evidence of importance of energy synthesis for normal embryo development [22]; here, these identified target genes were regulated by miRNAs, such as miR414, miR2661, and miR5021, also support the significance of energy metabolism for chrysanthemum embryo development.


## 5. Conclusion

Polyploid breeding will pay more attention from genomic research in the future as rapid advances in the next generation sequencing technology, which makes unprecedented opportunities to explore and understand the regulatory of genomic or transcriptomic changes. As a critical regulatory factor, the function of miRNAs has attracted a lot of attention during plant growth and development. In chrysanthemum distant hybridization, breeders always faced the barriers existed in embryo development. The present study provided some explanation about the embryo abortion and the miRNAs related to embryo development even their target genes. We propose that late embryonic miRNAs, especially miR164a, regulate NAC transcription factor and thereby affect the embryonic development. miR160 mediated the auxin response, and miR414, miR2661, and miR5021 regulate the genes involved in energy metabolism; together, they regulate the embryo development in chrysanthemum hybridization.



---

*Source: 6790478-2017-12-13.xml*