Transcriptome analysis around the onset of strawberry fruit ripening uncovers an important role of oxidative phosphorylation in ripening

Although much progress has been made towards understanding the ripening of non-climacteric fruit using the strawberry as a model plant, the defined molecular mechanisms remain unclear. Here, RNA-sequencing was performed using four cDNA libraries around the onset of ripening, and a total of 31,793 unigenes and 335 pathways were annotated including the top five pathways, which were involved in ribosome, spliceosome, protein processing, plant-pathogen interaction and plant hormone signaling, and the important DEGs related to ripening were annotated to be mainly involved in protein translation and processing, sugar metabolism, energy metabolism, phytohormones, antioxidation, pigment and softening, especially finding a decreased trend of oxidative phosphorylation during red-coloring. VIGS-mediated downregulation of the pyruvate dehydrogenase gene PDHE1α, a key gene for glycolysis-derived oxidative phosphorylation, could inhibit respiration and ATP biosynthesis, whilst promote the accumulation of sugar, ABA, ETH, and PA, ultimately accelerating the ripening. In conclusion, our results demonstrate that a set of metabolism transition occurred during green-to-white-to-red stages that are coupled with more-to-less DEGs, and the oxidative phosphorylation plays an important role in the regulation of ripening. On the basis of our results, we discuss an oxidative phosphorylation-based model underlying strawberry fruit ripening.

The regulation of strawberry fruit ripening involves many phytohormones, transcription factors, cis-elements, and metabolic enzymes. This suggests that the ripening-regulated molecular mechanisms of non-climacteric fruit, in contrast to those of climacteric fruit regulated by ethylene, is highly complex. Given that one obvious sign of the onset of strawberry fruit ripening is the development of red colouration, transcriptome-based RNA sequencing carried out at the reddening stage may provide insight into the molecular mechanisms of non-climacteric fruit ripening. We constructed and sequenced the cDNA libraries of receptacles at the large green, white, initial-red and partial-red stages using the Illumina HiSeq 2000 platform. Based on our results, we outline an oxidative phosphorylation-based model for the onset of strawberry fruit ripening.

Results
RNA sequencing analysis and de novo assembly. To maximize the transcriptional data obtained from ripening strawberries, we prepared four mixed cDNA libraries from ten uniformly sized fruits at each of four stages of ripeness for RNA sequencing using the Illumina HiSeq 2000 instrument 62 . Each sequenced sample yielded 100 bp reads from paired-end sequencing of cDNA fragments. After quality assessment and data clearance, 5.58-7.71 billion (G) reads with more than 90% Q20 bases (i.e., an average base quality greater than 20) were retained as high-quality reads for each library and used in subsequent analyses (Supplemental Table 1). The average 'G + C' content of the above was 48% [46%, 50%, 49%, and 48% for CM1 (large green), CM2 (white), CM3 (initial-red), and CM4 (partial-red) libraries, respectively] (Supplemental Table 1).
All valid reads were combined to perform de novo assembly by the paired-end method with Trinity software 63 (trinityrnaseq_r20131110 version). A total of 98,848 unigenes were obtained, among which 24,358 genes were longer than 1 kb. An overview of the assembled transcripts and unigenes are presented in Supplemental Table 2. The length distributions of unigenes are shown in Supplemental Figure 1. These results demonstrate the effectiveness of Illumina pyrosequencing in rapidly capturing a large portion of the transcriptome. Sequence annotation. We used several complementary approaches to annotate the assembled sequences.
The unigenes were compared against diverse protein databases, including the National Center for Biotechnology Information (NCBI) non-redundant protein (NR) database, SWISS-PROT, TrEMBL, Cdd, pfam, KOG database, Kyoto Encyclopedia of Genes and Genomes (KEGG), and Gene Ontology (GO). The resulting functional annotations are listed in Supplemental Table 3 and provide the best unigenes with an E-value of 1e −5 and identity of 30%. We performed a sequence similarity search against the NCBI NR and SWISS-PROT protein databases using the Basic Local Alignment Search Tool (BLAST) algorithm, specifying an E-value of less than 1e −5 . This analysis indicated that of the 98,848 unigenes, 41,282 (41.76%) matched the NR database and 25,312 (25.61%) had similarity to proteins in the SWISS-PROT database. A total of 14,267 unigenes classified into 25 KOG annotations are provided in Supplemental Figure 2. The top three functional categories were signal transduction mechanisms, posttranslational modification/protein turnover/chaperones, and general function prediction only. Altogether, 43,769 (44.3%) unigenes were successfully annotated in NR, SWISS-PORT, TrEMBL, Cdd, pfam, KOG, KEGG, and GO. The high percentage (55.7%) of unmapped unigenes that could be assigned a putative function may be due mainly to the short sequence reads generated by the sequencing technology and the relatively short sequences of the resulting unigenes, most of which probably lack the conserved functional domains. Another possible explanation may be that some unigenes were non-coding RNAs or unknown.
Based on the BLAST UniProt annotation, we carried out GO analysis, which led to the annotation of 31,793 unigenes. The results comprise three classifications: "biological process, " "cellular component", and "molecular function. " The most frequently associated GO terms for the biological process group were "metabolic process" (GO:) and "cellular process" (GO:); for the cellular component group, they were "cell" (GO:) and "cell part"; and for the molecular function group they were "binding" (GO:0005488) and "catalytic activity" (GO:0003824) (Supplemental Figure 3).
To analyze the KEGG pathway of unigenes, we used KAAS to obtain corresponding KO numbers. These KO numbers were then used to obtain corresponding KEGG pathways. By analysing the relationship between unigenes and enzymes in KEGG comment files and mapping them to the pathway, we created 335 pathway charts. This process predicted a total of 301 pathways represented by a total of 14,011 unigenes (Supplemental Table 4). The top five pathways were involved in ribosomes, spliceosomes, protein processing, plant-pathogen interaction, and plant hormone signal transduction (Supplemental Table 5).
Gene expression patterns throughout fruit development. We used the distribution of gene expression levels to evaluate the normality of the library data. Gene expression level was determined by calculating the number of reads and normalising to the RPKM. As shown in Fig. 1, the majority of mRNAs were expressed at low Scientific RepoRts | 7:41477 | DOI: 10.1038/srep41477 levels, whereas a small proportion was highly expressed. In CM1, most genes (10,307) have expression levels close to 1 (< 1, namely log 2 [2 RPKM]), while only 0.07% of the genes (27 of 39,086) were highly expressed at a level above 10. Compared with CM1, the distribution of the CM2, CM3, and CM4 libraries showed similar expression patterns; in other words, the later stages (white, initial-red, partial-red) were distinct from the large green stage, suggesting that the white stage might represent the onset of fruit ripening.
To analyze variation in gene expression throughout strawberry fruit development, the differentially expressed genes (DEGs) between each library pair (CM1-CM2, CM1-CM3, CM1-CM4, CM2-CM3, CM2-CM4, and CM3-CM4) were listed in Fig. 2.The results showed that more DEGs were detected when comparing CM2-CM4 to CM1, and among these DEGs, downrregulated genes were predominant, indicating that the onset of fruit ripening is coupled with many genes that were downregulated. Comparing CM3/CM4 with CM2, DEGs decreased sharply and upregulated genes were predominant, showing that these upregulated genes may therefore be involved in the onset of ripening. Notably, CM3 and CM4 had similar gene expression patterns. Some important DEGs related to ripening were listed in the Heatmap Fig. 3, mainly included oxidative phosphorylation (NADH, ATPase), plant hormones (ABA, IAA, GA, ETH, JA and PA), antioxidation (GST, Vc), protein translation and processing, sugar metabolism, pigment and softening. Taken together, these results suggest that a set of metabolism transition occurred during green-to-white-to-red stages, and the processes of fruit de-greening are more complex than those of fruit red-colouring, and that white fruit is a distinct stage.
Enrichment pathway analysis of DEGs. We performed KEGG pathway enrichment analysis to categorise the biological functions of DEGs. We mapped all the genes to terms in the KEGG database. Comparing CM1 with CM2, CM3, and CM4, 41, 52, and 54 pathways were upregulated while 42, 48, and 54 pathways were downregulated, respectively. Comparing CM2 with CM3 and CM4, 36 and 39 pathways were upregulated while 42 and 43 pathways were downregulated, respectively. Comparing CM3 with CM4, 27 pathways were upregulated while 30 pathways were downregulated. Some DEGs specific to each stage are listed in Fig. 4, compared CM2/CM3/ CM4 with CM1, 828, 930 and 1312 DEGs were found, respectively (Fig. 4A); compared CM3/CM4 with CM2, 284 and 483 DEGs were found, respectively (Fig. 4B); compared CM4 with CM3, only 192 DEGs were found.  The ten most DEG-rich pathways between library pairs are provided in Supplementary Table 6. Between CM1 and CM2, protein turnover and processing increased within pulp cells, as well as including oxidative phosphorylation, carbohydrate metabolism and the plant hormones (ABA, JA). The carbohydrate metabolism pathways included carbon fixation, glycolysis/gluconeogenesis, pentose and glucuronate interconversions, as well as starch degradation, sucrose accumulation, amino sugar, and nucleotide sugar metabolism related to fruit softening. Meanwhile, DNA turnover, photosynthesis, cell cycle, emergency reactions, and the plant hormones (IAA, GA, ETH and PA) signaling decreased within pulp cells.
Between CM1 and CM3/CM4, pulp cells maintained high levels of protein turnover and processing, oxidative phosphorylation, carbohydrate metabolism and the plant hormones (ABA, JA). Glutathione, phenylpropanoid/ flavonoid, cysteine, and methionine metabolism were also active. Notably, cysteine and methionine metabolism were involved in ETH and PA biosynthesis. Carbon fixation, starch degradation, plant-pathogen reaction, and Comparing CM2 to CM3/CM4, pulp cells sustained high levels of protein/DNA/RNA turnover and processing, as well as glutathione, flavonoid, sucrose metabolism and the plant hormones (ABA, JA, ETH and PA). Meanwhile, photosynthesis, carbon fixation, and the plant hormone (GA and IAA) signal transduction, as well as cysteine, methionine, and starch metabolism, decreased. Notably, oxidative phosphorylation metabolism also decreased.
Comparing CM3 to CM4, pulp cells sustained high levels of phenylpropanoid/flavonoid biosynthesis, carbohydrate metabolism (involved in pectate lyase, pectinesterase and beta-glucosidase), cysteine and methionine metabolism and the plant hormones (ABA, JA, ETH and PA). Porphyrin and chlorophyll metabolism, related to chlorophyll degradation, were also sustained at high levels. Photosynthesis, carbon fixation, and the plant hormone (GA and IAA) signal transduction decreased, as well as starch metabolism and oxidative phosphorylation.
Together, these results show that the processes of strawberry fruit ripening involve a decrease in oxidative phosphorylation and the plant hormone (GA and IAA) signaling, while an increase in carbohydrate, phenylpropanoid, flavonoid, glutathione, methionine and the plant hormone (ABA, JA, ethylene, and PA) signaling. The fact that the oxidative phosphorylation decreases at red-coloring, suggests a regulatory role of oxidative phosphorylation in ripening.

Downregulation of PDHE1α by VIGS promotes strawberry fruit red-colouring. Given that
the white stage represents the onset of fruit ripening and only several genes have an increased expression in this staged, including NADH dehydrogenase (ubiquinone) flavoprotein 1, pyruvate dehydrogenase E1 component alpha, heat shock 70 kDa protein 1/8 and EIN3-binding F-box protein (Fig. 3). It is known that pyruvate dehydrogenase contributes to linking the glycolysis metabolic pathway with the citric acid cycle and releasing energy via NADH dehydrogenase-mediated oxidative phosphorylation, and its E1 subunit is considered to be the rate-limiting step for the pyruvate dehydrogenase 64 .
To further explore the role of oxidative phosphorylation in strawberry fruit development, we used the transcriptome data to investigate the mRNA expression levels of pyruvate dehydrogenase. The results showed that the E1 component subunit alpha of pyruvate dehydrogenase gene, PDHE1α (GenBank NO. XM_004297634.2), had an increased expression in white stage, and then followed a high-to-low pattern during strawberry fruit ripening (Fig. 5A). In order to investigate a role of the PDHE1α in the fruit ripening, the PDHE1α gene expression was silenced using syringe-infiltrated a mixture of Agrobacterium strain GV3101 cultures containing pTRV1 and pTRV2 carrying a 589-bp fragment of the gene (pTRV2-PDHE1α ) into large green fruit in a 1:1 ratio. Control fruits were infiltrated only with TRVs alone. One week after infiltration, control fruit remained white (Fig. 5B); in contrast, RNAi fruit turned red (Fig. 3C), and PDHE1α transcripts were downregulated compared to the control (Fig. 5D). The results that downregulation of PDHE1a expression stimulated the fruit ripening, making PDHE1a a negative regulator of strawberry fruit ripening.

Alteration of PDHE1α expression affects the ripening-related physiological parameters.
To further elucidate the role of PDHE1α in the regulation of strawberry fruit ripening, we assayed ripening-related physiological parameters in both RNAi and control fruit, including plant hormones (ABA, GA, IAA, JA, ETH, and PA), soluble solids (sugar) concentration, biochemical measurements of oxidative phosphorylation (ATP), and respiration rate. The levels of ABA, ETH, PA, and soluble solids concentration were all remarkably high in RNAi fruit compared with control fruit (Table 1). In contrast, the levels of IAA, GA 3 , and JA showed non-significant changes between the RNAi and control fruit. Notably, both ATP content and respiration rate were markedly lower in RNAi fruit compared with control fruit. These results suggest that downregulation of PDHE1α expression inhibits respiration and ATP biosynthesis, resulting in an increase in soluble solids concentration and an accumulation of ABA, ETH, and PA, which ultimately promotes the ripening in RNAi fruit.

Discussion
While the molecular mechanisms of strawberry fruit ripening are highly complex, the external signs (de-greening and red-colouring) are very simple. On the basis of our transcriptome-based RNA-sequencing analysis around the onset of ripening using large green, white, initial-red, and partial-red fruit, the previous and our data indicate that the pyruvate dehydrogenase E1α (PDHA1) is the first component enzyme of the pyruvate dehydrogenase (PDH) complex that catalyzes the overall conversion of pyruvate to acetyl-CoA, which is subsequently used in cellular respiration by both the citric acid cycle and oxidative phosphorylation to generate ATP, so this E1α subunit is considered to be the rate-limiting step for the pyruvate dehydrogenase complex 64,65 . The downregulation of PDHE1α expression inhibits respiration, ATP biosynthesis, and ABA and ETH accumulation (Fig. 5 and Table 1), to some extent, suggesting that a link between oxidative phosphorylation and strawberry fruit ripening is existed. Further explain is as follows: It is known that oxidative phosphorylation play a central role in supplying the carbon skeleton and motive force for biochemical reactions, and that pyruvate dehydrogenase is a key step in the oxidative phosphorylation stage of sugar metabolism 66 . Our RNA-sequencing analysis indicated that the mRNA expression level of the E1 component subunit alpha of the pyruvate dehydrogenase gene PDHE1α follows a trend similar to that of oxidative phosphorylation activity at the onset of strawberry fruit ripening (Fig. 5). Interestingly, the downregulation of PDHE1α expression could inhibit respiration and ATP biosynthesis, whilst also promoting the accumulation of sugar, ABA, ETH, and PA, ultimately accelerating fruit ripening (Table 1). These results are consistent with previous reports on the roles of ABA, ETH, and sugar in the regulation of strawberry fruit ripening 2,7-12,67,68 .
On the one hand, glycolysis-derived oxidative phosphorylation provides the required carbon and energy fluxes for fruit development; on the other, cellular respiration through oxidative phosphorylation exhausts carbohydrates that provide the basis for fruit development. In the present study, the number of the differentially expressed genes (DEGs) is remarkable lower after white stage (comparing CM3/CM4 with CM2) in comparison to de-greening stages (comparing CM3/CM4/CM2 with CM1; Fig. 2, Supplemental Table 6), suggesting that the processes of early development are more complex (more DEGs) before white stages, once on the set of red-coloring in white fruit, the processes of ripening is relative simple (less DEGs). To a large extent, a metabolism transition occurred during green-to-white-to-red stages. In response to this change, fruit decreased oxidative phosphorylation intensity, which in turn promoted the accumulation of sugar, ABA, and ETH (Table 1), which can accelerate fruit ripening 7,67 .
On the basis of the available data, we have tentatively outlined the processes of strawberry fruit ripening at the molecular level. With the action of sun energy, CO 2 and H 2 O can be transformed into sugar by carbon fixation in the chloroplasts of leaves and fruit. Sugar can be turned into CO 2 and H 2 O by oxidative phosphorylation in order to produce the motive force for biochemical metabolic reactions in fruit. During fruit development, protein/ DNA/RNA turnover and oxidative phosphorylation follow a decreasing trend that corresponds to DEG variation ( Fig. 2 and Table 1). It is interesting to note that, most upregulated DEGs were involved in sugar accumulation, fruit softening, and pigment biosynthesis (phenylpropanoid/flavonoid metabolism), while the capacity of fruit for carbon fixation and resistance decreased (Fig. 3 and Table 1). In response to these changes, fruit cells accelerate lysosome action and evoke cellular antioxidant action through glutathione and methionine metabolism 20,69 . Undoubtedly, cellular metabolism is regulated by plant hormones in response to fruit developmental cues: coupled with the onset of the fruit ripening, the action of ABA, ethylene, and polyamine increased (Fig. 6).   The present study provides the first description of an oxidative phosphorylation-based model for strawberry fruit ripening based on transcriptome data. Given the complexity of non-climacteric fruit ripening, the roles of aerobic respiration and anaerobic respiration in the regulation of fruit ripening are subjects for future study.
RNA extraction, library construction, and RNA-seq. Three fruits were randomly taken out from the frozen samples for every RNA isolation and cDNA synthesis. RNA was extracted using the RNeasy plant mini kit (Qiagen, Dusseldorf, Germany) from the receptacles of CM1-CM4 fruit. DNase digestion with the RNase-Free DNase set (Qiagen) was performed to remove contaminating DNA, and the RNA samples were then processed using the RNA library prep kit (New England BioLbs, Ipswich, MA, USA) (NEB) and sequenced with the Illumina HiSeq2000 platform. The experiments were repeated three times.
Sequence data analysis and assembly. The raw reads were filtered with the FASTQ_Quality_Filter tool from the FASTX-toolkit. Reads with more than 35 bp having a quality score higher than 20 were kept 70 . All valid reads were combined to perform de novo splicing by the paired-end method with Trinity software 71 . The longest transcripts per locus were used as a unigene.
Functional annotation. Several complementary approaches were utilised to annotate the unigenes.
Searches were conducted using the Basic Local Alignment Search Tool (BLASTX) 72 . The unigenes were compared against the NCBI NR, SWISS-PROT, TrEMBL, Cdd, pfam, and KOG databases with an E-value of 1e-5 and identity of 30% 73 . Functional annotation by gene ontology terms (GO, http://www.geneontology.org) was carried out using the Blast2GO software 74,75 . To achieve KEGG orthology (KO) assignment 76 , the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways were assigned to sequences using the online KEGG Automatic Annotation Server (KAAS) (http://www.genome.jp/kegg/kaas/) 77,78 . The output of KEGG analysis includes KO assignments and KEGG pathways 77,78 .
The bowtie2-2.2.2 software was used to compare reads with unigenes using the single-end mapping method 71,79 . To compare the levels of unigene expression among the four libraries, the transcript level of each expressed unigene was calculated and normalised to the reads per kilobase of the exon model per million mapped reads 80 (RPKM). The significance of differential unigene expression was determined by using the Chi-squared test with a threshold of P = 0.05. P-values were adjusted to account for multiple testing by using the false discovery rate (FDR) and assigned error ratio Q-values (< 0.05) 81 . The unigenes with an adjusted P-value of < 0.05 and an absolute value of log 2 (expression fold change) > 1 were deemed to be differentially expressed, while the unigenes with an FDR-adjusted P-value of < 0.05 were considered statistically significant. PDHE1α was amplified using primers (sense 5′ -TCTAGAACCACTGCACCTTCCTCG -3′ , underlined is Xba І, and antisense 5′ -CTCGAGTCACGCTCCTG TCTCACG -3′ , underlined is Xho І) cloned into the virus vector Xba І-Xho І-cut pTRV2. Agrobacterium strain GV3101 containing the pTRV1, pTRV2, or pTRV2 derivatives pTRV2-PDHE1α 589 was used for RNAi (Fig. 5). Agrobacterium-mediated TRV infiltration by syringe injection with a needle into strawberry fruits was performed as described by Jia et al. 7 . (2011). For semi-quantitative RT-PCR (SqRT-PCR) analysis of PDHE1α transcripts, the cDNA was used as a template for PCR amplification through 25 cycles and actin was used as a reference gene (sense 5′ -CAGTTAGGAGAACTGGGTGC-3′ and antisense 5′ -TGGGTTTGCTGGAGATGAT-3′ ). The experiments were repeated three times.
Determination of fruit respiration, ethylene release rate, and soluble solid content. 200 g of fruit were sealed in a 500 mL glass bottle at 25 °C. After 1 h, 10 mL of gas sample were used to measure respiration rate with an infrared carbon dioxide analyser (AMETEKMode lCD-3A) and ethylene release rate with a HP P5890a gas chromatograph. The soluble solids content of flesh was measured using a hand-held sugar measurement instrument (MASTER-100H, ATAGO Master, Japan), onto which fruit juice was applied to obtain a reading. Ten uniformly sized fruit were used for analysis. The experiments were repeated three times. ATP assay. ATP synthesis rates were measured using a luminometer (Glomax 96 Microplate Luminometer; Promega). 20 mL of isolated pulp chromoplasts were added to microplate wells containing 80 mL of luciferase/ luciferin reagent (ENLITEN; Promega) and 80 mL of buffer (600 mM sorbitol, 10 mM TES, 2 mM MgCl2, 25 mM KH2PO4, and 0.33 mM EDTA, pH 7.4). Ten uniformly sized fruit were used for ATP synthesis measurement, and the experiments were repeated three times.
Hormone analysis. In order to determine ABA, IAA, GA, and JA contents, 0.5 g receptacle samples were ground in a mortar and homogenised in extraction solution (80% methanol, v/v). Extracts were centrifuged at 10,000 × g for 20 min. The supernatant liquid was eluted through a Sep-Pak C18 cartridge (Waters, Milford, MA, USA) to remove polar compounds, and then stored at − 20 °C for enzyme-linked immunosorbent assay (ELISA) according to the manufacturer's instructions (China Agricultural University, Beijing, China). Ten uniformly sized fruit were used for detection, and the experiments were repeated three times.
Scientific RepoRts | 7:41477 | DOI: 10.1038/srep41477 Polyamines assay. The polyamines assay was performed using high performance liquid chromatography (HPLC) with the Waters ™ 996 photodiode Array Detector, Millennium 32 Chromatography software, and Diamond C18 columns (250 × 4.6 mm, 5 μ m). The wavelength of detection was 230 nm, the flow rate was 1.2 mL/min, and and mobile phase was acetonitrile:methanol:water (3:4:3, v-v:v). 3 mL of precooled 5% perchloric acid (V/V) were added to 1 g samples, and the mixture was ground to homogeneity at 4 °C and centrifuged (15,000 rpm, 30 min, 4 °C) after 1 h extraction. 1 mL of supernatant, 7 μ L of benzoyl chloride, and 1 mL 2 M NaOH solution were then added to 10 mL plastic centrifuge tubes. After 20 s of vortexing, the mixtures were incubated at 37 °C for 20 min in a water bath. 2 mL of saturated NaCl solution were added to the mixtures, which were then extracted with 2 mL diethyl ether after mixing. The mixture was centrifuged at 1,500 rpm for 5 min and 1 mL of ether phase was collected and dried under vaccum. The residuum was dissolved with 200 μ L methanol and filtered through a 0.45 μ m membrane, and 10 μ L of filtrate were collected for injection. Ten uniformly sized fruit were used for analysis, and the experiments were repeated three times.