Fatty acid accumulation analysis based on transcriptome sequencing of developmental herbaceous peony ‘Hangshao’ ( Paeonia lactiflora ‘Hangshao’) seeds

Background: Paeonia lactiflora ‘Hangshao’ is widely cultivated in China because its root can be used to produce raw materials for traditional Chinese medicine ‘Radix Paeoniae Alba’. Due to the presence of abundant unsaturated fatty acids in its seed, it also can be regarded as a new oil plant. However, the process of the biosynthesis of unsaturated fatty acid in herbaceous peony ‘Hangshao’ remained largely unknown. Therefore, transcriptome analysis is helpful to better understand the molecular mechanisms. Results: Five main fatty acids, stearic acid, palmitic acid, oleic acid, linoleic acid and α-linolenic acid, were detected, and their absolute contents increased first and then decreased during seed development. A total of 150,156 Unigenes were obtained by transcriptome sequencing, with an average length of 1,030 bp. There were 1,550 Unigenes annotated in the seven functional databases including NR, NT, GO, KOG, KEGG, SwissProt and InterPro. Based on KEGG database, 1,766 Unigenes were annotated in the lipid metabolic pathways, among which 103, 74 and 70 Unigenes are annotated into fatty acid biosynthesis pathway, fatty acid elongation pathway and unsaturated fatty acid synthesis pathway; respectively. A total of 1480 DEGs were detected. Among them, 83 DEGs were enriched in the fatty acid metabolism pathway, including 12 DEGs involved in the fatty acid biosynthesis and 1 DEG involved in fatty acid elongation. Furthermore, qRT-PCR was used to analyze the expression patterns of nine fatty acid biosynthetic related genes including FBCP, BC, FabD, FabF, FATB, KCR, FAD2, FAD3 and FAD7, and it showed that they all highest expressed at 45 DAF. Conclusions: This study provides the first comprehensive genomic resources characterizing herbaceous peony seeds gene expression at the transcriptional level. These data lay the foundation for elucidating the molecular mechanism of the lipid biosynthesis and fatty acid accumulation for herbaceous peony 'Hangshao'. Abstract Background: Paeonia lactiflora ‘Hangshao’ is widely cultivated in China because its root can be used to produce raw materials for traditional Chinese medicine ‘Radix Paeoniae Alba’. Due to the presence of abundant unsaturated fatty acids in its seed, it also can be regarded as a new oil plant. However, the process of the biosynthesis of unsaturated fatty acid in herbaceous peony ‘Hangshao’ remained largely unknown. Therefore, transcriptome analysis is helpful to better understand the molecular mechanisms. Results: Five main fatty acids, stearic acid, palmitic acid, oleic acid, linoleic acid and α-linolenic acid, were detected, and their absolute contents increased first and then decreased during seed development. A total of 150,156 Unigenes were obtained by transcriptome sequencing, with an average length of 1,030 bp. There were 1,550 Unigenes annotated in the seven functional databases including NR, NT, GO, KOG, KEGG, SwissProt and InterPro. Based on KEGG database, 1,766 Unigenes were annotated in the lipid metabolic pathways, among which 103, 74 and 70 Unigenes are annotated into fatty acid biosynthesis pathway, fatty acid elongation pathway and unsaturated fatty acid synthesis pathway; respectively. A total of 1480 DEGs were detected. Among them, 83 DEGs were enriched in the fatty acid metabolism pathway, including 12 DEGs involved in the fatty acid biosynthesis and 1 DEG involved in fatty acid elongation. Furthermore, qRT-PCR was used to analyze the expression patterns of nine fatty acid biosynthetic related genes including FBCP , BC , FabD , FabF , FATB , KCR , FAD2 , FAD3 and FAD7 , and it showed that they all highest expressed at 45 DAF. Conclusions: This study provides the first comprehensive genomic resources characterizing herbaceous peony seeds gene expression at the transcriptional level. These data lay the foundation for elucidating the molecular mechanism of the lipid biosynthesis and fatty acid accumulation for herbaceous peony 'Hangshao'.

LA and ALA are belonged to PUFA which can reduce thrombosis, decrease cardiovascular and resist cancer [13]. According to the location of the last double carbon bond relative to the terminal methyl end of the molecule, PUFA can be divided into omega-3 (n-3) series and omega-6 (n-6) series. ALA, an important omega-3 PUFA, is the necessary fatty acid which cannot be synthesized in the human body and must be obtained from dietary intake and is the precursor of docosahexaenoic acid (DHA) and eicosapentaenoic acid (EPA), which are essential to the growth and development of the brain and retina [14]. As well, LA, an important omega-6 PUFA, is the necessary fatty acid which cannot be produced by human body and is the precursor to the omega-6 series of fatty acids such as γ-Linolenic acid (C18:3 Δ6, 9,12 , GLA) and Arachidonic acid (C18:4 Δ5, 8,11,14 , AA) [15]. Nevertheless, the high intake of omega-6 PUFA with low intake of omega-3 PUFA might promote the pathogenesis of many diseases, including cardiovascular disease, cancer, and inflammatory and autoimmune diseases [16]. Therefore, the intake of the ratio of omega-6/omega-3 was lower than 5:1 suggested by Food and Agriculture Organization (FAO) of the Untied Nations and the World Health Organzation (WHO) [17].
The previous researches showed that the ratio of omega-6/omega-3 is about 1.5:1 in herbaceous peony 'Hangshao' mature seed oil. However, the ratio of omega-6/omega-3 is increased during the seed development because that the relative content of ALA was decrease gradually and the relative content of LA is increased slightly [10,11].
Nevertheless, few research of the fatty acid accumulation at a molecular level for 6 herbaceous peony seed oil had been reported. Fatty acid biosynthesis in plant is a very complicated biological pathway regulated by various enzymes. RNA sequencing (RNA-Seq) has provided convenience that the key genes of fatty acid biosynthesis and unsaturated fatty acid biosynthesis will be rapidly found at the molecular level. De novo transcriptome sequencing had been applied in identifying and profiling the different expression gene of some oil plants such as soybean [18], peanut [19], rape [20], sunflower [21], olive [22], perilla [23], oil palm [24], tree peony [17], tea-oil camellia [25], Eucommia ulmoides [26], sea buckthorn [27], Plukenetia volubilis [28] and so on. However, herbaceous peony research was mainly concentrated on the thermo-tolerant related differently expressed genes [29], codon usage patterns [30] and anthocyanin biosynthetic genes [31] based on transcriptom sequencing. The fatty acid metabolism of herbaceous peony remained unknown. Meanwhile, the genetic control of the unsaturated fatty acid accumulation of herbaceous peony was unexplored. In this study, the absolute content of five main fatty acids for herbaceous peony 'Hangshao' developing seed was measured by gas chromatography-mass spectrometry (GC-MS), and the first herbaceous peony seed transcriptome is generated to uncover genes related to fatty acid biosynthesis using highthroughput Illumina sequencing technology. In total, 150156 unigenes were obtained from nine samples of herbaceous peony 'Hangshao' seeds transcriptome within three developing stages (30,60 and 90 days after flower, DAF) and most of genes in fatty acid biosynthesis were identified. Nine differentially expressed genes were selected for quantitative real-time PCR (qRT-PCR) analyses at five developmental stages (30, 45, 60, 75 and 90 DAF) which provided a large-scale transcriptome annotation of herbaceous peony 'Hangshao' and valuable gene resources related to fatty acid biosynthesis in seeds of herbaceous peony 'Hangshao'.

7
Fatty acid composition and their absolute content change Five main fatty acids, PA, SA, OA, LA and ALA, were detected through GC-MS analysis of the fatty acid components of the seeds from herbaceous peony 'Hangshao' (Fig. 1A).
Through quantitative analysis of the contents of five fatty acids (Fig. 1B), it was found that the content of ALA, LA and PA showed a trend of first rising and then falling during seed development, and the maximum value appeared in the S4 period. While the content of OA and SA increased all the time during seed development, and the maximum value appeared in S5 period. However, the absolute contents of the five fatty acids in the five periods showed that ALA content was the highest and SA content was the lowest, and the contents were successively ALA > LA > OA > PA > SA. The total contents of five fatty acids from S1 to S5 seeds were 44.9 mg·g -1 , 121.70 mg·g -1 , 211.9 mg·g -1 , 261.1 mg·g -1 and 247.5 mg·g -1 , respectively. Meanwhile, the content of ALA was exceeded 100 mg·g -1 in both S4 and S5 stages.
Sequencing reads filtering and de novo assembly Nine cDNA libraries containing three biological repeats were created from three stages of developing herbaceous peony 'Hangshao' seeds (the early, middle and late stages of seed development, that is 30, 60, and 90 days after flower respectively) collected from the same strain and sequenced using Illumina high-throughput sequencing platform in order to explore the molecular mechanism of fatty acid biosynthesis and accumulation for herbaceous peony 'Hangshao'. The sequencing reads which containing low-quality, adaptor-polluted and high content of unknown base (N) reads, should be processed to be removed before downstream analyses. The ratio of the amount of clean reads for nine samples was more than 82% (Additional file 1). Due to without reference genome of herbaceous peony, clean reads had to be assembled after sequencing to get reference sequence for subsequent analysis. After reads filtering, we used Trinity to perform de novo assembly with clean reads [32]. Then, Tgicl were used twice on cluster transcripts to remove abundance and get the final unigenes list for downstream analyses [33]. After assembling all samples together and filtering the abundance, we got 150,156 unigenes, the total length, average length, N50, and GC content of unigenes were 154,674,160 bp, 1,030 bp, 1,825 bp and 39.97% respectively ( Table 1). The size of unigenes' sequences mainly ranged from 300 to 3000 nt, and unigene number gradually decreased without obvious separation as the size of sequences increased. These results mean that a best continuity and high quality of RNA sequencing was conducted (Additional file 2). Of the unigenes, 58,645 (39.1%) were short sequences between 300 nt to 400 nt, whereas only 8927 (6.9%) were longer than 3,000 nt. These results showed that the short sequences advantaged over other sequences of the unigenes of herbaceous peony. All clean reads were uploaded to National Center for Biotechnology Information (NCBI) Sequence Read  Table 2).
The similarity analysis between herbaceous peony 'Hangshao' seed unigenes and NT Nuleic acid database and NR protein database of NCBI were conducted using Blastn [34].
Based on the NR database function annotation result, the ratio of different species on unigene annotataion was calculated and the distributions map was drew (Fig 2A). The unigenes of herbaceous peony 'Hangshao' seed were matched to those of four species such as Vitis vinifera (32.79%), Nelumbo nucifera (5.84%), Theobroma cacao (5.45%) and Jatropha curcas (3.66%).
Based on the 25 functional groups in KOG database, the unigenes were annotated and the unigene distribution was calculated ( Fig 2B). Among them, unigenes mapped to 'General function prediction only' category was the most with the number of 17,454 (25.89%), followed by 'Signal transduction mechanisms' category with the number of 6,039 (8.96%), and 'Posttranslational modification, protein turnover, chaperones' category with the number of 5,185 (7.69%), whereas the unigenes mapped to 'Cell motility' category (88, 0.13%) were the least. The number of unigenes mapped to 'Lipid transport and metabolism' category was 1,868 (2.77%).
Among the 21 sub-categories, 'global and overview maps' was the maximum group with 12,390 unigenes, followed by 'carbohydrate metabolism' (4,905 unigenes) and 'translation' (4,321 unigenes), and the smallest group containing only 17 unigenes was 'drug resisistance: antimicrobial'. Whereas, the number of the unigenes for 'lipid metabolism' was 1,766.

Unigenes related to fatty acid biosynthesis
Based on the analysis of the unigenes in the KEGG pathway and functional annotation, the critical enzymes related to lipid metabolism pathways had been classified and listed in Table 3. According to the previous studies and references, the map of fatty acid biosynthesis pathway was constructed on the basis of these identified enzymes, including de novo synthesis of fatty acid, fatty acid elongation and biosynthesis of unsaturated fatty acid pathway (Fig 3). A total of 247 unigenes related to fatty acid biosynthesis were In addition, a total of 1,480 DEGs were obtained from the intersection analysis of the three groups of 30 vs 60 DAF, 60 vs 90 DAF and 30 vs 90 DAF (Fig 5A). And the total of 1480 DEGs was annotated to 17 biological processes, 14 cellular components and 12 molecular functions in the GO database ( Fig 5B). In biological process, the largest group was 'metabolic process' which was annotated by DEGs, followed by 'cellular process' and 'single-organism process'. And in cellular component, the largest group was 'membrane' which annotated by DEGs, followed by 'cell' and 'cell part'. In addition, in molecular function, the two largest group were 'catalytic activity' and 'binding'. Moreover, the total of 1480 DEGs was annotated to 19 KEGG pathways (Fig 5C). The largest group was 'global and overview maps' which was annotated by DEGS, followed by 'carboydrate metabolism' and 'lipid metabolism'.  (2) and TER (2), and 20 down-regulated DEGs were indentified to fadD(1), KCS (12), KCR(4), HCD (2) and ECR (1), respectively. Moreover, 23 DEGs were indentified in the KEGG pathway of biosynthesis of unsaturated fatty acid, which were all down-regulated DEGs, and they were indentified to SAD(5), FAD2 (9), FAD3(2), FAD7 (5) and FAD8 (2).  (Fig 6). Overall, the expression patterns of 9 DEGs verified by qrt-pcr were highly consistent with the results of transcriptome sequencing, indicating that the results of DEGs analysis were reliable.

Discussion
Lipid, which was produced by fatty acid and glycerinum as substrate, was stored in the seed in the form of triacylglycerol (TAG). The accumulation of fatty acids is closely related to the accumulation of lipids, that is, the increase of fatty acid content will accelerate the accumulation of lipids [36]. The proportion of fatty acids, especially the proportion of unsaturated fatty acids, is one of the evaluation criteria for oil quality [37]. However, the content of LA and ALA constituted different proportions of essential fatty acids have different effects on human physiological functions [38,39]. At present, many studies have shown that the intake ratio of LA/ALA has the most significant physiological function at 1:1 [40]. In this study, five main fatty acids were detected by GC-MS in the seed of herbaceous peony 'Hangshao', that is PA, SA, OA, LA and ALA, and which was consistent with the study on the component of other peony [11,17]. The absolute contents of the five fatty acid components during the seed development were all the highest for ALA and the lowest for SA. The specific manifestations were ALA>LA>OA>PA>SA. The ratio of LA/ALA was lower during the developmental seeds of herbaceous peony 'Hangshao', maintained between 0.55 and 0.71, close to the optimal ratio of 1:1.
At present, the studies on fatty acid of herbaceous peony are limited to the extraction of seed oil and its component, the related gene cloning and expression of ALA synthesis [11], and few research of the fatty acid biosynthesis for herbaceous peony seed had been reported. With the development of sequencing technology, species without reference genome can use RNA-seq to directly sequence the total mRNA of samples. Transcriptome sequencing is also sensitive to changes of gene expression, which has important reference significance for species with complex genomic backgrounds such as herbaceous peony [17]. In this study, Illumina HiSeq TM 2000 sequencing platform was used to carry out transcriptional sequencing on the seeds of herbaceous peony 'Hanghao' at different development stages, and 100.07 Gb data were obtained. After assembling and removing the redundancy, 150,156 unigenes were obtained, and 67,096 unigenes were annotated in the NR database. In order to obtain comprehensive transcriptome information of herbaceous peony 'Hangshao' seeds, the following measures were adopted in this study.
Firstly, the seeds were collected from the same plants at early (30 DAF), middle (60 DAF) and mature (90 DAF) stages, and cDNA libraries were constructed from the total RNA of the seeds in 3 follicles at each stages, and transcriptome was sequenced in parallel.
Secondly, for unigenes annotation, seven databases (NR, NT, SwissProt, KOG, KEGG, GO and InterPro) were used to get the most comprehensive annotation information possible.
According to the annotation information of unigene from transcriptome database, two important pathway related to fatty acid biosynthesis were found in the developmental seed of herbaceous peony 'Hangshao'. One pathway of fatty acid biosynthesis is in the plasmid, and the other is in the endoplasmic reticulum. In terms of fatty acid biosynthesis in the plasmids, the genes catalyzing the generation of malonyl-CoA including 11 BC,  [41], and catalyze LA to form ALA by Δ15 fatty acid dehydrogenase including FAD3, FAD7 and FAD8 [42]. In this process, FAD genes play a key role in the formation of ALA [43]. After the synthesis of fatty acids, the TAG was assembled in the endoplasmic reticulum through Kennedy pathway, and finally formed into oil body by budding for storage [44].
At present, 17 lipid metabolism pathway and a series of related protein have been annotated in the KEGG database, including 33, 23 and 29 related genes in the pathway for de novo synthesis of fatty acid, fatty acid elongation and biosynthesis of unsaturated fatty acid, respectively [45]. In the de novo synthesis of fatty acids, acetyl CoA is condensed with CO 2 to form malonyl CoA, and then malonyl CoA is catalyzed by FabD to form malonyl-ACP. Then, malonyl-ACP is reacted with acetyl CoA to form butyryl-ACP (C4:0-ACP) through condensation, reduction, dehydration, reduction by a series of enzymes. In this study, the unigenes of 1 BC (down), 1 α-CT (down), 2 BCCP (down), 1 FabD (down), 1 FabG FAD8 (down) were found in the stage of 60 d vs 90 d. The above DEGs expression showed a general downward trend, which may be due to the fact that the peak of lipid accumulation was 75DAF, while the maximum rate of unsaturated fatty acid accumulation was 45 DAF (Fig 1). Although the accumulation amount of fatty acid increased slowly, the accumulation rate of fatty acid decreased gradually. According to transcriptome, the expression of desaturase (SAD, FAD2, FAD3, FAD6, FAD7 and FAD8) showed no significant change at the stage of 30 d vs 60 d, possibly due to the maximum accumulation was happened at 45 DAF (Fig 6). However, with the decrease of fatty acid accumulation rate, the expression level of FAD gene family was gradually reduced from 60 d to 90 d. In the genes of BC, α-CT, BCCP, FabD, FabH and FabG that catalyzed the synthesis of malonyl-CoA in the plastid, the expression of FabG (CL2376.Contig2_All) showed a gradually increasing trend during seed development, however, the expression of other genes showed a decreasing trend, which was consistent with the gradually decreasing rate of fatty acid synthesis.
In the transcriptome of herbaceous peony 'Hangshao' seed, the unigene annotated to FabD, FabH and FATA was only 1, respectively. FabD is malonyl CoA-acyl carrier protein transacylase (MACT), which gene encodes malonyl CoA-ACP transallase that binds to malonyl-CoA to form a stable malonyl-MACT intermediate. Then, the acyl protein binds to the malonyl-MACT intermediate, and the malonyl group is transferred to the mercaptan at the end of the acyl carrier protein by hydrophobic interaction, so as to synthesize malonyl-ACP [46]. Malonyl ACP is a carbon chain donor that extends the carbon chain in the de novo synthesis of fatty acids. Under the catalysis of the complex enzyme system of FAS, it eventually forms C16:0-ACP and C18:0-ACP through multiple cycles. So, FabD is a key enzyme in the fatty acid synthesis pathway, and its expression level is closely related to the lipid accumulation. At present, the overexpression of FadD gene in rape and other plants has a good effect on promoting the accumulation of C16 and C18 fatty acid [47].
FabD belongs to DEGs in the transcriptome of herbaceous peony 'Hangshao' seed, and shows a downward trend in the seed development, which is consistent with the change trend of fatty acid accumulation rate. FabH is ketoacyl-ACP synthase III (KASIII), a key enzyme that initiates de novo synthesis of fatty acid. It acts on the condensation in the first cycle of fatty acid synthesis, which catalyzes condensation of malonyl ACP with acetyl CoA to produce β-ketobutyryl-ACP. Then after reduction, dehydration and reduction, βketobutyryl-ACP was catalyzed to butyryl-ACP (C4:0-ACP) [48]. At present, FabH gene has been cloned from plants such as tung tree [49], sunflower [50] and jatropha [51].

Condensation reaction by FabH was conducted that acetyl CoA receptor substrate release
CoA and bind C2 unit offered by malonyl-ACP after decarboxylation to form β-ketobutyryl-ACP. In addition, some studies showed that FabH enzymes can use other receptor substrate to generate different final products of fatty acid synthesis, which indicate that FabH enzyme is one of the main factors determining the carbon chain structure of fatty acid [52,53]. FabH, as the initator of de novo synthesis of fatty acid, plays a second role in detemining the length of the carbon chain of fatty acid [54]. Moreover, some studies believe that the synthesis rate of fatty acid is affected by FabH enzymatic activity, indicating that it may be one of the rate-limiting enzymes in the fatty acid synthesis pathway [55]. The FPKM value of FabH gene in the transcriptome of herbaceous peony 'Hangshao' seeds decreased with the development of seeds, but had no significant difference, which was consistent with the accumulation rate of fatty acid in herbaceous peony 'Hangshao' seeds. FATA is acyl-ACP thioesterase A, a type of acyl-ACP thioesterase.
De novo synthesis of fatty acid was terminated by acyl-ACP thioesterase (FAT), which hydrolyzes the thiol-lipid bond located between acyl and ACP to form fatty acid. Then, fatty acid are transported from the plastids and then are esterified to form acly-CoA in the plasmid membrane [56]. The specificity of FAT largely determines the chain length and unsaturated degree of vegetable fatty acid in glycerolipid and TAG [57]. FAT is a soluble enzyme encoded by nuclear genes, which can be classified into FATA and FATB according Quantitative analysis of main fatty acid composition Sample hydrolysis: 0.5g seeds were ground into granules and put into a flask. Then, 100 μL methyl nonadecanoate, 100 mg pyrogallic acid, several zeolites and 2 mL 95% ethanol were added and mixed. Subsequently, 10 mL hydrochloric acid solution were added and mixed well. Next, the flask was put into 70℃ water bath to hydrolyze for 40 min. After hydrolysis, the flask was cooled to room temperature.
Fat extraction: the hydrolyzed sample was added with 10 mL 95% ethanol and mixed. And the hydrolysate in the flask was transferred to the separation funnel. Then, the flask and the stopper were washed with 50 mL ethyl ether petroleum ether mixture. Next, the rinse solution was incorporated into the separation funnel, capped, shaken for 5 min, and set aside for 10 min. Subsequently, the extracting solution was collected into 250 mL flask.
The hydrolysate was extracted 3 times as described in the previous steps. Lastly, separatory funnel was flushed by petroleum ether with ether mixture (volume ratio 1:1) and the flushing fluid were collected into a constant weight flash which would be boiled in water bath and dried in 100℃ for 2 h. Fatty acid esterification: 2 mL of 2% sodium hydroxide methanol solution was added into the fat extract, and heated in 80℃ for 30 min. Then, 3 mL 14% boron trifluoride methanol solution was added and heated in 80℃ for 30 min. Next, the mixture was cooled to room temperature and put into a centrifugal tube and 1 mL n-hexane was added. And it was extracted by shock for 2 min, then left for 1 h, and then stratified. 100 μL supernatant was volumed to 1 mL with n-hexane and filtered by 0.45 μm organic phase filter to 2 mL chromatographic sample bottle.
Chromatographic and mass spectrometric conditions: the conditions was referred to the reference [10].
Calculation formula of fatty acid absolute content: W (mg•g -1 ) =C×V×N×K/ (1000×M), where W is the absolute content of fatty acid, C is the concentration of fatty acid methyl Kit With gDNA Eraser (TaKaRa, Japan). Actin (JN105299) was used as an internal control in P. lactiflora [61].All gene-specific primers in this study were shown in Additional file 8, and synthesized by Shanghai Sangon Biological Engineering Technology and Services Co., Ltd.
(Shanghai, China). qRT-PCR was performed using SYBR ® Premix Ex TaqTM (Perfect Real Time) (TaKaRa, Japan). The amplification was carried out under the following conditions: 55 ℃ for 2 min, followed by an initial denaturation step at 95℃ for 30 s, and 40 cycles at 95 ℃ for 5 s, 55 ℃ for 15 s, and 72℃ for 30 s. Relative expression levels of target genes were calculated by the 2 -△△Ct comparative threshold cycle method [62].

Declarations
Ethics approval and consent to participate Not applicable.

Consent for publication
Not applicable.

Availability of data and materials
All data generated or analyzed during this study are included in this published article and its supplementary information files.

Competing interests
The authors declare that they have no competing interests.   Table 1 Quality metrics of transcripts of seeds of Paeonia lactiflora 'Hangshao' Table 2 Functional annotation of uingenes for seeds of Paeonia lactiflora 'Hangshao' in public functional database Table 3 Number of Unigenes related to fatty acid biosynthesis in seeds of Paeonia lactiflora