Ho Wen, Wenqing Wng, Xi Jing, Minyu Wu, Hongjin Bi ,Cuiyun Wu,*, Lirong Shen,*

a Department of Food Science and Nutrition, Zhejiang Key Laboratory for Agro-Food Processing, Zhejiang University, Hangzhou 310058, China

b College of Horticulture Forestry, Tarim University of Land-Reclamation, Aral 843300, China

Keywords:

Korla fragrant pear

Development period

Chlorogenic acid content

De novo transcriptome assembly

Phenylpropanoid biosynthesis pathway

A B S T R A C T

Korla fragrant pear (KFP) with special fragrance is a unique cultivar in Xinjiang, China. In order to explore the biosynthesis molecular mechanism of chlorogenic acid (CGA) in KFP, the samples at different development periods were collected for transcriptome analysis. High performance liquid chromatography analysis showed that CGA contents of KFP at 88, 118 and 163 days after full bloom were (20.96 ± 1.84), (12.01 ± 0.91) and(7.16 ± 0.41) mg/100 g, respectively, and decreased with the fruit development. Pears from these typical 3 periods were selected for de novo transcriptome assemble and 68 059 unigenes were assembled from 444 037 960 clean reads. One ‘phenylpropanoid biosynthesis’ pathway including 57 unigenes, 11 PALs,1 PTAL, 6 4CLs, 9 C4Hs, 25 HCTs and 5 C3’Hs related to CGA biosynthesis was determined. It was found that the expression levels of 11 differentially expressed genes including 1 PAL, 2 C4Hs, 3 4CLs and 5 HCTs were consistent with the change of CGA content. Quantitative polymerase chain reaction analysis further showed that 8 unigenes involved in CGA biosynthesis were consistent with the RNA-seq data. These findings will provide a comprehensive understanding and valuable information on the genetic engineering and molecular breeding in KFP.

1. Introduction

Xinjiang is the main pear production area in China due to its huge temperature difference between day and night with long sunshine time which is very beneficial to fruit growth and nutrient accumulation.Korla fragrant pear (KFP) belongs toPyrus sinkiangensisYü, a unique species grown at southern area of Xinjiang. Because of its crisp taste and rich fragrance, KFP has become an important pear cultivar in Xinjiang.

During growth and development of pear fruits, the physical and chemical properties, titratable acid content, soluble solids content varied significantly with the increase of pear volume and color change [1]. Meanwhile, the content of secondary metabolites,such as phytohormones including abscisic acid and polyphenolic antioxidants in fruits change sharply [2]. Chlorogenic acid (CGA),a main secondary metabolite in some plants such asChromolaena odorataleaves [3],Lonicera macranthoides[4]and potatos [5]has important biological activities such as antioxidant, anti-free radical and anti-inflammatory [6]. Irondi et al. [7]confirmed that the fruit extract ofTetrapleura tetrapteracontaining CGA could effectively inhibit xanthine oxidase and Fe2+-induced lipid peroxidation in tissues of rats. Therefore, a growing number of researchers believe that CGA could prevent various diseases related to oxidative stress, and CCA complex could induce cell death in human colon cancer HCT-116 cells [8]. CGA was also detected in pear fruits. Arbutin and CGA were the main phenolic compounds in peel of 16 pear cultivars [9].It was found that the contents of CGA and arbutin in core was higher than those inflesh which was peel [10]. Previous studies have shown that CGA is synthesized in plants via 3 approaches. The first one is produced by caffeoyl-CoA under the catalysis of hydroxycinnamoyl transferase (HCT). The second one is produced byp-coumaroyl quinic acid which is catalyzed byp-coumarate 3-hydroxylase(C3’H). The third one is produced by caffeoylD-glucose with the catalysis of hydroxycinnamoylD-glucose: quinate hydroxycinnamoyl transferase (HCGQT). The remaining key enzymes involved in CGA biosynthesis are phenylalanine ammonia-lyase (PAL), phenylalanine/tyrosine ammonia-lyase (PTAL), cinnamic acid 4-hydroxylase (C4H)and 4-coumarate-CoA ligase (4CL). The overexpression related genes encoding HCT led to the amassment of CGA indicating that HCT was a rate-limiting enzyme in CGA biosynthesis [11]. Therefore, the first and the second approaches are mainly responsible for the synthesis of CGA. The third approach is only found in a few plants such as tomato leaves [12]. In order to explore the biosynthesis molecular mechanism of CGA in KFP, it is necessary to analyze the gene expression profiles and the key pathways for CGA biosynthesis.

Transcriptome analysis is an effective method to study the relationship between plant growth and gene expression [13,14].Currently, the transcriptome researches such as the regulation of epidermal wax [15], lignification of stone cells [16,17]and anthocyanin biosynthesis [18]of pear fruit have received extensive attention. Studies on polyphenols including changes in CGA content and related gene expressions in pears have been implemented [19,20].In 2013, the draft genome of the pear (P. bretschneideri) was reported for the first time with a 512.0 Mb sequence and 42 812 protein-coding genes [21]. Although some studies have used transcriptomics to study the inheritance and characteristics of KFP in Xinjiang, there were limited studies focusing on the development and bioactive substances of KFP. Pei et al. [22]found that the calyx persistence in KFP was regulated by many genes related to cell wall degradation and plant hormones. The expression ofkfpMYBseemed to be correlated with calyx persistence in KFP [23]. In this study, fruits from KFP during the development were selected forde novotranscriptome assemble to explore the key enzyme genes in the biosynthesis of CGA of KFP,which will provide a theoretical basis for pear species selection and breeding in Xinjiang.

2. Materials and methods

2.1 Plant materials

KFP pear fruit samples were collected on June 7 (58 days after full bloom (DAFB)), June 22 (73 DAFB), July 7 (88 DAFB), July 22(103 DAFB), Aug. 6 (118 DAFB), Aug. 21 (133 DAFB), Sept. 5 (148 DAFB), Sept. 20 (163 DAFB) at Luntai plant Germplasm Resources Garden (41°N, 84°E), Xinjiang Academy of Agricultural Sciences in 2019, respectively. For each sampling time, 5 fruits each tree in different locations (i.e., north, south, center, west, and east) for a total of 25 fruits from 5 trees were collected. Three pears from each period were immediately cut into the mixtures containing both the peel and the flesh,then the samples were stored at –80 °C after frozen in liquid nitrogen for transcriptome analysis. The rest of fruits were sliced and dried to a constant weight at 70 °C for using in CGA content determination.

2.2 Determination of CGA and arbutin content

KFP fruit samples containing both the peel and the flesh were subjected to the following operations, drying for 24 h, cooling for 12 h, passing through a 100 mesh after crushing, weighing 1.000 ×gpowder accurately, diluting with 25 mL methanol(50%,V/V) solution and the extracting in ultrasonic water bath for 30 min. Subsequently, the samples were centrifuged for 30 min(4 °C, 6 000 ×g). Finally, 1 mL supernatant was filtered by 0.22 μm microporous membrane for CGA analysis.

High performance liquid chromatography (HPLC) analysis was performed on a Shimadzu LC-20AT HPLC analyzer (Shimadzu Corporation, Kyoto, Japan) which was equipped with a LC-20AT pump, SIL-20AC autosampler, RF-10AXL fluorescence detector and CBM-20A system controller. Separation was carried out using an Agilent ZORBAX 300SB-C18column (250 mm × 4.6 mm,5 μm). The mobile phase was methanol - 0.4% (V/V) phosphoric acid in elution, the detection wavelength was 280 nm, and the flow rate was 1 mL/min. A volume of 20 μL was injected and the HPLC chromatogram was monitored for 25 min. The content of CGA and arbutin was calculated according to the peak area of 3 parallel experiments.

2.3 RNA extraction and sequencing

Samples collected at 88 DAFB (T1), 118 DAFB (T2) and 163 DAFB (T3) were chosen for RNA-seq according to the development periods and developmental change of CGA content. Total RNA was extracted using Trizol Reagent according the manufacturer’s instructions (Invitrogen) and the genomic DNA was removed by DNase I (TaKara). RNA quality was determined by Agilent 2100 Bioanalyzer and quantified using the ND-2000 (Nanodrop Technologies). The RNA integrity number (RIN) was examined using Agilent 2100 Nanodrop. The transcriptome library was prepared following TruSeq RNA Sample Preparation Kit from Illumina (San Diego, CA, US) using 5 μg of total RNA, which was then sequenced with the Illumina NovaSeq 6000 System supplied by Shanghai Majorbio Bio-pharm Technology Co., Ltd.

2.4 De novo assembly and function annotation

To get clean reads, SeqPrep and Sickle was used to remove the low-quality reads with connectors from the original sequencing data.Trinity (k-mer = 25) was utilized to assemble all clean reads to obtain ade novotranscriptome assembly. The longest transcript of each gene was taken as unigene for further gene function analysis.

Unigenes were compared by Diamond, HMMER3, Blast2Go and KOBAS with six public databases including Non-Redundant Protein Sequence Database (NR), an annotated protein sequence database(SWISS-PROT), evolutionary genealogy of genes: Non-supervised Orthologous Groups (eggNOG), a database of protein families (Pfam),Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes(KEGG).

2.5 Analysis of differentially expressed genes

Expectation-Maximization (RSEM) was used to compare the RNA-Seq data after quality control with thede novotranscriptome assembly sequences to estimate the gene expression levels and obtain the TPM (Transcripts Per Million reads) value of each unigene. The differentially expressed genes (DEGs) between different development periods were obtained with certain screening criteria (P-adjust < 0.05,|log2Fold Change| ≥ 1) using DESeq2 based on negative binomial distribution. TheP-value obtained from the statistical test was corrected by BH (FDR correction with Benjamini/Hochberg) to obtainP-adjust.

2.6 Validation of RNA-seq by quantitative PCR

In order to verify the accuracy of RNA-Seq analysis data,four DEGs and four other unigenes involved in chlorogenic acid biosynthesis were selected and subjected to quantitative PCR (qPCR)detection. Total RNA were reverse-transcribed to synthesize cDNAs using MonScript RT All-in-One Mix (Monad, Zhuhai, China). Primer of selected unigenes was designed using Primer Premier 5 software(Supplementary Table 1). qPCR was performed using a MonAmp?SYBR?Green qPCR Mix (Monad, Suzhou) on ViiA 7 (Applied Biosystems) with Low ROX (Monad). All qPCR experiments were performed with three technical and 3 biological repeats.Tubulin alpha-3was taken as the reference gene and the relative expression level was calculated using the 2-ΔΔCtmethod.

Table 1De novo transcriptome assembly information.

2.7 Statistical analysis

Data analysis was performed using the software, Statistical Product and Service Solutions (SPSS) 22.0. The data from biochemical assays were presented as mean ± SD (standard deviation). The statistical differences between the groups were shown by One-way analysis of variance (ANOVA) test. Duncan’s multiple comparison test was used to compare all the groups in pairs, andP-value ＜ 0.05 was considered to have a significant difference.

3. Results

3.1 Determination of CGA content at different development periods

With the development of fruits, the volume and weight of KFP increased continuously, and the peel color changed from green to yellow and red (Fig. 1B). According to the single KFP fruit weight, it could be inferred that 58-88 DAFB represented in young fruit period,103-133 DAFB represented fruit enlarging period, and 148-163 DAFB represented fruit mature period with the weight of a mature fruit being 149.25 ± 3.30 (Fig. 1A). HPLC analysis results of KFP at 88, 103, 118, 133, 148 and 163 DAFB (Fig. 2) showed that the contents of CGA and arbutin decreased sharply with the fruit growth and development. The CGA contents at the three typical periods at 88, 118 and 163 DAFB were (20.96 ± 1.84), (12.01 ± 0.91) and(7.16 ± 0.41) mg/100 g (dry weight, DW), respectively, which were different significantly (P-value ＜ 0.05) (Fig. 2A). Thus these samples at the three typical periods were selected for transcriptome analysis.

Fig. 1 KFP development periods. (A) The appearance of KFP. (B) The single fruit weight of KFP. Error bars represent the standard error of ten independent experiments. June 7 (58 DAFB), June 22 (73 DAFB), July 7 (88 DAFB), July 22 (103 DAFB), Aug. 6 (118 DAFB), Aug. 21 (133 DAFB), Sept. 5 (148 DAFB),Sept. 20 (163 DAFB).

Fig. 2 (A) The CGA contents of KFP during development (88, 103, 118, 133, 148 and 163 DAFB). (B) The arbutin contents of KFP during development (88,103, 118, 133, 148 and 163 DAFB). Error bars represent the standard error of three independent experiments. Different lowercase letters (a-d) indicate a significant difference from each other at P < 0.05 by Duncan’s multiple comparison test.

3.2 Sequencing and de novo assembly

A total of 444 037 960 clean reads from 448 107 462 raw reads with an average Q30 of 94.21% and an average GC content of 47.32% were obtained after filtration. 68 059 unigenes with an average length of 861 bp and N50 length of 1 519 bp were assembled by Trinity to obtain ade novotranscriptome assembly of KFP. The clean reads were mapped to the assembled unigenes and the mean alignment rate reached 83.50% (Table 1). Nearly 25.25% of unigenes were over 1 000 bp in length, indicating the result ofde novoassembly was effective (Supplementary Fig. 1).

3.3 Function annotation and classification

The functional annotation results of unigenes were obtained by comparing with six databases. Among them, unigenes annotated to NR database were the most, accounting for 56.10% (Table 2).

Table 2Functional annotation of unigenes.

GO database is the world’s largest source of information on gene function and classification, consisting of three categories: molecular function, biological process and cellular component. A total of 25 691 unigenes were assigned into GO database, of which there were 30 204 annotations for ‘cell components’, followed by 29 925 for ‘molecular function’. The least one was ‘biosynthesis process’ with 25 039 annotations. Among all level 2 terms, unigenes annotated to ‘binding’represented the majority (51.99%), followed by ‘catalytic activity’(49.39%) (Fig. 3). It was shown that the same unigenes could be annotated into multiple level terms. Therefore, the total annotation numbers in the same category were greater than that of unigene,which was the same as KEGG annotations. After GO annotation, the functional classification of all unigenes can be roughly understood.

Fig. 3 GO classification of unigenes in KFP.

As an extension of the COG database of NCBI, eggNOG database is a hierarchical, functionally and phylogenetically annotated orthology resource, providing direct homologous grouping of proteins at different classification levels. A total of 24 308 unigenes were annotated and distributed in 23 COG categories. Although the function of 54.13% of unigenes was unknown, the groups‘posttranslational modification, protein turnover, chaperone’ and‘transcription’ were also annotated to 7.88% and 6.32% of unigenes(Supplementary Fig. 2), indicating that the metabolic and biosynthesis processes of pears were active.

KEGG is a great knowledge base for systematically analyzing gene function, associating genomic information and functional information. By comparing with the KEGG database, a total of 15 780 unigenes were annotated and 8 829 of them were assigned to the 144 pathways. Among all first pathways categories, unigenes annotated to ‘metabolism’ were the most, accounting for 63.84%. There were 22 pathways related to secondary metabolism, among which the most unigenes were annotated to ‘phenylpropanoid biosynthesis’ pathway(Table 3). According to further mining of annotation information, we found that there were 57 unigenes encoding 6 enzymes involved in CGA biosynthesis pathway in KFP.

Table 3The top ten KEGG pathways related to secondary metabolism.

3.4 Analysis of differential expression genes

The comparative transcriptome analysis of KPF at the three typical development periods was performed in order to discover DEGs related to CGA biosynthesis. A total of 4 728, 7 712, 2 174 DEGs were found in the ‘T2vsT1’, ‘T3vsT1’ and ‘T3vsT2’ group as the thresholds ofP-adjust < 0.05, FC > 2 and FC < 0.5 (Fig. 4).Among unigenes in ‘T3vsT1’, 3 910 unigenes were up-regulated and 3 802 unigenes were down-regulated, indicating that the gene expression patterns of KFP varied greatly during development, but tended to stabilize in fruit mature period.

Fig. 4 DEGs analysis of KFP. (A) Volcano plots of the DEGs between KFP at T1 and T3. (B) Venn diagram of the groups between different development periods. The cross regions of the circle represent unigene numbers shared by each group. ‘FC’ means fold change. T1, T2, T3 represent 88, 118 and 163 DAFB, respectively.

GO annotation and enrichment analysis were conducted on the DEGs in ‘T3vsT1’ group. The analysis result showed that a total of 47 GO terms (Level 2) including 20 for ‘biological process’ (BP),14 for ‘molecular function’ (MF), and 13 for ‘cellular component’(CC) were annotated. Further enrichment analysis showed the top 20 GO terms in -lg(P-adjust) (Fig. 5). The ‘P-adjust’ represents whether the enrichment results were statistically significant, The bigger the -lg(P-adjust), the more significantly enriched the term was. And the GO terms belonging to MF were the most, with the number reaching 17. Among these 20 terms, ‘integral component of membrane’ was annotated with 1 649 unigenes, followed by ‘oxidoreductase activity’ and ‘transporter activity’ with the unigenes number of 482 and 316.

Fig. 5 Top 20 terms of DEGs GO annotation by the enrichment levels. ‘Rich Factor’ represents the ratio of DEGs number enriched in the GO term to totally annotated unigenes number in the GO term. The larger the rich factor, the greater the enrichment degree was.

KEGG annotation and enrichment analysis were conducted on the DEGs in ‘T3vsT1’ group. The analysis results showed that a total of 137 pathways were annotated, of which the -lg(P-adjust) of ‘plant hormone signal transduction’ and ‘phenylpropanoid biosynthesis’pathway was significantly higher than other pathways, indicating that the DEGs in ‘T3vsT1’ group were significantly enriched in the two pathways. The ‘phenylpropanoid biosynthesis’ pathway was annotated to 76 DEGs, which was the secondly most only to the‘plant hormone signal transduction’ pathway related to the growth and development of KFP. The other pathways with higher enrichment significance were ‘photosynthesis - antenna proteins’, ‘cutin, suberine and wax biosynthesis’, ‘fatty acid biosynthesis’ and so on (Fig. 6).

Fig. 6 Enrichment significance and annotation number of DEGs in the pathway based on the KEGG.

3.5 Analysis of key genes involved in CGA biosynthesis

In the ‘phenylpropanoid biosynthesis’ pathway, we found the key enzymes involved in CGA biosynthesis route 1 and route 2(Fig. 7). A total of 57 unigenes encoding 6 enzymes including 11PALs, 1PTAL, 64CLs, 9C4Hs, 25HCTsand 5C3’Hswere found in the assembled unigenes database (Table 4). However, the UGCT enzyme and HCGQT enzyme in route 3 were not annotated.

Fig. 7 Schematic representation of the CGA biosynthesis pathway coupled with heatmaps of relevant genes participating in the pathway. (Red indicates that the unigene has a higher expression level in the sample, and blue indicates a lower expression level).

Table 4Genes encoding key enzymes involved in KFP CGA biosynthesis.

Previous analysis showed that the CGA content of KFP was significantly different between T1 and T3, therefore the 76 DEGs annotated in ‘phenylpropanoid biosynthesis’ pathway were excavated.A total of 20 DEGs including 1PAL, 3C4H, 54CLsand 11HCTsinvolved in CGA biosynthesis were found. Among the 20 DEGs, 11 down-regulated and 9 up-regulated unigenes were observed. Since the CGA content decreased significantly during KFP development, these 11 down-regulated unigenes including 1PAL, 2C4Hs, 34CLsand 5HCTswere worthy of attention (Table 5). Especially, the expression levels ofDN15799_c0_glencoding 4CL,DN8181_c0_glandDN28670_c0_g3encoding HCT at T3 period were 0.066, 0.059 and 0.088 times of that at T1 period. These results are in agreement with the trend of CGA content at the T3 and T1 of KFP. Therefore, we believed that these 11 unigenes, particularly the unigenes encoding HCT played a key role in the synthesis of KFP CGA.

Table 5Down-regulated DEGs involved in CGA biosynthesis.

3.6 Validation of RNA-Seq data by qPCR

Eight unigenes encoding key enzymes were randomly selected for qPCR to verify the accuracy of RNA-seq. The relative expression levels of these genes were different in different development periods, and the expression patterns were consistent with the data of RNA-seq (Fig. 8). Therefore, our RNA-seq results provided reliable data for future research on key enzyme genes of CGA biosynthesis in KFP.

Fig. 8 Relative expression levels of eight DEGs in three development periods of KFP. We compared expression levels determined by RNA-seg and qPCR analyses for each DEG.

Fig. 8 (Continued)

4. Discussion

KFP with persistent calyx, spindle shape and red green peel is differing from popular pear species,P. pyrifoliaNakai andP. bretschneideriRehd. As well known, KFP is a typically late maturing pear cultivar, which usually matures in middle and late September such as Yali pear [24]. While some early maturing cultivars such as Huangguan pear were usually mature before September [25].Compared with Yali and Niitaka [26], KFP with its lower hardness due to low content of stone cells inflesh is more suitable for people to chew and eat. Another distinctive popular feature of KFP is the rich flavor and volatile compounds including hexanal, ethyl hexanoate,ethyl butanoate, ethyl acetate, hexyl acetate, ethanol,α-farnesene,butyl acetate, and ethyl (E,Z)-2,4-decadienoate [27]. The content of the volatile compounds as above were significantly higher than that of Yali [28]. In addition, KFP peel possessed a relatively higher wax concentration than other popular cultivars such as Nanguoli, Xuehuali and Zaosu. Its higher wax concentration makes it an excellent parent for pear breeding against fruit splitting, pathogenic or insect attacks and mechanical damage [29]. The traditional breeding purpose mainly focused on the cultivation characteristics of resistance against to insect attacks, early or late maturing and increasing yield. While with the improvement of people’s living level, to enhance functional and nutrition values are becoming a new breeding tends. To increase the content of polyphenolic compounds in apples, investigating candidate genes had been realized [30]. Pear fruit is rich in CGA which is one of the main polyphenols, and may becoming an important direction of pear breeding in future.

Our result of the CGA content variation was consistent with the previous report onP. pyrifolia[31]. It was found that the contents of CGA and arbutin in immature pear for this cultivar decreased gradually during the development in fruit. The polyphenolic compounds decrease during fruit development is due to a cultivar of factors, such as the increase of fruit volume and the hydrolysis of polyphenols [32]. Other opinion believed that the reduction of hydroxyl cinnamate and flavanol monomers such as CGA was consumed for the synthesis offlavanols and proanthocyanins [33].

CGA is a major phenolic acid produced by many plants through secondary metabolism. At present, three CGA biosynthesis pathways including a series of key enzymes have been identified in various plants (Fig. 7). Previous studies inLonicera japonica,Saussurea involucrata,Erigeron breviscapusandMorus albashowed that PAL,PTAL, 4CL and C4H had been verified as key enzymes upstream of the CGA biosynthesis [34-37]. Coumarin-CoA was generated from phenylalanine via catalysis by these enzymes. The present study identified two possible pathways for CGA biosynthesis and candidate genes related to CGA biosynthesis in KFP. These pathways were also identified inEucommia ulmoides, tobacco, tomato and other plants [38-40]. However, the HCGQT in the third pathway was only identified in sweet potato root, indicating that this pathway is not the right main pathway for CGA biosynthesis in plants [41].

This study found 57 unigenes encoding 6 enzymes involved in CGA biosynthesis pathway in KFP and 20 of them including 1PAL, 3C4H, 54CLsand 11HCTswere identified as DEGs. Among these DEGs, we found that only onePALwas identified as DEG with relatively low TPM value, though 11 unigenes were annotated asPALand one unigene was annotated asPTAL. In addition, there were 9 unigenes annotated asC4Hwith twoC4HDEGs expressed higher in T1 and lower in T3. In particular, one of the twoC4HDEGs (DN19674_c0_g1) had a pretty higher expression level with the TPM value over 50. Coumaric acid generated from cinnamic acid via catalysis by C4H is a key intermediate product upstream of the CGA synthesis pathway [42]. Therefore, the high expression level of this unigene was the basis of CGA biosynthesis. In addition to participating in CGA biosynthesis, coumaric acid was also an important upstream intermediate product of the phenylpropane biosynthesis pathway, participating in the biosynthesis of various secondary metabolites, such as catechins and lignin [43,44]. So this may also serve as the reason for the high level ofC4Hexpression.These results suggested that those different homologous genes encoding the same enzyme may function differently, which further suggesting that these unigenes upstream of the CGA biosynthetic pathway may be shared by the phenylpropane pathway and other metabolic pathways.

In addition, the expression level of three4CLDEGs (DN15799_c0_g1,DN15799_c0_g2,DN23802_c0_g3) were corresponded well with the CGA content, indicating their essential role in CGA biosynthesis. The other two4CLsbelonging to DEGs had the opposite trend to the change of CGA content, one of which had a pretty higher TEM value, indicating that they may participate in the biosynthesis of other secondary metabolites. Previous studies had shown that the CGA content of nitrogen-deficient tobacco increased significantly,and a group of genes including4CL2were up-regulated required for CGA biosynthesis, indicating that4CL2is a key gene regulating CGA content [45]. Besides, the4CL2gene from tobacco was introduced into the yeastSaccharomyces cerevisiae, demonstrating that yeast can producep-coumaroyl shikimate, which is the precursor of CGA biosynthesis [46].

The downstream route of CGA biosynthesis involved the CGA production from coumaroyl-CoA under the catalysis of C3’H and HCT enzyme. In this study, we also found that 5 unigenes and 25 unigenes were annotated asC3’HandHCT, respectively. C3’H belongs to CYP98A subfamily and is a key enzyme that directly catalyzes the production of CGA byp-coumarinyl quinic acid in rote 2.Moreover, we found that none of theseC3’Hwas identified as DEG, but one of them namedDN1383_c0_g1had a relative higher expression level, and its expression level at T3 was 0.52 times that at T1, consisting with the change trend of CGA content. In the same way, such research found thatC3’Hcorrelated with the high production of CGA in Tartary buckwheat was identified to be in high expression in flowers [47]. These results indicated thatC3’Hwas associated with the accumulation of CGA. Unlike route 2, as the last key enzyme in route 1, HCT could use shikimate acid and quinic acid as substrates to generate CGA. Further research found such a fact thatMaHCT4cloned from mulberry leaves favored quinic acid over shikimic acid as its acyl acceptor [37]. In present study, 25 unigene were annotated asHCT, but only 5 unigene transcripts were associated with CGA content in different development periods. Among them, the sequences of protein encoded by theDN19828_c0_g1andDN14620_c0_g2were highly similar to the shikimateO-hydroxycinnamoyltransferase ofP. bretschneideri, which were considered to be involved in the biosynthesis of CGA. Apart from this, the sequence of proteins encoded by theDN8181_c0_g1was highly similar to the BAHD acyltransferase At5g47980 ofP. bretschneideri. It is well known that HCT is a member of the BAHD superfamily, which is pervasive in land plants likeP. appendiculatumandMarchantia paleacea[48]. However, we found that not allHCTspositively regulated CGA biosynthesis in this study. Therefore, the function of key downstream genes in CGA biosynthesis in KFP need to be verified in future studies.

5. Conclusion

In summary, we conducted RNA-seq andde novotranscriptome analysis on KFP grown in Xinjiang, studied the CGA accumulation model during the development period of KFP, and explored the candidate pathway and genes involved in KFP CGA biosynthesis.Transcriptome analysis showed that there were two potential routes for CGA biosynthesis in KFP, which was involved in‘phenylpropanoid biosynthesis’ pathway. These routes involved 6 key enzymes: PAL, PTAL, 4CL, C4H, HCT, and C3’H. In addition,we found that some unigenes including 1PALs, 2C4Hs, 34CLsand 5HCTswere closely related to the content of CGA, indicating that they played an important role in the biosynthesis of CGA of KFP.

These results revealed the potential mechanism of CGA biosynthesis, and provided a wealth of transcriptome data as the guidance for future breeding and functional component research of pears, especially for KFP grown in Xinjiang.

Conflict of interest

The authors declare there is no conflict of interest.

Acknowledgement

We would like to thank Prof. Liang Liu from Department of Statistics, The University of Georgia, United States for his help on writing and revising the manuscript, Prof. Qikang Gao of Center of Analysis and Measurement, Faculty of Agriculture, Life and Environment Science, Zhejiang University for his assistance in experiment of qPCR analysis. This work was supported by Major scientific and technological projects of XPCC (2020KWZ-012).

Appendix A. Supplementary data

Supplementary data associated with this article can be found, in the online version, at http://doi.org/10.1016/j.fshw.2022.03.007.