RNA-Seq reveals novel genes and pathways associated with hypoxia duration and tolerance in tomato root

Due to climate change, economically important crop plants will encounter flooding periods causing hypoxic stress more frequently. This may lead to reduced yields and endanger food security. As roots are the first organ to be affected by hypoxia, the ability to sense and respond to hypoxic stress is crucial. At the molecular level, therefore, fine-tuning the regulation of gene expression in the root is essential for hypoxia tolerance. Using an RNA-Seq approach, we investigated transcriptome modulation in tomato roots of the cultivar ‘Moneymaker’, in response to short- (6 h) and long-term (48 h) hypoxia. Hypoxia duration appeared to have a significant impact on gene expression such that the roots of five weeks old tomato plants showed a distinct time-dependent transcriptome response. We observed expression changes in 267 and 1421 genes under short- and long-term hypoxia, respectively. Among these, 243 genes experienced changed expression at both time points. We identified tomato genes with a potential role in aerenchyma formation which facilitates oxygen transport and may act as an escape mechanism enabling hypoxia tolerance. Moreover, we identified differentially regulated genes related to carbon and amino acid metabolism and redox homeostasis. Of particular interest were the differentially regulated transcription factors, which act as master regulators of downstream target genes involved in responses to short and/or long-term hypoxia. Our data suggest a temporal metabolic and anatomic adjustment to hypoxia in tomato root which requires further investigation. We propose that the regulated genes identified in this study are good candidates for further studies regarding hypoxia tolerance in tomato or other crops.

Global climate change has been associated with frequent flooding periods over the past 60 years. Flooding has negative consequences for crop yields, with reductions of up to 10% or 40% in severe cases causing a decrease in food supply for an increasing human population 1,2 . Therefore, the development of flood tolerant cultivars is crucial for sustainable agriculture. Although some improvements in flood tolerance have been achieved through natural genetic variation, e.g. Submergence 1 (SUB1) gene in rice, flood tolerance in other crops requires further investigation 3 .
One component of flood stress is oxygen deprivation, or hypoxia, which is due to the ca. 10 4 fold lower diffusion rate of oxygen in water than in air 3 . Oxygen deficiency can also occur in tissues with high metabolic activity as well as bulky or dense tissues such as phloem, meristems, seeds and fruits 1 . Hypoxia is associated with disturbances in energy supply due to the inhibition of mitochondrial respiration, pH alteration, changes in redox state, reactive oxygen (ROS) accumulation and Ca +2 spiking 4,5 . A set of 49 hypoxia-induced genes encoding transcription factors (TF) as well as metabolic activities such as glycolysis and fermentation has been identified in Arabidopsis thaliana 6 . Moreover, it has been shown that hypoxia is important for general plant development such as shoot meristem activity and root architecture 7,8 .
Hypoxia sensing and signalling in Arabidopsis thaliana is conducted via the oxygen-and NO-dependent N-degron pathway which regulates the proteolysis or stability of subgroup VII ETHYLENE-RESPONSE FACTOR (ERFVII) TFs, such as HYPOXIA RESPONSIVE1 and 2 (HRE1 and HRE2) and RAP2.12. During the No significant differences were observed between short term hypoxia treatments and controls in fresh and dry weight differences of roots (Fig. 1a,b), water content (%) and SPAD values (Fig. 1c,d). As anticipated, these results indicate that 6 h hypoxia is not long enough to affect the above-mentioned parameters. However, hypoxia resulted in a significantly lower root fresh weight after 48 h (Fig. 1a). Nonetheless, no statistically significant differences  were observed in response to hypoxia (48 h) for root dry weight and root water content (%) (Fig. 1b,c). Moreover, after 48 h, relative chlorophyll levels, represented by SPAD values, showed a significant (P < 0.05) reduction in response to hypoxia (Fig. 1d). This indicated that long-term hypoxia at the root could have an impact on leaf performance. We observed that short-term hypoxia (6 h) resulted in regulation changes in only 267 genes (38 down-regulated and 229 up-regulated) while long-term hypoxia (48 h) resulted in regulation changes of 1421 genes (524 down-regulated and 897 up-regulated) (Fig. 2a). The fact that long-term hypoxia led to a fivefold higher number of regulated genes indicated that hypoxic stress duration has a strong impact on adaptive responses induced by gene expression.

Functional classification of differentially regulated genes using GO terms and MapMan categories.
To identify the most highly regulated hypoxia-responsive pathways in response to short-and long-term stress, a Gene Ontology (GO) term enrichment analysis (molecular function, biological process and cellular compartment) was performed on the up-and down-regulated genes. To visualize substantially regulated GO terms after 6 h and 48 h hypoxia, all significantly enriched (Padj < 0.05) GO terms describing molecular function for the upand down-regulated genes are displayed (Fig. 2b). The list of all the GO terms associated with the regulated genes is provided in Supplementary Table 1. It must be noted that GO terms refer to the proteins encoded by the genes and therefore in some cases, the word "activity" is used in GO term results.
GO term analysis revealed that down-regulated genes encoded proteins with peroxidase and oxidoreductase activity, at both time points. Among those GO terms, which were associated only with down-regulated genes Enriched GO terms (Padj < 0.05), describing molecular function, among down-and up-regulated genes in response to 6 h and 48 h hypoxia. The regulated genes in all samples were analysed for enriched GO terms using the online tool PANTHER 14.0 and Solanum lycopersicum as a reference organism. The light and dark green bars represent all significantly enriched GO terms associated with down-regulated genes in response to 6 h and 48 h hypoxia, respectively. The light and dark red bars represent all significantly enriched GO terms associated with up-regulated genes in 6 h and 48 h hypoxia samples, respectively.
at 48 h, we observed that cell wall, water channel and transmembrane transporter activity showed the highest enrichment. Moreover, the up-regulated genes encoded for proteins that were involved in cation binding, metal ion binding and vitamin binding at both time points. Thiamine pyrophosphate binding and sulfur compound binding showed the highest enrichment at 6 h only while thioredoxin activity and phenylalanine ammonia-lyase activity represented the maximum enrichment in up-regulated genes at 48 h. It was noteworthy that genes up-regulated at 48 h only and not 6 h, were enriched in oxidoreductase activity and transcription factor activity, although these did not represent the highest enrichment levels found.
These results indicate that the short-and long-term hypoxia had distinct and also conserved effects on cellular activities associated with different stress response pathways such as cell wall, redox homeostasis and water transport. Moreover, regulation of transcription appears to be more active at later time points (>6 h hypoxia) resulting in modulation of hypoxia responsive genes.
MapMan categories 29 based on ITAG2.3 annotations, were used to provide independent verification of results and to gain a more in-depth view of the biological pathways related to differentially regulated genes (Supplementary Table 2). We compared the result of both GO terms (Supplementary Table 1) and MapMan categories (Supplementary Table 2) related to regulated genes in our study. We observed similarities between GO terms and MapMan categories related to cell wall, antioxidant activity, oxidoreductase activity, organic acid metabolic process, glycolytic process (glycolysis) and reactive oxygen species metabolic processes. However, only MapMan categories were used for the description of the regulated genes in our study.
Validation of differentially expressed genes in response to hypoxia using qPCR. Validation of RNA-Seq data was performed using qPCR on 18 hypoxia responsive genes at both 6 h and 48 h, excluding two genes, which were up-regulated only at 48 h (PGB1 and PGB3) (Supplementary Table 4). We observed a positive correlation between fold-change values (hypoxia/control) of RNA-Seq and qPCR data for both 6 h (r 2 = 0.90) and 48 h (r 2 = 0.90) time points (Fig. 3a,b). These data confirm the reliability of the RNA-Seq results used in this study. The expression fold-changes of the 18 selected genes for RNA-Seq and qPCR are provided in Supplementary  Table 3.
Short-and long-term hypoxia regulates expression of genes involved in carbon (C) flux and amino acid metabolism. We observed that three tomato genes annotated as SUS4 showed regulation changes (one down-regulated and two up-regulated) at 6 and 48 h hypoxia (Fig. 4a). Moreover, several genes coding for enzymes, which are involved in glycolysis (ENO2, PPC4 and TPI) and fermentation (PDC2, ADH1), were up-regulated at both time points in response to hypoxic stress (Fig. 4a). These data indicate that induction of the genes encoding enzymes involved in alternative routes for ATP production is one of the earliest responses under hypoxia.
Alanine (Ala) and 2-oxoglutarate (2OG) as well as gamma-aminobutyric acid (GABA) shunts are amino acid related pathways leading to accumulation of Ala, GABA and succinate upon hypoxia 30 . Among different regulated genes involved in amino acid synthesis and degradation (Fig. 4b), two transcripts encoding AlaAT showed up-regulation in response to hypoxia. However, one transcript (solyc0g123610.2.1) was up-regulated only at 6 h hypoxia while the other one (solyc0g123600.2.1) showed up-regulation at both time points. Furthermore, a transcript encoding glutamate decarboxylase 4 (GAD4) was up-regulated at 6 and 48 h hypoxia in our study. Further investigations are required to confirm the role of Ala-and GABA shunts in the acclimation of tomato roots to hypoxia.

Hypoxia affects the link between N metabolism and NO formation/scavenging. NO formation
under normoxic conditions is very limited due to the high ratio of NO 3 − to NO 2 − (50 to 100 fold) and NR preference for NO 3 − rather than NO 2 − . However, hypoxia results in a higher accumulation of NO 2 − and therefore a higher level of NO can be produced via the NR pathway and mitochondrial cytochrome c oxidase [31][32][33] .
We investigated the effect of hypoxia on N metabolism related processes such as induction of genes encoding NO 3 − transporters, nitrate and nitrite reductase as well as nitric oxide (NO) scavenging via phytoglobin (PGB). www.nature.com/scientificreports www.nature.com/scientificreports/ We observed that a transcript encoding NO 2 − reductase 1 (NIR1) was down-regulated after 6 h hypoxia, leading to accumulation of NO 2 − (Fig. 5). It has been suggested that NO and its interaction with PGB might play a role in hypoxic stress adaption 34 . We observed that two PGB encoding transcripts, PGB1 and PGB3, showed up-regulation (3.5 fold and 5 fold, respectively) after 48 h hypoxia (Fig. 5).
Seven high-and low affinity NO 3 − transporter encoding (NRTs) transcripts showed expression changes, mostly after 48 hypoxia. The only exception was NRT1.5, which was down-regulated at 6 h and 48 h post-hypoxia. Out of three transcripts annotated as NRT2.4, two were up-regulated. The other NRT2.4 as well as NRT1.1, NRT1.2 and an uncharacterized NO 3 − transporter gene were down-regulated after 48 h hypoxia (Fig. 5). Altogether, these results suggest that N metabolism might be involved in NO production and scavenging in tomato root adaption to hypoxic stress with a gene regulation specific response to the duration of hypoxia (Figs. 5 and 6). However, a precise measurement of NO level and N metabolites, in particular NO 2 − under hypoxia is necessary to confirm their role in hypoxia responses in tomato roots.
Redox related gene expression changes appear to be responsive to long-term hypoxic stress. ROS production is associated with low oxygen signalling and adaption in plants 35 . In the current study, genes encoding members of several redox-related enzyme families showed time-dependent regulation changes in response to hypoxia: catalase (CAT, one gene), peroxidase (PRX, 33 genes), monodehydroascorbate reductase (MDAR, 4 genes), glutathione reductase (GRX, 4 genes), thioredoxin (TRX, 5 genes). Eight genes encoding members of PRX, MDAR and TRX, were regulated at 6 h (3 up-regulated and 5 down-regulated). It was noteworthy that 47 redox-related genes were regulated after 48 h hypoxia with 22 and 25 genes being down-and up-regulated, respectively (Fig. 7). Moreover, it was observed that out of 5 regulated TRX transcripts, 4 showed up-regulation (Fig. 7).
Overall, our results indicate that there was a fine-tuning of redox homeostasis, which was dependent upon the duration of hypoxic stress. transcription factors as key regulators of the response to short-and long-term hypoxic stress. Ethylene response factor (ERF) TFs, in particular ERFVII, belonging to the APETALA2/ethylene response factor (AP2/EREBP), are among the most studied TFs known to be involved in hypoxic response. However, there is evidence that the huge transcriptome modulation observed under low oxygen is not limited to ERF TF targets 36 . Therefore, we investigated gene expression changes of members of different transcription factor families in response to short-and long-term hypoxic stress. A total of 122 tomato transcripts annotated as TF  Table 3).
These data indicate a fine-tune regulation of TFs which was sensitive to the duration of hypoxic stress. The hypoxic responsive TF encoding genes are particularly interesting candidates for hypoxia tolerance studies due to their crucial role as master regulators of stress response.
Hypoxia resulted in regulation of the genes linked to ethylene biosynthesis, ROS production and pcD. Aerenchyma formation through programmed cell death (PCD) in root cortical cells is a coping strategy under hypoxic conditions. PCD is the consequence of ethylene accumulation and up-regulation of respiratory burst oxidase (RBOH) genes involved in ROS production 37 .
We investigated the expression changes of rice orthologue genes in tomato root with potential roles in aerenchyma formation.
Our data analysis showed that genes encoding ethylene biosynthesis enzymes underwent upregulation in response to 48 h hypoxia e.g. ACS9 (>13 fold) and ACO1 (ca. 45 fold). We observed that out of eight members of the RBOH gene family in tomato (RBOHA to RBOHH), RBOHB was strongly up-regulated (ca. 18 and 53 fold after 6 h and 48 h hypoxia, respectively) (Fig. 9). A methallothionein encoding gene (MT2B, encoding metallothionein-like protein 2B), involved in ROS scavenging, was down-regulated after 48 h of hypoxia (Fig. 9). MT2B down-regulation might lead to ROS accumulation and subsequent aerenchyma formation, as it has been reported in rice 38 .
Among diverse regulated protein kinase encoding genes, we observed up-regulation of several CPKs transcripts: CPK9, CPK17 and CPK21. This up-regulation was observed at 6 h as well as 48 h after hypoxia (Fig. 9).
Among regulated transcripts encoding cell wall proteins, we observed 4 EXP encoding genes (EXPA18, EXPA3, EXPA6 and EXPB1) to be down-regulated after 48 h hypoxia. In total, eight XTH encoding genes showed expression changes in the current study with six (XTH24, XTH31, XTH32, XTH5, XTH8 and XTH9) being down-regulated and two (XTH16 and XTH22) being up-regulated (Fig. 9). Among cell wall related up-regulated www.nature.com/scientificreports www.nature.com/scientificreports/ genes, an expansin-like encoding gene, EXLB1, showed the highest up-regulation at 6 h (>24 fold) with a slight decrease (ca. 12 fold) at 48 h. Moreover, the cellulose-synthase-like CSLA09, was another highly up-regulated gene (>20 fold) at 48 h. Discussion transcriptome reprogramming in response to short-and long-term hypoxic stress. The root is the first organ to experience hypoxic stress under flood conditions. Hypoxia, due to lowered oxygen transport under flood stress, is particularly damaging to the root because of the lack of photosynthetic oxygen production 39 . The tomato is a widely used crop plant. However, little is known about the underlying molecular mechanism associated with hypoxia tolerance in tomato roots.
The current study utilises comparative transcriptomics to elucidate the differential gene expression taking place in tomato (cv. Money maker) roots in response to short-(6 h) and long-term (48 h) hypoxia. Long-term hypoxia resulted in significantly (P<0.05) lower fresh weight and chlorophyll content (SPAD values), (Fig. 1a,d) while short-term hypoxia did not show any significant effect on the above-mentioned parameters.
Overall, we observed a ~50 times higher number of regulated genes at 48 h compared to 6 h hypoxia (Fig. 2a). This substantial discrepancy indicates that a major shift occurs in transcriptome response after2 days of hypoxia, www.nature.com/scientificreports www.nature.com/scientificreports/ highlighting the importance of hypoxic stress duration. It was also apparent that the number of up-regulated genes in response to hypoxia was higher than down-regulated genes at both time points.
To validate the RNA-Seq data, we confirmed the expression changes in 18 differentially regulated genes using qPCR. We observed a high level of consistency between the qPCR and RNA-Seq expression data (Fig. 3a,b) showing a strong positive correlation at both time points (r 2 = 0.90), indicating the reliability of our transcriptome data. C flux and amino acid metabolism reconfiguration in response to hypoxic stress. The regulation of gene expression involved in primary metabolism and energy homeostasis is crucial to avoid the detrimental effects of mitochondrial oxidative phosphorylation inhibition under hypoxic stress 40 . The core hypoxia responsive genes related to primary carbon metabolism reconfiguration such as alcohol dehydrogenase (ALDH), pyruvate decarboxylase (PDC2) and sucrose synthase (SUS4) showed up-regulation in our study (Fig. 4a). Induction of the above-mentioned gene expression in different plant species under hypoxic stress has been reported. Additionally, it has been reported that ADH, PDC as well as SUS mutant plants exhibit reduced hypoxia tolerance 6,40,41 .
Aside from carbon metabolism related gene expression changes, we observed that genes encoding enzymes with pivotal roles in N metabolism showed up-regulation in response to hypoxia, e.g. alanine aminotransferase (AlaAT1, two transcripts), a key enzyme of amino acid metabolism and GAD4, encoding another key enzyme involved in GABA shunt, (Figs. 4b and 10). GAD4 and one of the two transcripts annotated as AlaAT were up-regulated at both time points while the other AlaAT encoding transcript showed up-regulation only at 48 h.
The AlaAT enzyme reversibly transfers an amino group from glutamate to pyruvate leading to the formation of 2-oxoglutarate and Ala 40,42 . Ala accumulation and up-regulation of the AlaAT gene has been observed in different plant species under hypoxic stress 43,44 . However, despite proline and betaine, Ala may not function as a protectant. It has been proposed that under hypoxia, the AlaAT/NADH-GOGAT cycle uses glycolysis produced pyruvate to store carbon in the form of Ala synthesized by AlaAT. In this process, glutamate provides the amino group while 2-oxoglutarate is consumed by NADH-GOGAT leading to glutamate and NAD + regeneration. Ala accumulation under hypoxic conditions will be beneficial for post-hypoxic conditions by providing pyruvate and glutamate via the AlaAT reverse reaction 40 .
The finding of up-regulation of glutamate decarboxylase (GAD4) in our study, one of the key enzymes in the GABA shunt, is in line with former studies in Arabidopsis thaliana which showed among the five homologous GAD, only GAD4 is induced upon hypoxia 43 . GABA production and its subsequent transamination to succinic semialdehyde might be beneficial for cytosolic pH stabilization as well as providing a route for conversion of pyruvate to Ala 3 . However, whether Ala and GABA accumulate under hypoxia and if they are involved in tomato root responses to hypoxia needs further investigation.
Hypoxia-induced RoS production, redox associated gene expression and ethylene related genes. Oxidative stress response and redox signalling, including production of ROS and NO as well as other reactive nitrogen species (RNS), is associated with low oxygen stress 45,46 . Production of ROS signalling is associated with the enhanced activity of PM bound NADPH oxidases (RBOH gene), through Ca +2 signalling and phosphorylation, in response to stress perception 47,48 .
Apart from GRX and TRX related genes which showed more up-regulation, down-regulation was more frequent among the members of antioxidant related families such as CAT, PRX and MDAR at both time points in response to hypoxia, i.e.3 out of 4 and 4 out of 5 regulated genes, respectively (Fig. 7b).
TRXs have been shown to be involved in several processes associated with oxidative proteins such as DNA damage related proteins, inducing activity in antioxidant protecting enzymes as well as in the regulation of scavenging mechanisms or signalling pathways in the antioxidant network. Moreover, TRXs have been shown to be involved in cold tolerance in rice 49 . However, the precise role of TRX members in response to hypoxia needs further investigation.
The expression of antioxidant genes (SOD1, AOX1A, APX and MnSOD) has been shown to play a role in aerenchyma development in wheat. It was shown that after 2 h hypoxia, antioxidant gene expression was up-regulated but down-regulated after 24 h 50 . www.nature.com/scientificreports www.nature.com/scientificreports/ Among cell wall related up-regulated genes, we observed several members of the EXP and XTH families (Fig. 9). Different classes of cell wall proteins such as xyloglucan endotrans glucosylase (XET) and expansin (EXP) are involved in cell wall loosening at lower pHs, e.g. upon hypoxia during waterlogging. However, many cell wall related genes, among them XTH and EXP encoding genes, showed down-regulation, particularly at 48 h hypoxia (Fig. 9). This indicates that hypoxia did not strongly induce genes involved in cell wall degradation in our study.
NO production during hypoxia increases through different pathways including coordinated activity of the root-specific plasma membrane-bound nitrate reductase (PM-NR) and nitrite:NO reductase (NI-NOR) 51 , nitrate reductase (NR) and the mitochondrial electron transport chain 31,52 . The induction of NO after hypoxia in the root occurs rapidly, within minutes, and therefore it has been suggested to be involved in hypoxic stress response in the root 53 . Some of the functions, which have been attributed to hypoxia-induced NO, include a role in the production of ATP 54 , interaction with alternative oxidase (AOX) 55 We observed up-regulation of two PGB encoding transcripts, PGB1 and PGB3, after 48 h hypoxia (Figs. 5 and 10). Up-regulation of both PGB encoding genes only after 48 h hypoxia might indicate that PGB expression in tomato roots is more a consequence of ATP level rather than as a direct response to O 2 level, as has been suggested in former studies using respiratory chain-and ATP production inhibitors 58 . Non-symbiotic hemoglobins are induced in root exposed to hypoxic stress and are involved in NO scavenging leading to the formation of nitrosyl-phytoglobins after interaction with NO 58,59 . Moreover, it has been shown that overexpression of PGB1 in A. thaliana enhances hypoxia tolerance and post-hypoxia survival 60 . www.nature.com/scientificreports www.nature.com/scientificreports/ Apart from cell death-associated aerenchyma formation in cortex cells, the emergence of adventitious roots in rice is associated with epidermal cell death induced by ethylene as well as H 2 O 2 and O 2 •− which act as signalling molecules in this process 61 . In contrast with findings in Arabidopsis and rice, ethylene is not involved in early flood recognition and adaption to hypoxia in the tomato cultivar Moneymaker 14 . Concordantly, we observed that only long term hypoxia led to the induction of genes encoding ethylene biosynthesis enzymes ACS9 and ACO1 (Fig. 9). However, further studies are required to confirm ethylene accumulation in tomato roots under long term hypoxia. It has been reported that in rice roots, induction of ACS1 and ACO5 genes upon hypoxia results in ethylene biosynthesis, which in turn increases expression of RBOHH. Moreover, Ca 2+ release from apoplast to cytosol results in RBOHH protein activation via protein kinases (CDPKs) such as CDPK5 and/or CDPK13.

Regulated transcription factors serve as potential targets for improved hypoxia tolerance.
Several TF families such as AP2/EREBP, MYB, WRKY, NAC and bZIP have been shown to be associated with different abiotic stresses 64 . Members of these families showed regulation changes in response to hypoxia in our study (Fig. 8 and 1 (MYB84, 48 h, (>15 fold)), solyc06g054620.2.1 (ZFN1, 6 h (>15 fold) and 48 h (>15 fold)). We suggest that these are good candidate genes for further investigations regarding hypoxia tolerance in tomato and other crops.
In the following section, we discuss the most promising TF candidates among the top 10 highly up-regulated TFs in our study: NAC060, RAP2.6 L and MYB84.
NAC060, has been shown to be induced by the sugar-ABA signalling pathway. However, NAC060 availability in the nucleus hinders sugar-ABA signalling leading to seedlings' insensitivity to sugars 65 . A former study in rice showed that the link between low oxygen signals and the sugar-sensing cascade (via SnRK1) is important for regulation of sugar and energy production and therefore rice growth under flood conditions 66 . These data indicate that NAC060 is a promising candidate for further investigation with respect to its role in sugar sensing and hypoxia tolerance in tomato.
Induction of RAP2.6 L in response to waterlogging and ABA has already been reported in Arabidopsis 67 . Overexpression of RAP2.6 L resulted in a reduction of water loss as well as membrane leakage. Moreover, overexpression lines showed higher expression levels of genes encoding several antioxidant enzymes, the ABA biosynthesis gene and one of the hypoxia marker genes, ADH1. The authors concluded that RAP2.6 L overexpression can delay premature senescence under waterlogging stress mainly through stomatal closure. Moreover, RAP2.6 L function might be linked to the ABI1-mediated ABA signalling pathway. Further studies are required to confirm the role of RAP2.6 L in tomato roots under hypoxia 67 .
Different stresses such as drought and salt stress as well as H 2 O 2 and ABA have been shown to induce MYB84 expression in soybean (Glycine max). MYB84 has been shown to be involved in drought tolerance in soybean and its overexpression lines showed a higher resistance to drought including enhanced survival rate, longer root length, elevated proline and ROS level, higher activity of antioxidant enzymes and rehydration rate as well as reduced malondialdehyde (MDA) content. Moreover, GmMYB84 was shown to bind directly to the GmRBOHB-1 and GmRBOHB-2 promoters, suggesting a possible link between GmMYB84, H 2 O 2 and antioxidant enzymes in controlling root growth in response to drought stress as well as under control conditions 68 . Feng et al., (2004), reported that MYB84 is the closest homologue of MYB68 and shows an overlapping expression pattern in root pericycle cells which suggests a redundant role in root development 69 . The role of MYB84 in hypoxia resistance has not been previously studied. However, since RBOHB and MYB84 showed up-regulation in our study, a similar model could be proposed for tomato roots in response to hypoxic stress. However, more studies are required to confirm the link between MYB84 and RBOHB in response to hypoxia. MYB84 appears to be a potential target for improving multi-stress resistance in plants.
Group VII members of the AP2/ERF (apetala2/ethylene response factor) transcriptions of Arabidopsis thaliana and rice have been shown to be involved in low oxygen sensing via the N-degron pathway 41,70 , such that they are stabilized under low oxygen conditions and induce the expression of hypoxia core genes. RAP2.12, RAP2.2 and RAP2.3 have been shown to accumulate in response to hypoxia; however, their accumulation has www.nature.com/scientificreports www.nature.com/scientificreports/ also been observed under normoxia 71,72 . On the other hand, two other ERF-VIIs members, HRE1 and HRE2, are hypoxia inducible 73 . RAP2.2 and RAP2.3 are RAP2.12 homologues and have a redundant function in response to various stresses 72 . RAP2.2, a member of the same subfamily as the rice (Oryza sativa) submergence tolerance gene, SUB1A, has been shown to be involved in ethylene signal transduction as well as hypoxic stress. Overexpression of RAP2.2 in Arabidopsis thaliana resulted in better hypoxia tolerance while its knockout lines showed an impaired resistance to hypoxia compared to wild type. Therefore, it has been concluded that RAP2.2 plays an important role in plant hypoxia resistance 71 . In the current study, RAP2.2 showed>7 and>8 fold up-regulation after 6 h and 48 h hypoxic stress, respectively. Up-regulation of RAP2.2 has also been shown to be involved in activation of ethylene biosynthesis genes such as ACS7 and ACO1 which resulted in ethylene production and regulation of genes responsive to waterlogging in Arabidopsis 70,71 . In accordance with this, we observed up-regulation of ACS9 and ACO1 after 48 h hypoxia. However, ethylene accumulation and its role during hypoxia in tomato root need further investigation.
We noted that, the RAP2.3 gene showed down-regulation (>3 fold) at both time points in our study. Moreover, a member of the AP2/EREBP family, Hypoxia Response Attenuator 1 (HRA1), was the third highest up-regulated TF gene, at both 6 h (>18 fold) and 48 h (>30 fold) hypoxia. HRA1 has been reported to be highly up-regulated in response to hypoxia and to restrain the induction of core hypoxia-responsive genes in Arabidopsis thaliana. Since the constitutive accumulation of RAP2.12 hinders a submergence and hypoxia tolerance response, it is supposed that a constant induction of stress response genes has a negative effect on stress response. HRA1 has a negative regulatory role (negative feedback mechanism) on RAP2.12 level and activity resulting in higher flexibility of hypoxia response under fluctuating O 2 concentration 74,75 .
We compared differentially regulated genes (FDR <0.05) in response to 6 h hypoxia in our study with data presented in a very recent study based on the roots of submerged tomato seedlings 18 . We observed that ca. 12% (55 out of 450) of differentially regulated genes (FDR <0.05) in response to 6 h hypoxia in our study were also differentially regulated in the roots of submerged tomato seedlings (Supplementary Table 5). The relatively low overlap between regulated genes in these two studies is possibly due to: 1) the differences in the root samples which were analysed (whole root in this study vs. 1 cm of the root tip). Since root tips are exposed to more hypoxia developmentally, it is possible that their response to transient hypoxia differs from the upper parts of the root. 2) Plant age and developmental phase can affect the root response to hypoxia. We used 5 week old plants in our study while the above-mentioned study was conducted at the seedling stage. 3) Technical differences in hypoxia application can partially be responsible for the observed discrepancy. In the current study, a hydroponic system was used and therefore only roots were exposed to low oxygen while in the other study the whole seedlings were submerged.
Our data suggest that hypoxia tolerance in tomato roots is dependent on a precise modulation of transcription via members of various TF families, apart from members of AP2/EREBP, to maintain metabolic balance within the cell. Highly regulated TF genes have the potential to be considered as early and/or late biomarkers of hypoxia response in tomato root and possibly other crop plants. We propose a model representing cellular processes potentially involved in hypoxia tolerance with an emphasize on the genes identified in this study (Fig. 10).

conclusion
We found evidence that the tomato root is sensitive to the duration of hypoxia which is reflected in the strong difference between the number and the function of the regulated genes after 6 h and 48 h hypoxia. Short-term hypoxia (6 h) resulted in molecular changes which in part were transient, implying on acclimation to the stress. However, long-term hypoxia (48 h) led to a stronger reprogramming of the transcriptome which, to a large extent, was not observed under short-term hypoxia. The long-term hypoxia tolerance mechanism of the tomato root seems to be an escape mechanism to avoid low oxygen concentration through aerenchyma and adventitious root formation. This is associated with down-regulation of antioxidant related genes, such as CAT and MT2B, and up-regulation of ethylene biosynthesis genes (ACS9 and ACO1) as well as up-regulation of RBOHB which is involved in ROS production.
We identified transcription factor genes belonging to different families that were specifically regulated under short-or long-term hypoxia and some which were regulated during both time points. We hypothesize that some of these genes may have specific functions in regulating the target genes involved in hypoxia tolerance and may contribute to the differences between sensitive and tolerant varieties. Further studies are required to address this question and investigate the target genes for improving the hypoxia tolerance and selection of tomato cultivars with stronger adaptive responses to stress conditions. In summary, our data suggest a transcriptional acclimation to short term hypoxia. However, hypoxia progression results in a transcriptional reprogramming to support an escape mechanism probably through aerenchyma and adventitious root formation. This indicates the ability of a cultivated crop such as tomato to temporally adjust its response mechanism to hypoxia, both metabolically and anatomically. However, further investigations are required to confirm the precise mechanism of the escape strategy in tomato root.
Methods plant material and growth conditions. Tomato plants (Solanum lycopersicum L. cv. Moneymaker) were cultivated on sand in the greenhouse at 500 μmol photons/m2/s and 25 °C under a 14/10-h light/dark regime. All plants were treated with modified Hoagland nutrient solution containing 5 mM nitrate (NO 3 − ) as described previously 76 . Three week old plants were transferred to hydroponic conditions such that the roots were submerged in plastic pots containing ~6 L of nutrient solution and aerated by mild bubbling using aquarist air pumps (Hailea ACO-9620, Raoping, Guangdong, China) and air outlets (Tetratec, Osnabrück, Germany). The hypoxia treatment was conducted on five weeks old plant roots using N 2 gas (≥ 99.99 Vol. %) (Air Liquide, Germany). We harvested the roots after 6 h and 48 h of root exposure to N 2 gas. Oxygen concentration was measured using CellOx 325 DO electrodes (WTW, Germany) following manufacturer's instructions. For each root sample (after 6 h or 48 h hypoxia) three biological replicates, each containing a pool of three technical replicates, were harvested and immediately frozen in liquid nitrogen and stored in −80 °C (Fig. 10). We harvested the roots 6 h and 48 h after the initiation of root exposure to N2 gas. Hypoxia treatment was initiated at 8:00 a.m. The 6 h hypoxia and control samples were collected at 2:00 p.m. 48 h hypoxia and control samples were collected two days later at 8.00 a.m. The hypoxia and control samples at each time point were harvested and compared separately to remove the impact of other unknown or known additive factors between time points such as circadian rhythm 1 .
Relative chlorophyll levels of leaf #3 (the third leaf above the cotyledon) were determined using a Konica Minolta SPAD-502 chlorophyll meter. For each sample, SPAD values from three similar positions on the leaf #3 were recorded and their averages were calculated.
RnA isolation and cDnA synthesis. For RNA-Sequencing (RNA-Seq) and qRT-PCR analysis, total RNA was extracted from 250 mg frozen root tissue using a phenol-chloroform extraction method 77 . RNA concentration was quantified photometrically using a NanoDrop (ND-1000, Thermo Scientific, Wilmington, DE, USA) and RNA integrity was tested on 1.2% agarose gel. 2 μg DNaseI-digested total RNA was used for cDNA synthesis using oligo-(dT) 18  Read mapping and identification of differentially expressed genes. Adaptor clipped reads obtained from the NextSeq. 500 Illumina platform (LGC Biosearch Technologies, Berlin, Germany), were processed to omit short fragments (final length <20 bases) and low quality reads. The filtering of rRNA sequences was conducted using RiboPicker 0.4.3. The remaining reads were mapped onto the tomato reference genome (ITAG 2.4) (The Tomato Genome Consortium, 2012) using CLC Genomics Workbench (Qiagen, V. 7.5.5). Quantification of the gene expression levels were conducted using the CLC Genomics Workbench (Qiagen, V. 7.5.5). Sequencing data are deposited in the Sequence Read Archive (SRA) database (bioproject accession PRJNA553994) at the National Centre for Biotechnology Information (NCBI). The bioproject's metadata are available at https://dataview.ncbi.nlm.nih.gov/object/PRJNA553994?reviewer=684sto9a948tin240f0tlt1o1h.
For normalization and estimation of P-values the TMM (trimmed means of M values 78 ) and the edger algorithm 79 were used respectively. The log fold change values are provided as calculated by the edgeR algorithm. The P-values were adjusted for multiple testing 80 . All calculations were performed with the CLC Genomics Workbench software (Qiagen, V. 7.5.5). Differentially expressed genes (DEGs) with ≥2-fold expression change and P adj <0.05 (Supplementary Table 6) were selected for subsequent analysis. The FDR threshold was used for the P-value in multiple tests (P adj ). GO term enrichment analysis was performed using the Panther database 81,82 to identify significantly enriched GO terms (P adj <0.05) associated with DEGs in response to 6 h and 48 h hypoxia. For biological pathway analysis of differentially regulated genes, MapMan categories based on ITAG 2.3 annotations 29 were used. qpcR primer design and assay. Primers for quantitative real-time PCR (qPCR) were designed using QuantPrime 83 (Supplementary Table 7). qPCR reactions were performed in 5 μl total reaction volumes including 2.5 μl Power SYBR Green Master Mix (ThermoFisher Scientific), 0.5 μM forward and reverse primers and 0.5 μl cDNA. ACTIN was used as reference gene 84

Data availability
All materials and data sets represented in the current study are available in the main text or the supplementary materials. RNA-Seq data are deposited in the Sequence Read Archive (SRA) database (bioproject accession PRJNA553994) at the National Centre for Biotechnology Information (NCBI). The bioproject's metadata are available at https://dataview.ncbi.nlm.nih.gov/object/PRJNA553994?reviewer=684sto9a948tin240f0tlt1o1h.