Comparative transcriptome analysis of two contrasting kiwifruit (Actinidia) genotypes under waterlogging stress

Kiwifruit vines are generally sensitive to waterlogging stress. So far, molecular responses of different kiwifruit genotypes for waterlogging stress are less well-explored. In this study, using RNA-sequencing, we examined transcriptional regulation in the roots of a waterlogging-tolerant genotype KR5 (Actinidia valvata), and a sensitive genotype ‘Hayward’ (Actinidia deliciosa) subjected to 0, 12, 24, and 72 h of waterlogging. Compared with 0 h, transcriptional adjustments of these two genotypes occurred as early as 12 h and became notably pronounced 72 h after waterlogging. Waterlogging stress for 72 h promoted the expression of genes involved in ethylene biosynthesis, sucrose and hexose transport, anaerobic fermentation, nitrate reduction, alanine accumulation, and reactive oxygen scavenging in both genotypes. The differential regulation of genes encoding 9-cis-epoxycarotenoid dioxygenase, phosphoglucomutase, alanine-glyoxylate transaminase, and other enzymes pointed to their diverse strategies upon waterlogging in these two genotypes. In addition, more sucrose and trehalose contents, as well as a higher activity of alcohol dehydrogenase and manganese superoxide dismutases were stimulated in KR5 roots after 72h of waterlogging than that in ‘Hayward’. Overall, our results provided more insights into the molecular basis of the waterlogging response in kiwifruit.


Introduction
The flooding events have occurred with increasing frequency under the influence of global climate change [1]. The sustainability of irrigated agriculture has been threatened by this stress for centuries [2]. One direct challenge induced by waterlogging is the decrease of oxygen concentration in the rhizosphere [3]. Under long-term waterlogging, mitochondrial respiration of root is restricted followed by an energy crisis that can ultimately kill the plant [4]. The plants adapt themselves to flooding stress by altering their gene expression and metabolism. Especially, some tolerant species or their genotypes have developed 'better' strategies to respond to different water conditions than sensitive ones [5]. In this regard, comparative transcriptome has been widely used to examine the molecular responses to waterlogging stress in many plant species [6][7][8][9].
Oxygen sensing is responsible for the acclimation responses thereafter. Under hypoxia, coupled ethylene and nitric oxide (NO) signalings mediate the stabilization of the subgroup VII ETHYLENE RESPONSE FACTORs (ERF-Ⅶ ) proteins [10,11]. These proteins are principle activators of many hypoxia-responsive genes [12]. In many plants, genes involved in signal transduction pathways, transcription, fermentation, and sucrose-related pathways are conservatively modulated in

De Novo Transcriptome Sequencing and Read Assembly
In total, 24 cDNA libraries from root tissues of KR5 and 'Hayward' were generated using Illumina HiSeq. The sequencing resulted in 20 ~ 26 million clean reads with an average GC content of 45.03% ~ 47.62% per library. After assembling the clean reads, 496,799 transcripts were yielded for KR5 with an N50 value of 744 bp, as well as 227,946 unigenes with an N50 value of 822 bp. While, 487,169 transcripts were produced in 'Hayward' with an N50 value of 704 bp, as well as 236,597 unigenes with an N50 value of 772 bp. All unigenes were longer than 200 bp and the length distribution was shown in Figure S1.

Unigene Annotation and Functional Classification
According to the Nr (NCBI non-redundant protein sequences), 17.22% ~ 18.55% of Vitis vinifera sequences matched with unigenes of these two kiwifruit genotypes (Figure 2A). The KOG analysis showed that "general function prediction only" represented the largest group followed by two categories related to "posttranslational modifications, protein turnover, and chaperones" and "signal transduction mechanisms" (categories O and T; Figure 2D). The other enriched metabolism categories were partly involved in amino acid, carbohydrate, and lipid transport and metabolism (categories E, J and I; Figure 2D).

The W 72h Induced a Large Number of DEGs in Roots of KR5 and 'Hayward'
Compared with W 0h, the number of differentially expressed genes (DEGs) increased gradually as the waterlogging treatment proceeded in KR5 ( Figure S2). While, the number of active genes was very low after W 12h and increased during W 24h to W 72h in 'Hayward'. By comparing the DEGs between each neighbor time-point, the number of DEGs reached a maximum between W 0h -W 12h and between W 12h -W 24h in KR5 and 'Hayward', respectively. As W 72h induced a great number of DEGs than other treatments, transcript adjustments at this time point might be important for these two genotypes.

The GO Functional Analysis of Waterlogging-Responsive DEGs
The GO analysis of the common 6113 waterlogging-responsive DEGs in KR5 showed that the go terms involved in phosphorylation and phosphotransferase activity were highly enriched ( Figure  3). Analysis of the common 369 waterlogging-responsive DEGs in 'Hayward' indicated the go terms "carbohydrate metabolism process", "steroid biosynthetic process", and "metal ion binding" were dramatically enriched ( Figure 4). In comparison, "carbohydrate metabolism", "the phosphorylation processes" and "oxidoreductase activity" were activated in both genotypes, though they were enriched at different time points during waterlogging. Specifically, the go terms "nucleotide-sugar metabolic process", "protein metabolic process", and "primary metabolic process" were solely enriched in KR5. While, go terms like "hydrolase activity, hydrolyzing O-glycosyl compounds" and "ammonia ligase activity" were specifically enriched in 'Hayward'.

K-Means and KEGG Analysis of Waterlogging-Responsive DEGs
Among the common 6113 waterlogging-responsive DEGs in KR5, 60.1% of them showed up-regulation trend of transcript abundance over the time course analyzed (subcluster 1 and subcluster 2, Figure 5A). The expression levels of these DEGs started to increase as early as after W 12h. The KEGG analysis of subcluster 1 showed that carbohydrate, amino acid, fatty acid, and alkaloid metabolism pathways were highly enriched ( Figure 5B). Several down-regulated genes were especially enriched in subcluster 4. The KEGG analysis of subcluster 4 indicated that the pathways like "sulfur metabolism", "phenylpropanoid biosynthesis", "pentose and glucuronate interconversions", "cyanoamino acid metabolism", and "ABC transporters" were enriched.
Among the common 369 waterlogging-responsive DEGs in 'Hayward', 95.7% of them were all along up-regulated across the whole waterlogging period (subcluster1, 2, and 3, Figure 6A). An early induction of these up-regulated DEGs also occurred after W 12h. The KEGG analysis of the highly expressed genes in both subcluster 2 and subcluster 3 showed that "glycolysis", "carbon fixation in photosynthetic organisms", "carbon metabolism", "biosynthesis of secondary metabolites", and "biosynthesis of amino acids" were enriched pathways in response to waterlogging ( Figure 6B). Besides, the pathways named "fatty acid degradation", "tyrosine metabolism", and "alpha-linolenic acid metabolism" were enriched with down-regulated DEGs.

Oxygen-Sensing Responses in KR5 and 'Hayward' After W 12h
In this study, transcript patterns of some oxygen-sensing genes after W 12h in kiwifruit were found similar to that in Arabidopsis thaliana [10,11]. Our results showed that the expression of most PHYTOGLOBINs (PGBs) was significantly up-regulated after W 12h in both KR5 ( Figure 7) and 'Hayward' (Table S1) under the possible effects of ethylene accumulation. Then, to examine the expression of ERF-VIIs in kiwifruit, we searched the assembled unigenes using the HMM profile of the Apetalla2 (AP2) / ERF domain. In total, 186 and 200 putative AP2 / ERF unigenes were identified in KR5 and 'Hayward', respectively ( Figure S3). By comparing Arabidopsis ERF-VIIs and kiwifruit ERFs, one and two putative kiwifruit ERF-VIIs were found in KR5 and 'Hayward', respectively. The transcript levels of TRINITY_DN56147_c1_g1 from KR5 and TRINITY_DN80164_c1_g2 from 'Hayward' were increased after W 12h, but neither was found significant at p ≤ 0.05. A significant increase of PCOs was found after W 12h in KR5 (Figure 7), but not in 'Hayward'. Unigenes encoding the fermentative enzymes ALCOHOL DEHYDROGENASE (ADH) were activated after W 12h in both genotypes.

Differential Regulation of Genes Related to Plant Hormone Signal Transduction After W 72h
The ethylene plays an important role in plant waterlogging responses. The rate-limiting biosynthetic genes of it including 1-aminocyclopropane-1-carboxylate synthase (ACS) and ACC oxidase (ACO) were mostly up-regulated after W 72h in both genotypes (Table S2). The W 72h treatment also promoted the expression of ethylene receptor genes (ETRs) in KR5. In some species, the accumulated ethylene negatively regulated the abscisic acid (ABA) level [38]. In KR5, unigenes encoding 9-cis-epoxycarotenoid dioxygenase (NCED), a rate-limiting enzyme in abscisic acid (ABA) biosynthesis, were inhibited after W 72h (Table S2). Moreover, 8 out of 16 unigenes encoding ABA 8'-hydroxylase, an enzyme which hydroxylates ABA, were enhanced after W 72h in KR5. Both NCEDs and ABA 8'-hydroxylase encoded unigenes were mostly promoted after W 72h in 'Hayward'. The expression levels of PYR/PYL (ABA receptors) were mostly enhanced in the roots of these two genotypes. DELLA proteins, repressors of gibberellin (GA), showed decreased levels of transcript in both genotypes. In contrast, transcript levels of BZR1BRASSINAZOLE-RESISTANT1 (BZR1), a regulator of a host of BR-responsive genes, significantly increased after W 72h in 'Hayward' while declined in KR5.

KR5 Accumulated More Sucrose and Trehalose Than 'Hayward' After W 72h
Sucrose content significantly increased after W 72h in KR5 while decreased in 'Hayward', compared with that of W 0h ( Figure 8A). Suc-phosphate synthase (SPS) and Suc-phosphate phosphatase (SPP) are responsible for sucrose synthesis. Transcript levels of these two enzymes declined in KR5 while no significant changes for them were found in 'Hayward' ( Figure 8B, Table  S3). Interestingly, expressions of unigenes encoding sucrose transporters (SUT) and hexose transporters (STP) were mostly increased in both genotypes. Cell wall invertase (CWI) is also involved with the phloem unloading of sugars [39,40]. All eight DEGs encoding CWI were significantly increased in KR5 while no significant change was found in 'Hayward'. Sucrose synthase (SUS) and sucrose invertase (INV) are involved in sucrose cleavage. The W 72h treatment promoted the expression of SUSs in both genotypes. Among the INVs, transcriptions of vacuolar invertase (VIN) encoded unigenes were inhibited while most neutral invertase (NI) encoded unigenes were enhanced in these two genotypes.
The starch was also accumulated after W 72h in waterlogged roots of both genotypes ( Figure  8A), but 'Hayward' contained more starch content than that of KR5. The ADP-glucose pyrophosphorylase (AGPase) catalyzes the first committed step in starch biosynthesis. However, unigenes encoding AGPase were mostly down-regulated in these two genotypes ( Figure 8B, Table  S3). Furthermore, the transcripts associtated with starch synthesis e.g. granule-bound starch synthase (GBSS), soluble starch synthase (SSS), and starch branching enzyme (SBE) were also decreased after W 72h in KR5. In 'Hayward', two out of three SBE encoded unigenes were found down-regulated. The Glc-1-P which is an important precursor for starch biosynthesis is synthesized via plastidial phosphoglucomutase (PGM) in plant roots. The W 72h treatment inhibited the expression of most PGM unigenes in KR5, while it promoted all PGM unigenes in 'Hayward'. In addition, most of the unigenes involved in starch cleavage (e.g., α-amylase and isoamylase) were down-regulated while 11 out of 19 unigenes encoded β-amylase were up-regulated after W 72h in KR5. On the other hand, one differentially expressed unigene encoding α-amylase and seven β-amylase encoded unigenes were down-regulated in 'Hayward'. While, one differentially expressed unigene encoding isoamylase was up-regulated after W 72h.
The trehalose biosynthesis pathway shares similarities to sucrose and starch synthesis. The accumulation of trehalose content was stimulated after W 72h ( Figure 8A). In KR5, 23 out of 27 trehalose-phosphate synthase (TPS) unigenes were up-regulated at the transcription level ( Figure  8B, Table S3). Among the other four unigenes, three unigenes annotated as TPS5 were all down-regulated. In 'Hayward', the only one significantly down-regulated TPS gene was also annotated as TPS5. Most of the unigenes encoding trehalose-6-phosphate phosphatase (TPP), which catalyzes the conversion of T6P to trehalose, showed promoted expression in both genotypes.

KR5 Showed a Stronger Capacity of Anaerobic Fermentation Than 'Hayward' After W 72h
Most unigenes encoding the three key regulatory enzymes of glycolysis, viz. [hexokinase (HXK), phosphofructokinase (PFK), and pyruvate kinase (PK)] were up-regulated after W 72h in KR5 (Table S4). Similarly, most of HXK encoded unigenes and one PK encoded unigene were also elevated in 'Hayward'. In contrast, no significant change of PFK encoded unigenes were found in 'Hayward'. The ATP-dependent PFKs catalyze the phosphorylation of D-fructose 6-phosphate. Additionally, PPi-PFKs (PFP) can render the action reversible by using inorganic phosphate (PPi). A big number of PFPs showed an increased level after W 72h in both genotypes. Glyceraldehyde phosphate dehydrogenase (GAPDH) catalyzes a phosphorylation step, by using ppi rather than ATP. Most of differentially expressed GAPDHs were up-regulated after W 72h in both KR5 and 'Hayward'.
The expression of ethanol fermentative genes [e.g., pyruvate decarboxylase (PDC), alcohol dehydrogenase (ADH), and lactate dehydrogenase (LDH)] were significantly promoted after W 72h in these two genotypes at the transcript level (Table S4). Compared with W 0h, the W 72h treatment significantly inhibited the PDC and LDH enzyme activity in roots of KR5, while promoted their activities in 'Hayward' roots ( Figure 9A). Anyway, enzyme activities of PDCs and LDHs in waterlogged KR5 roots were still higher than that in 'Hayward'. The ADH activity was enhanced by waterlogging treatment in both genotypes. A higher increase in ADH was induced in KR5 than in 'Hayward'.

Differential Regulation of Unigenes Associated With Nitrogen and Amino Acids Metabolism After W 72h
The nitrate transporters / channels belong to NPF, NRT, and ChLoride Channel (CLC) families [41]. In this study, expression levels of the NPFs and NRTs decreased after W 72h in roots of both KR5 and 'Hayward', While enhanced expression for most of CLCs was found only in KR5 roots (Table S5). Both nitrate reductase (NR) and nitrite reductase (NIR) are involved in the nitrate assimilation. The NRs and NIRs gene expression was shown to increase significantly after W 72h in KR5, while no significant change was found in 'Hayward'.
In the roots of both KR5 and 'Hayward', the expression of alanine aminotransferase (AlaAT) was induced after W 72h (Table S5). The synthesis of alanine by AlaAT is accompanied by the formation of 2-oxoglutarate (2-OG) at the expense of glutamate and pyruvate. In recycle, pyruvate can be regenerated by alanine-glyoxylate transaminase (AGT). The expression of AGTs was enhanced in KR5, while it was inhibited in 'Hayward'. The metabolite 2-OG can then be metabolized to the GABA shunt or TCA cycle. Among the reactions in the GABA shunt, ammonium is assimilated by glutamate dehydrogenase (GDH), and GABA production is then catalyzed by glutamate decarboxylase (GAD). Expression of GDHs and GADs were found mostly promoted by W 72h treatment in roots of KR5 and 'Hayward' (Table S5). Besides, aspartate and 2-OG can be converted to glutamate and oxaloacetate by aspartate aminotransferase (AspAT). The up-regulation of AspATs were observed after W 72h in both genotypes.

A Larger Increase in Transcript Level and Enzyme Activity of Manganese Superoxide Dismutases After W 72h in roots of KR5
Transcript levels of respiratory burst oxidase (RBOH), which catalyzes the production of superoxide anion, increased after W 72h in both genotypes. Superoxide dismutases (SODs) catalyze the dismutation of superoxide anion to hydrogen peroxide (H2O2). They are classified into copper-zinc (Cu/Zn-SOD), iron (Fe-SOD), or manganese SODs (Mn-SOD) according to their metal cofactors [42]. The expression of Mn-SOD (mitochondria) was found highly promoted while the other two types of SOD were generally repressed in KR5 at the transcript level (Table S6). The W 72h treatment also promoted the expression of Mn-SODs in 'Hayward', while didn't significantly influence the expression of the two types of SODs. At the protein level, the KR5 showed a higher activity of Mn-SOD than 'Hayward' under normoxia ( Figure 9B). Furthermore, the Mn-SOD activity was significantly elevated in KR5 roots, while a significant decline in activity occurred in 'Hayward' roots after W 72h.
The peroxidase (POD), catalase (CAT), and glutathione peroxidase (GPX) reduce H2O2 to water and oxygen. We found that 34 out of 60 and 16 out of 22 POD encoded unigenes were down-regulated after W 72h in KR5 and 'Hayward', respectively (Table S6). In KR5, all the four CAT encoded unigenes were significantly down-regulated while GPX encoded unigenes were increased. Expression of CATs and GPXs showed no significant change after W 72h in 'Hayward'. Ascorbate is one of the non-enzymatic antioxidants. Monodehydroascorbate reductase (MDAR) and dehydroascorbate reductase (DHAR) maintain the reduced states of antioxidants. Transcript levels of most MDAR and DHAR encoded unigenes increased after W 72h in KR5 and 'Hayward'.

Differential Regulation of Genes Encoding Transcription Factors
We analyzed the waterlogging-regulated TFs after W 72h. In total, 444 and 125 TFs from ten chosen protein families were identified from KR5 and 'Hayward', respectively ( Figure S4). About half of the total 444 TFs from KR5 were up-regulated while 64.0% of the total 125 TFs from 'Hayward' were down-regulated at the transcript level. Both the ERFs and bHLHs were pronounced in these two genotypes. In this study, the only ERF-VII unigene identified in KR5 was significantly promoted after W 72h. While, one out of two ERF-VII unigene (TRINITY_DN80164_c1_g2) in 'Hayward' was also elevated but the results were not significant at p ≤ 0.05. Besides, most of WRKYs were noticeably up-regulated in KR5 after W 72h.

Other Genes Predicted to Highly Expressed in Response to Waterlogging
As discussed above, W 72h was an important time point for transcript adjustments in these two genotypes. So, the top 1% up-regulated DEGs in the comparison W 0h VS W 72h were chosen for further analysis. The K-means results showed that nearly all of them were up-regulated as early as after W 12h ( Figure S5B, Figure S6B). According to the KEGG analysis results, some pathways were mutually enriched in both KR5 and 'Hayward', e.g. fatty acids degradation, tyrosine metabolism, glycolysis, and alpha-linolenic acid metabolism. Besides, stilbenoid, diarylheptanoid, gingerol biosynthesis and carbon fixation in photosynthetic organisms were specifically enriched in KR5 ( Figure S5A). Whilst, carotenoid biosynthesis and glutathione metabolism were only enriched in 'Hayward' (Figure S6A).
In Arabidopsis thaliana, about 51 core hypoxia genes (CRGs) were found up-regulated in response to submergence, regardless of organ or cell type [8]. The orthoFinder was applied to find the orthologues in kiwifruit genotypes of these core hypoxia response markers. Hence, we identified the clusters that were homologous to 7 CRGs among kiwifruit and Arabidopsis thaliana (Table S7). By comparing the expression level, genes annotated as alanine aminotransferase 2 and hypoxia response attenuator were found up-regulated in response to waterlogging in kiwifruit ( Figure S7).

Validation of DEGs By qRT-PCR
Each eight genes from KR5 and 'Hayward' were respectively selected and used for qRT-PCR analysis to verify the transcriptome data. These genes encode TFs or enzymes involved in different pathways. We found that the change trends identified by qRT-PCR and RNA-seq were similar ( Figure. S8) and, therefore, the results of RNA sequencing could be used to quantify DEG expression.

Transcriptional Adjustments Under Waterlogging Stress in Two Kiwifruit Genotypes
Waterlogging always limits the growth of kiwifruit vines [34]. In this study, the tolerant kiwifruit genotype KR5 was found more able to keep them from harm by waterlogging than the sensitive genotype 'Hayward' (Figure 1), as reported previously [35]. Plant adaptation to waterlogging relied on the adjustment of the global gene expression [23,43]. The RNA-sequencing result showed that KR5 shaped its transcription in response to waterlogging earlier and vigorously than 'Hayward' (Figure S2). Ethylene and NO signaling seemed to play an important role in kiwifruit hypoxia sensing (Figure 7, Table S1), like other plant species [11,44]. Induced expression of PCOs in KR5 indicated the involvement of ERF-VIIs in the process. The manipulation of kiwifruit ERF-VIIs could improve the waterlogging tolerance in tobacco [45]. For this reason, the role of kiwifruit ERF-VIIs in transcript regulation under waterlogging stress needs more studies.
Some waterlogging responsive genes expressed themselves differently according to the duration of different waterlogging treatments. These genes reflected the conserved transcript patterns against waterlogging in both genotypes. Among these common DEGs, 60.1% of them in KR5 and 95.7% of them in 'Hayward' showed up-regulation trend of transcript abundance during W 12h to W 72h (Figure5, Figure 6). The KR5 might take an advantage over 'Hayward' to adapt to waterlogging stress by adjusting a larger set of down-regulated genes. Usually, down-regulation of metabolism like protein and DNA synthesis could decrease the energy consumption [46]. We also paid attention to the top 1% DEGs up-regulated in the comparison of W 0h VS W 72h. Nearly all of these 1% DEGs were kept enhanced during W 12h to W 72h in both genotypes ( Figure S5, Figure S6). The W 12h was proved to be an important turning point for the transcript adjustments.
The results for GO analysis showed that "carbohydrate metabolism", "phosphorylation", "oxidoreductase activity", and "steroid metabolisms" were enriched in both two genotypes against waterlogging (Figure 3, Figure 4). This revealed the basic metabolism responses in kiwifruit upon waterlogging, similar to other plant species [26,28,47]. The glycolysis and fatty acid degradation pathways were mutually activated by these two kiwifruit genotypes. Thus, glycolysis helps the waterlogged vines to endure the energy crisis [23]. Acceleration of fatty acid degradation might decrease the fatty acid level [48], and produce more energy and metabolites. In KR5, numerous pathways involved in the metabolism of the amino acids, such as arginine, alanine, aspartate, and glutamate, were uniquely enriched ( Figure 5, Figure S5). Perhaps, a release of free amino acids contributes to maintain the osmotic potential of its roots under waterlogging stress [49]. In 'Hayward', pathways involved with carbohydrate metabolism and carotenoid biosynthesis were enriched in the up-regulated DEGs (Figure 6, Figure S6). The above results pointed to the mutual and unique transcript patterns of both genotypes under waterlogging stress.

Role of Hormones in the Acclimation of Waterlogging
Ethylene typically tend to accumulate in the tissues of flooded plants [50]. In this study, both ACSs and ACOs were promoted in response to waterlogging. Usually, ethylene acts in concert with other hormones [51]. Our results indicated that GA, ABA and BRs signaling possibly regulated the kiwifruit responses to waterlogging together with ethylene (Table S2). Down-regulation of repressors of GA signaling indicated the enhanced GA activity of these two genotypes. The catabolism of ABA is mainly regulated at the transcriptional level under stress [52]. Expression of key ABA biosynthetic genes i.e. NCEDs, followed a decreasing trend in KR5, whereas most of them were promoted in 'Hayward'. Transcript levels of BZR1, regulator of a host of BR-responsive genes, significantly declined in KR5 after W 72h while increased in 'Hayward'. Changes in transcript level of these hormone-related genes possibly brought about the different adaptive responses to waterlogging in these two kiwifruit genotypes.

Adaptation of Carbohydrate Metabolism, Glycolysis and Fermentation to Waterlogging
Carbohydrate metabolism plays an important role in orchestrating the plant response to O2 deficiency [53]. At the physiological level, W 72h stimulated a higher accumulation of sucrose and trehalose content in KR5 and a higher starch content in 'Hayward' (Figure 8, Table S3). The presence of higher soluble carbohydrate reserves in KR5 possibly helped them to survive longer during waterlogging. At the transcript level, the expression of biosynthetic genes of sucrose declined and transporters genes of sucrose and hexose were mostly increased in both genotypes. The up-regulated expression of CWINs in KR5 contributed more to an appropriate source/sink balance within the entire plant [54]. The activated expression of these transporter genes might drive the accumulation of sucrose in the KR5 roots. The sucrose cleavage via SUS is more energy-efficient than that of invertase [55]. The enhanced expression of SUSs indicated their important function during waterlogging in both genotypes. The trehalose 6-phosphate regulates carbon allocation and utilization in plants under diverse stresses [56]. A higher transcript level of TPSs linked to the higher content of trehalose in KR5. The inhibition of TPS5, opposite to other class II TPSs, possibly played a special role during waterlogging in these two genotypes [57]. Starch hydrolysis can produce large amounts of neutral products such as maltose and trehalose. The distinct transcript patterns of starch biosynthesis-related PGMs and starch hydrolysis-related BAMs might influence the content of hexose in both genotypes ( Figure 8, Table S3).
Under waterlogging stress, root cells altered their metabolism to increase glycolysis and the activity of fermentative enzymes [58]. A significant induction of PFP and GAPDH, two PPi-dependent enzymes, were found after W 72h at the transcript level in these two kiwifruit genotypes. Activation of enzymes that use pyrophosphate (PPi) as the energy source rather than ATP is beneficial for saving energy. Alcohol fermentation is necessary for sustaining glycolysis and the re-oxidization of NADH. Based on the enzyme activity, we can infer that KR5 had a stronger capacity of anaerobic fermentation than 'Hayward' (Figure 9). Especially, the higher ADH enzyme activity in KR5 could account for more regeneration of NAD + and the less acidosis of the cytoplasm [1].

Adaptation of Nitrogen Metabolism in Response to Waterlogging
NRs and NiRs expression increased significantly after W 72h in both genotypes, similar to the observations in other plant species [59,60]. Nitrate reduction, catalyzed by NR, plays a role in the acclimation to O2 deprivation by regenerating NAD + from NADH. The nitrite participates to produce the signal molecule NO and contributes to ATP synthesis during oxygen deprivation [20,61]. Further, NO emissions may activate ROS production and release Ca 2+ to regulate the downstream responses [31].
Modulation of alanine metabolism is important for the adaptation to oxygen deprivation [62]. In both KR5 and 'Hayward', the expression of AlaATs and GADs were mostly enhanced after W 72h. It seemed that changes in AlaAT expression and GABA shunt were associated with the alanine accumulation in kiwifruit [26]. Given that GAD catalyzes is a proton-consuming reaction, the up-regulation of GADs were beneficial for the acid amelioration under waterlogging [63]. Different expression patterns of AGTs in KR5 and 'Hayward' was directed towards their distinctness in alanine metabolisms. The up-regulation of the expression of GDHs and AspATs after W 72h pointed out the roles of glutamate and aspartate metabolisms in kiwifruit response to waterlogging.

Adaptation of Reactive Oxygen Species Scavenging in Response to Waterlogging
Reactive oxygen species (ROS) are overproduced as a consequence of hypoxia [35]. Significant increases in the mRNAs encoding enzymes involved in ROS generation were observed after W 72h in both genotypes, similar to other findings [60]. ROS are detoxified by various cellular enzymatic and nonenzymatic mechanisms [30,31]. The changed transcript levels of SOD, POD, CAT, GPX, MDAR and DHAR indicated the roles of them in ROS scavenging. The Mn-SOD has been found in many plant species [65,66], and the expression of mitochondrial Mn-SOD was induced by leave senescene, salt stress, and drought [67][68][69]. In this study, Mn-SODs were especially induced by waterlogging in both genotypes. Moreover, the enzyme activity of Mn-SODs in KR5 was higher than that in 'Hayward' after W 72h. The Mn-SODs possibly played a unique role in the ROS metabolisms in kiwifruit under waterlogging.

Plant Materials and Stress Conditions
Clonally propagated vines of 'KR5' from Actinidia valvata (a waterlogging-tolerant genotype) and of 'Hayward' from Actinidia deliciosa (a waterlogging-sensitive genotype) were used in this study. Vines were transplanted to 18-cm-diameter pots containing a mixture of peat moss, perlite, and sand (1:1:1 v/v) and were normally managed. At five to six leaf-stage, the vines were selected for waterlogging treatments. To perform waterlogging, three potted vines of one genotype were placed in a plastic container (45 cm × 35 cm × 16 cm) filled with tap water. The level of water was maintained at 2 cm above the soil surface. The vines were subjected to waterlogging for time periods of 0, 12, 24, and 72 h. The waterlogging time period of 0 h was used as the control.

RNA-Seq sampling, RNA Extraction, cDNA Library Construction, and Illumina Sequencing
Root samples were taken at 0, 12, 24, and 72 h after waterlogging treatment. Roots of nine vines per genotype were harvested at each time point. The excised roots were immediately placed in liquid nitrogen and then stored at -80 ℃. In total, 24 samples (two genotypes, four time points, three biological replications) were prepared for the RNA sequencing. Total RNA was extracted from kiwifruit roots using TRIzol Reagent (Invitrogen). An amount of one microgram for total RNA was used to construct the library using the NEBNext® Ultra™ RNA Library Prep Kit for Illumina®, according to the manufacturer's protocol. The library was constructed based on a 2×150-bp paired-end configuration, and sequencing was conducted using an Illumina HiSeq2000 instrument.

Quality Control, Assembly, and Annotation
Raw data were processed using Cutadapt (version 1.9.1) to remove technical sequences, and the trimmed reads were assembled using Trinity [70]. The duplicated contigs were removed using cd-hit to obtain the unigene sequence file. The BLAST software was used to annotate the unigene sequences. The annotation databases included GenBank's non-redundant (Nr) protein, eukaryotic ortholog groups (KOG), swissprot, kyoto encyclopedia of genes and genomes (KEGG), and gene ontology (GO).

Differential Gene Expression, K-Means Clustering, Go and Kegg Enrichment Analysis
Gene and isoform expression levels from the pair-end clean data were estimated using RSEM (v1.2.6) with the unigene sequence file being used as a reference. The DESeq Bioconductor package was used for differential expression analysis. After adjustments using the Benjamini and Hochberg approach for control of the false discovery rate (FDR), we set a P-value of ≤ 0.05 and the absolute value of log2Ratio ≥ 1 as the criterion for determining the differential expression of genes (DEGs).
Then the mean log2-expression ratios of DEGs were used for K-means clustering. The steps of the K-means algorithm are described in this study [71]. All the DEGs were grouped into clusters using the sum of squared error. We used GO-Term Finder to identify GO terms to which enriched genes with a significant P-value of less than 0.05 were annotated. The significant DEGs in KEGG pathways were enriched using in-house scripts.

qRT-PCR Analysis of Gene Expression
According to the manufacturer's protocol, cDNA was synthesized using a reverse Tra Ace qPCR RT kit (TOYOBO, Japan). The qRT-PCR reactions were performed with a LightCycler 480 II machine (Roche) using a LightCycler 480 SYBR Green I Master Mix in a 20-μL reaction volume. Melting curve analysis was conducted to validate the amplification specificity of each run. The expression was normalized against the Actin gene [72]. The sequences of primers used for amplification are listed in the table under the supplementary files section (Table S8).

Measurements of Root Activity, Relative Water Content, Carbohydrate Content and Enzyme Activities
Six vines per genotype at each time point were collected. The fresh roots were used for the root activity measurements by the triphenyl tetrazolium chloride (TTC) method [73]. The above-ground parts of vines were freshly weighed using an electronic banlance and then were oven-dried at 65℃ for 48 h. The relative water content (RWC) of the aboveground tissues was calculated by the following formula. RWC = (FW -DW) / FW ×100. The roots stored at -80 ℃ were used to measure the carbohydrate contents and enzyme activities. The sucrose contents were measured using the anthrone-sulphuric method [73]. The starch and trehalose contents were measured using their corresponding kits from Solarbio [74,75]. The activities of pyruvate decarboxylase (PDC), alcohol dehydrogenase (ADH), and lactate dehydrogenase (LDH) were assayed by the kits from Suzhou Comin Biotechnology Company. Activities of different types of superoxide dismutases were assessed by the corresponding kits from Nanjing Jiancheng Bioengineering Institute. Protein concentrations of extractions were determined by the bicinchoninic acid (BCA) method. The sampling weights, dilution ratios, operation steps, and calculation methods were all according to the manufacturer's instructions for each kit. All the tests were repeated at three biological replicates. The statistical analysis were conducted using SPSS Statistics, version 22.0. Significant differences by analysis of variance (ANOVA) were assayed by Duncan's multiple range test (P < 0.05).

Conclusions
The transcriptional adjustments of genes involved in hormone signaling transduction, carbohydrate metabolism, nitrogen metabolism, and reactive oxygen species scavenging occurred in kiwifruit roots upon waterlogging stress. The waterlogging-tolerant genotype KR5 shaped its transcript more vigorously than the sensitive genotype 'Hayward' under waterlogging. The differential regulation of genes encoding 9-cis-epoxycarotenoid dioxygenase, phosphoglucomutase, alanine-glyoxylate transaminases, and other enzymes might characterize the associated tolerance mechanisms of these two kiwifruit genotypes. The KR5 showed a stronger capacity of anaerobic fermentation than 'Hayward' along with enhanced accumulation of sucrose and trehalose in the root. The promoted transcript levels of manganese superoxide dismutases indicated their special roles in kiwifruit responses to waterlogging.

Supplementary Materials
Supplementary materials can be found at www.mdpi.com/xxx/s1.