Comparative Transcriptomic Analyses Provide Insights into the Enzymatic Browning Mechanism of Fresh-Cut Sand Pear Fruit

Pear (Pyrus spp.) is one of the most commonly consumed temperate fruits, having considerable economic and health importance. Fresh-cut or processed pear fruits are prone to browning because of the abundant phenolic compounds; however, little is known about the molecular mechanisms underlying enzymatic browning of fresh-cut sand pear fruit. In this study, fruits of two sand pear genotypes (low browning cultivar ‘Eli No.2′ and high browning cultivar ‘Weiningdahuangli’) were used to analyze the molecular mechanism of enzymatic browning by SMRT-seq and RNA-seq. The results generated 69,122 consensus isoforms, 21,336 new transcripts, 7105 alternative splicing events, and 254 long non-coding RNAs (lncRNAs). Furthermore, five genes related to enzymatic browning were predicted to be targets of six lncRNAs, and 9930 differentially expressed genes (DEGs) were identified between two different flesh browning cultivars. Meanwhile, most DEGs (e.g., PAL, 4CL, CAD, CCR, CHS, and LAR) involved in the phenylpropanoid biosynthesis pathway were up-regulated, and the expression of PPO and POD were highly expressed in the high-browning cultivar. Interestingly, the transcript level of PbrPPO4 (Pbr000321.4) was significantly higher than other PPO and POD genes, and a high level of total polyphenol and PPO activity were observed in the high browning cultivar. We found that the expression of lncRNA PB.156.1 was significantly positively correlated with the target gene PbrPPO4 (Pbr000321.4). The results suggest that PbrPPO4 might act as a major contributor and a key enzyme encoding gene in regulating fresh-cut sand pear fruit enzymatic browning; the expression of PbrPPO4 was probably regulated by lncRNA PB.156.1. Altogether, the transcriptomic and physiological analyses expand the knowledge of sand pear flesh enzymatic browning at the molecular level and provide a foundation for germplasm resources for molecular breeding of high polyphenol and low browning cultivars in sand pears.


Introduction
Browning of fruits and vegetables during cutting, peeling, packaging, storing, and processing is a common phenomenon, which affects produce quality. This phenomenon, known as enzymatic browning, refers to the oxidation of phenolic compounds into quinones by polyphenol oxidase (PPO) in fruits and vegetables [1]. Pear (Pyrus spp.), one of the most commonly consumed and widespread fruits in the world, has considerable economic and health importance. However, pear fruits are prone to browning during fresh-cut, processing, and postharvest storage, and the browning seriously reduces their appearance, shelf life, taste, flavor, and nutrition [2][3][4].
PPO and peroxidase (POD) are the key enzymes engaged in the polyphenol degradation and browning discoloration in horticultural products [5]. PPO acts on phenols

Plant Materials
For this study, 10-year-old healthy and vigorously grown plants, grafted on P. calleryana, were selected from National Sand Pear Germplasm Repository in Wuhan, Hubei, China (N30. 29, E114.15). Fifteen immature and mature fruits of two sand pear cultivars (P. pyrifolia), the high browning (HB) cultivar (Weiningdahuangli), and the low browning (LB) cultivar (Eli No.2) were obtained. The two samples were both full bloom on 20 March 2017. Immature fruit samples of both cultivars were collected at 90 days after flowering, while mature fruit samples of LB and HB were collected on 28 July 2017 and 14 September 2017, respectively. The mature fruits of HB are oblate, with an average weight of 266 g and an average total soluble solids (TSS) of 11.2%; the mature fruits of LB are obovate, with an average weight of 208 g and an average TSS of 12.3%. Fresh peel and pulp tissues were dissected from immature and mature fruits. In addition, some fruits were cut and kept at room temperature and after 1 h mature pulp tissues were collected. All samples were immediately frozen in liquid nitrogen and stored at −80 • C for subsequent analysis.

Browning Measurement in Sand Pear Flesh Slices
The color changes were measured with a Konica-Minolta CR-5 chromameter (Osaka, Japan) in mature sand pear flesh slices. Three mature fruits and six replicate flesh slices per mature fruit were analyzed. The CIELAB colour parameters (L*, a*, b*) were measured immediately after the sand pear fruit was cut (L 0 *, a 0 *, b 0 *), and then the second reading was taken after 1 h (L 1 *, a 1 *, b 1 *). The browning index (BI Equation (1)) and total color changes (∆E Equation (3)) were calculated as follows. BI = (x − 0.31)/0.172 × 100 (1) where x was calculated as in Equation (2).

RNA Extraction
Total RNA was extracted using a TaKaRa Plant RNA extraction kit (Takara, Kyoto, Japan) following the provided protocol. RNA quality was checked and quantified using a Nanodrop ND-1000 (Thermo Fisher Scientific, Waltham, MA, USA) and the Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA). For SMRT sequencing, equal amounts of RNA from ten different tissues (immature peel, immature flesh, mature peel, mature flesh, and mature fresh-cut at 1 h of HB and LB cultivars) were combined. For Illumina sequencing, the total RNA from mature fruit fresh-cut at 0 h and 1 h of two cultivars was used.

SMRT, Illumina RNA Sequencing, and Data Analysis
Full-length (FL) transcriptome sequencing of SMRT sequencing was carried out at Biomarker Technologies Co., Ltd. (Beijing, China), following the manufacturer's protocols. Size fractionation and selection (1-2 kb, 2-3 kb, and >3 kb) were performed using the BluePippin™ Size Selection System (Sage Science, Beverly, MA, USA). Three SMRT libraries were generated by using Pacific Biosciences' 1.0 template prep kit. The Illumina cDNA libraries were generated by the NEBNext Ultra™ RNA Library Prep Kit (NEB, Beverly, MA, USA) according to the manufacturer's instructions. These processes were repeated for each sample three times.
FL isoforms, open reading frames (ORFs), and alternative splicing (AS) were generated from data assembly. In long noncoding RNA (lncRNA) identification from PacBio sequences, four pervasive coding potential assessment approaches, namely, Coding Potential Calculator (CPC), Coding-Non-coding Index (CNCI), Coding Potential Assessment Tool (CPAT), and Pfam database 27.0 analyses, were combined to sort the non-protein coding RNA candidates from putative protein-coding RNAs in the transcripts.
BLASTX alignment was performed between these unigene sequences and public protein databases: COG (Clusters of Orthologous Groups), eggNOG (evolutionary genealogy of genes: Non-supervised Orthologous Groups), NR (NCBI non-redundant), Pfam (a database of conserved protein families or domains), SWISS-PROT (a manually annotated and reviewed protein sequence database), GO (Gene Ontology), and KEGG (Kyoto Encyclopedia of Genes and Genomes). Transcript expression levels were quantified by the FPKM method (fragments per kilobase of transcript per million mapped reads). Differentially expressed genes (DEGs) with |Fold Change| ≥ 2 and FDR (false discovery rate) < 0.01 were defined after digitally normalized filtering.
KEGG pathway enrichment analysis was used as the unit of significance enrichment analysis. Hypergeometric tests were performed to identify pathways that were significantly enriched in DEGs compared to the whole gene background.

Quantitative Real-Time PCR (qRT-PCR) Analysis
The first-strand cDNA was performed using the PrimeScript™ RT Reagent Kit with gDNA Eraser (TaKaRa) according to the manufacturer's instructions. The qRT-PCR was performed on a CFX real-time PCR system (BioRad, Hercules, CA, USA) using Luna ® Universal qPCR Master Mix (NEB). The relative expression level was calculated using the comparative Ct method (2 −∆∆Ct ). The primers used in this study are shown in Table S1, and the YLS8 gene was used as a reference gene [34]. Each treatment process was duplicated three times.

Measurement of Total Phenolic Content and Monomeric Compounds
The Folin-Ciocalteu spectrophotometric method was used to measure the total phenolics (µg gallic acid/g fresh weight) of the sand pear flesh as previously described by Zheng [21]. The total phenolics content was quantified at 760 nm using a spectrophotometer (Hanon, Jinan, China). Folin-Ciocalteu reagent was supplied from Beijing Solarbio Science & Technology Co. (Beijing, China).The extraction and quantification of chlorogenic acid, catechin, and epicatechin were performed by using the high-throughput and multiplexed LC/MS/MRM method [16,17]. Monomeric phenolic compounds profiling was analyzed using an LC-MS system (1290 Infinity LC and Agilent 6420 QQQ, Agilent Technologies, USA). The detailed information is shown in the Supplementary Materials.

Assessment of Enzyme Activity
The activity of PPO and POD enzymes were detected using the operation manual of Ruixinbio (Quanzhou Ruixin Biological Technology, Quanzhou, China). PPO and POD enzyme activities were measured at 530 nm and 470 nm with a spectrophotometer system, respectively. One unit of enzyme activity (1 UEA) was defined as an increase in absorbance of 0.005 min −1 per g fresh tissue per mL reaction. Five repetitions were performed for each sample.

Flesh Browning Intensity
The flesh from 'Weiningdahuangli' (HB) showed a higher browning intensity compared with 'Eli No.2 (LB) ( Figure 1A). The L* of the HB flesh slices largely decreased more over 1 h than that of the LB slices ( Figure 1B); the HB flesh slices had higher red (a*) and yellow (b*) components than the LB slices ( Figure 1C,D); and a pronounced increase in BI was distinctive in the HB flesh slices ( Figure 1E). To compare the flesh browning variation between two cultivars over time, total color changes (∆E) were determined. HB flesh slices had a more pronounced total color difference (∆E = 11.57 ± 1.95) compared with the LB slices (∆E = 0.69 ± 0.31).

SMRT-and Illumina RNA-seq
To better understand the potential enzymatic browning mechanism of sand pear fruit, we obtained a representative full-length transcriptome to three libraries with different size fractions (1-2 Kb, 2-3 Kb, and 3-6 Kb) on six SMRT cells using the PacBio sequencing platform ( Figure S1). A total of 351,260 reads of insert (ROI) with a mean length of 2702 bp and 175,633 full-length non-chimeric (FLNC) with an average length of 2087 bp were screened (Table S2). In total, 69,122 consensus isoforms were identified with SMRT analysis using ICE (Iterative Clustering for Error Correction) algorithms (Table S3). The 47,797 high-quality isoforms were corrected with the Quiver program in combination with FLNC. We obtained 314 million clean reads of twelve sand pear flesh cDNA libraries with Illumina RNA-seq, and 21,325 low-quality isoforms were rectified using Illumina RNA-seq data (Table S4). A total of 27,610 full-length transcripts were obtained by comparing the consensus isoforms with the reference genome using GMAP. There were 21,336 new transcripts which were found by alternative splicing analysis of full-length transcripts. Using TransDecoder, 20,571 ORFs were identified, of which 16,552 (80.46%) were complete ORFs, and the length of the ORFs is shown in Figure S2. A total of 21,074 new transcripts were blasted to eight public databases, of which 8255 transcripts,12,584 transcripts, 9115 transcripts, 12,346 transcripts, 17,425 transcripts, 14,586 transcripts, 20,562 transcripts, and 17,531 transcripts were annotated in COG, GO, KEGG, KOG, Pfam, SWISS-PROT, eggnog, and NR, respectively (Table S5).

Identification of AS and lncRNAs
AS plays a significant role in controlling plant processes and gene expression [35]. In this study, 7105 AS events were identified in sand pear fruit, of which intron retention was the most common form (IR, 52.81%), followed by alternative 3 splice site (A3SS, 23.04%), exon skipping (ES, 12.17%), alternative 5 splice site (A5SS, 11.18%), and mutually exclusive exon (MXE, 0.8%) (Figure 2A). For example, we discovered that PbrAroD (3-dehydroquinate dehydratase)/SkDH (shikimate dehydrogenase), which was located upstream of phenylpropanoid biosynthesis, had two types of AS events with IR and A3SS to form 12 isomers.
LncRNAs are implicated in the regulation of a diversity of plant developmental processes such as flowering time, reproductive growth, fruit ripening, and abiotic stresses [36,37]. A total of 254 lncRNAs were found from sand pear fruit among all the CPC, CNCI, CPAT and Pfam datasets ( Figure 2B). Of these lncRNAs, 153 lincRNAs, 30 antisense-lncRNAs, 9 intronic-lncRNAs, and 48 sense-lncRNAs were identified.
LncRNA is often used as a signal molecule to regulate the expression of adjacent target genes. The target genes of lncRNA were predicted based on the lncRNA cis target and sequence complementarity. We concentrated on the target genes of lncRNAs involved in enzymatic browning. Five genes related to enzymatic browning were predicted to be targets of six lncRNAs ( Figure 2C

Analysis of DEGs and qRT-PCR Confirmation
Transcriptome profiling identified 9930 DEGs in sand pear flesh through comparison between two cultivars at 0 h and 1 h. Among the DEGs, 4738 and 5192 DEGs were found between LB-0h vs. HB-0h, and LB-1h vs. HB-1h, respectively ( Figure 3A). Twenty-six DEGs were chosen for the qRT-PCR analysis to proof the RNA-seq data; the results showed a significant correlation (r = 0.9566, p < 0.05) and a consistent expression trend between qRT-PCR and RNA-seq ( Figure 3C). A total of 3564 DEGs were expressed in two cultivars at 0 h and 1 h ( Figure 3B), of which nine candidate genes related to phenylpropanoid and flavonoid biosynthesis were up-regulated in HB compared to LB at 1 h ( Figure 3D). For instance, five DEGs (Pbr008363.1, Pbr017290.1, Pbr024635.1, Pbr006709.1, and Pbr022326.1) related to phenylpropanoid biosynthesis were found to encode PAL, C4H, 4CL, shikimate O-hydroxycinnamoyl transferase (HCT), and POD, respectively.

Annotation and Enrichment Analysis of DEGs
In order to understand the function of the DEGs and identify the enzymatic browning related genes in sand pear flesh, all DEGs from the LB and HB cultivar RNA-seq data were annotated in the COG, eggNOG, NR, Pfam, SWISS-PROT, GO, and KEGG databases. In total, 4738 and 5192 DEGs in both cultivars at 0 h and 1 h were successfully annotated using the above mentioned seven databases (Table S6).
According to the GO annotation, 2614 and 2845 DEGs at 0 h and 1 h were assigned to the three main GO functional categories ( Figure S3). These DEGs in the biological classification were mainly involved in the metabolic process, the cellular process, and the single-organism process at both time points. For the cellular components, cell part, cell, and organelle accounted for the majority. In molecular function, the dominant terms were catalytic activity, binding, and transporter activity. To evaluate the metabolic differences between the two cultivars, the biological pathways were analyzed using the KEGG database. Pathway enrichment of DEGs is used to analyze whether there is over-presentation on a particular pathway. The pathway enrichment analysis based on the KEGG database showed that the top three significantly enriched pathways were phenylalanine, tyrosine, and tryptophan biosynthesis (ko00400); protein processing in the endoplasmic reticulum (ko04141) and galactose metabolism (ko00052) in the LB vs. HB cultivars at 0 h (Table S7), and galactose metabolism (ko00052), phenylalanine, tyrosine, and tryptophan biosynthe-sis (ko00400) and ubiquinone and other terpenoid-quinone biosynthesis (ko00130) were significantly enriched in the LB vs. HB cultivars at 1 h ( Figure 4, Table S7).

Phenylpropanoid and Flavonoid Biosynthesis Related Gene Dynamics and Polyphenol Quantification
The gene expression in the RNA-seq data was normalized to FPKM (fragments per kilobase of exon model per million mapped reads) in different browning sand pear fruits flesh, and the expression levels of these metabolic pathway-related genes were calculated using FPKM value.
In the shikimate pathway, the results indicated that most of the phenylalanine biosynthesisrelated genes, one Arod/SkdH gene (Pbr017091.1), two chorismate mutase (CM) genes (Pbr019022.1 and Pbr019996.1), and four arogenate dehydratase (ADT) genes, were upregulated during flesh browning in the HB cultivar compared to the LB cultivar ( Figure 5). In the branch of the shikimate pathway, isochorismate synthase (ICS) was involved in salicylic acid biosynthesis, and most DEGs from shikimic acid to salicylic acid, tryptophan, and tyrosine were down-regulated in the HB cultivar compared to the LB cultivar ( Figure 5). In the phenylpropanoid pathway, one PAL (Pbr008363.1) gene, three 4CL (Pbr012851.1, Pbr027219.1, and Pbr036926.1) genes, two cinnamyl alcohol dehydrogenase (CAD) genes (Pbr026287.1 and Pbr005802.1), one cinnamoyl-CoA reductase (CCR) gene (Pbr022402.1), one chalcone synthase (CHS) gene (Pbr019531.1), one flavanone 3-hydroxylase (F3H) gene (Pbr033501.1), one dihydroflavonol 4-reductase (DFR) gene (Pbr020145.1), one anthocyanidin reductase (ANR) gene (Pbr039269.1), and two leucoanthocyanidin reductase (LAR) genes (Pbr013248.1 and Pbr027485.1) were significantly up-regulated during flesh browning in the HB cultivar compared to the LB cultivar ( Figure 5). Except for these DEGs, phenylpropanoid biosynthesis-related genes such as C4H and HCT were unremarkably up-regulated in the HB cultivar compared to the LB cultivar. To understand the dynamic changes in phenolics metabolism during fruit browning, we examined total and monomeric polyphenol contents. The total polyphenol content in the HB flesh was higher than that of the LB flesh before browning (at 0 h), and the content of total phenol and monomer phenol of the HB cultivar decreased more than that of the LB cultivar in sand pear flesh during browning (Table 1).

PPO and POD Gene Expression and Enzyme Activities
There were two PPO genes and six POD genes that were identified as DEGs, except for one POD gene, in which PPO and POD genes were up-regulated in the HB cultivar. In particular, the expression of one PPO gene (Pbr000321.4), which was annotated as PbrPPO4, was significantly higher than that of the other PPO and POD genes in the HB cultivar at 0 h and 1 h ( Figure 6A). Several PPO and POD genes were subsequently verified by qRT-PCR; this result also suggested that the relative expression of PbrPPO4 (Pbr000321.4) in the HB cultivar was~25-fold higher than that in the LB cultivar at 0 h and~1200-fold higher than another PPO gene (Pbr000322.1) in the HB cultivar ( Figure 6B).
PPO had significantly decreased and POD had significantly increased activity at 1 h in both cultivars compared with 0 h (Figure 6C,D). The results show that PPO activity was higher before flesh browning than after browning in the HB cultivar, and this finding was consistent with the expression trend of PbrPPO4 (Pbr000321.4). Compared to LB fruits, a significantly higher PPO and POD activity was detected in HB fruits ( Figure 6C,D).

Discussion
Full length (FL) transcripts can significantly upgrade the accuracy of genome annotation and assembly, which help to identify full length isoforms of a gene, AS events, lncRNAs, and alternative polyadenylation (APA) sites [38]. Recent advances in comprehensive information with SMRT-seq and Illumina RNA-seq have been applied to detect more gene isoforms and expose functional variety [39,40]. SMRT sequencing was used to obtain diverse FL cDNA sequences and enhance gene annotation in birch-leaf pear (P. betulifolia Bunge) [41], but the error rate of these data was not overcome by correction of Illumina RNA-Seq. In this study, 69,122 consensus isoforms were generated by SMRT from sand pear fruit and corrected by Illumina RNA-Seq. This complete transcriptome was used to detect 7105 AS events and 254 lncRNAs (Figure 2).
Plant lncRNAs are involved in induced responses to various abiotic and biotic stresses [37]. In birch-leaf pear, 251 lncRNAs were related to drought stress responses [42]. In ginkgo, four genes involved in phenylpropanoid biosynthesis targeted by five lncRNAs were identified [43]. We firstly identified 254 lncRNAs in P. prifolia during fruit fresh-cut browning. Five genes related to enzymatic browning were predicted to be targets of six lncRNAs (Figure 2). The target genes were involved in PPO, POD, and phenylpropanoid pathwayrelated genes. It is worth mentioning that two PPO (Pbr000321.4 and Pbr000322.1) targets of lncRNA PB.156.1 were predicted (Figure 2). PPOs have been concerned in several functional processes, including plant defense, regulation of plastidic oxygen levels, and the phenylpropanoid pathway [8]. In bananas, non-coding RNA miR528, which targets PPO genes, plays an important role in the peel browning during cold-stress [13]. Our findings on lncRNAs in sand pears indicate that the expression of lncRNA PB.156.1 was significantly positively correlated with the target gene PbrPPO4 (Pbr000321.4) during flesh browning.
The type and content of phenolic compounds significantly determine the browning characteristics [44]. In pears, cultivars with less phenolic contents were less susceptible to browning than cultivars with higher phenolic contents [45], and the degree of browning was closely correlated with the initial amount of phenolic compounds [46]. The content of phenolic compounds declined during pear slice browning [29], pear wine browning [47], pear internal browning [48], and apple flesh browning [11]. Similarly, we found that the content of total phenolics was higher in the HB cultivar than that examined in the LB cultivar at 0 h before flesh browning. After flesh browning, the concentration of total phenolics was reduced in the HB cultivar at 1 h ( Table 1). The HB flesh probably provided sufficient phenolic substances to PPO after wounding, which caused high browning.
The metabolism of phenolic compounds is also related to the flow of the phenylpropanoid and flavonoid biosynthesis pathways, and these pathways are downstream pathway of the shikimate pathway [27]. From glucose to phenylalanine, most DEGs were up-regulated at the HB cultivar compared to the LB cultivar ( Figure 5). This means that the HB cultivar might synthesize more L-Phenylalanine precursors during sand pear flesh browning in phenylpropanoid and flavonoid biosynthesis pathways, which contain a diversity of key enzymes, such as PAL, C4H, CCR, CHS, HCT, and LAR. In 'Yali' pears, the PbPAL1 and PbPAL2 might regulate core browning under modified atmosphere storage [22]. The expression of MdPAL, MdCHS, MdCHI, MdANS, MdFLS, and MdANR was positively correlated with Fuji apple flesh browning [49]. PAL CHS and 4CL genes were significantly expressed in the brown husk of walnuts [50]. After lettuce cutting, the expression of LsPALs was significantly higher after wounding at 4 h [51]. The expression of PAL encoded genes was largely up-regulated in lotus root pericarp unpacked after 30 days, but the expression was low in those packed with anti-browning solution such as 0.1% ascorbic acid solution, 1% onion solution, and 0.1% citrus powder solution [52]. In our study, most DEGs (e.g., PAL, 4CL, CAD, CCR CHS, and LAR) involved in the phenylpropanoid biosynthesis pathway were up-regulated ( Figure 5). This indicates that phenylpropanoid biosynthesis-related genes plays an important role in protecting plants against wounding-induced abiotic stresses and enhancing the relevant enzymes activity, causing an increase in phenolic compounds during sand pear flesh browning. Actually, in this study, a reduction of the total phenolics was observed during HB pear browning. This suggests that sand pear flesh high browning might depend on the initial number of phenolic compounds and the activation of enzymes.
PPO and POD play an important role in the browning of fresh-cut fruits and vegetables. The methods of reducing enzymatic browning are mainly achieved by inhibiting PPO and POD activity. A decreased expression of enzymatic browning related genes (POD, PPO1, PPO5, and LOX1) and reduced PPO activity were found in pear flesh cut with melatonin treatment at 0.1 mM [21]. The expression of PbPPO1 was correlated with the changes of core browning, phenolic contents, and PPO activity, but the expression of PbPPO4, PbPPO5, and PbPPO6 had no significant relation between core tissue and core browning [22]. PbPPO1 and PbPPO5 were involved in the core browning in 'Yali' pears during cold storage, but the expression of PbPPO4 and PbPPO6 were not correlated with the core browning [3].
In this study, two PPO genes and six POD genes were differentially regulated in LB and HB cultivars. The expression of most PPO and POD in the HB cultivar was significantly higher than that in the LB cultivar, and the expression of PbrPPO4 (Pbr000321.4) was~25fold higher in the HB cultivar than in the LB cultivar at 0 h ( Figure 6). In a transcriptomic analysis of 'Whangkeumbae' pear core browning, the DEGs of PPO4 (Pbr000321.4) were up-regulated under low-temperature storage [53]. The expression level of PPO4 showed discordance between 'Yali' pear [3,22], 'Whangkeumbae' pear [53], and two cultivars of this study with core browning and flesh browning. Therefore, this suggests that PbrPPO4 in different pear cultivars with different types of browning might be regulated independently and have a variety functions. In potato tuber tissue browning, inhibition of PPO by artificial microRNAs (amiRNAs) technology found that the different PPO genes had different contributions to the total PPO protein content, with StuPPO2 as the major contributor. Moreover, the accounting of StuPPO2 and StuPPO1 was~55% and~25-30%, respectively, while the accounting of StuPPO3 and StuPPO4 together was less than 15% [9]. Our results showed that the expression PbrPPO4 (Pbr000321.4) was significantly higher than that of other PPO and POD genes, and the expression trend of PbrPPO4 (Pbr000321.4) was consistent with the trend of PPO activity in the HB cultivar at 0 h and 1 h ( Figure 6). According to these results, PbrPPO4 (Pbr000321.4) might act as a major contributor and a key enzyme encoding gene in regulating sand pear flesh enzymatic browning.

Conclusions
A combination of single molecular real-time (SMRT) sequencing and Illumina RNA sequencing was used to detect FL isoforms, ORFs, AS, lncRNAs, DEGs, and functional annotations. A total of 9930 DEGs were found in sand pear flesh of the high-browning and low-browning cultivars at 0 h and 1 h. The expression profiles of all members of the phenylpropanoid and flavonoid biosynthesis, and the PPO and POD gene families were also analyzed, and many of them had higher expression in the HB cultivar than in the LB cultivar. PbrPPO4 (Pbr000321.4) had a significantly higher transcript level than other PPO and POD genes. A high level of total polyphenol and PPO activity were observed in the high browning cultivar, and the expression of PbrPPO4 declined during sand pear fresh-cut browning in the high browning cultivar as the PPO activity decreased and the intensity of sand pear fresh-cut browning increased. This work provides a transcriptome database and candidate genes for future molecular studies of enzymatic browning and lays a theoretical foundation for molecular breeding of high polyphenol and low browning germplasm resources in sand pears.

Supplementary Materials:
The following are available online at https://www.mdpi.com/article/ 10.3390/horticulturae7110502/s1, Figure S1: Summary of PacBio RS II single-molecule real-time (SMRT) sequencing, Figure S2: Length distribution of ORFs with new transcripts, Figure S3: Gene ontology (GO) categorization of sand pear flesh DEGs. (A) GO analysis of DEGs between low browning (LB) and high browning (HB) at 0 h in three main categories. (B) GO analysis of DEGs between low browning (LB) and high browning (HB) at 1 h in three main categories, Table S1: Primers designed for qRT-PCR validation, Table S2: The PacBio SMRT sequencing information of sand pear fruit, Table S3: ICE (Iterative clustering for error correction) cluster result statistics, Table S4: Samples sequencing data evaluation, Table S5:Annotation of new transcript datasets to eight public databases, Table S6: Summary of annotated differentially expressed genes between LB and HB sand pear flesh at 0 h and 1 h, Table S7: KEGG pathway enrichment of LB vs. HB at 0 h and 1 h in RNA-Seq.