Transcriptome Profiling of Watermelon Root in Response to Short-Term Osmotic Stress

Osmotic stress adversely affects the growth, fruit quality and yield of watermelon (Citrullus lanatus (Thunb.) Matsum. & Nakai). Increasing the tolerance of watermelon to osmotic stress caused by factors such as high salt and water deficit is an effective way to improve crop survival in osmotic stress environments. Roots are important organs in water absorption and are involved in the initial response to osmosis stress; however, few studies have examined the underlying mechanism of tolerance to osmotic stress in watermelon roots. For better understanding of this mechanism, the inbred watermelon accession M08, which exhibits relatively high tolerance to water deficits, was treated with 20% polyethylene glycol (PEG) 6000. The root samples were harvested at 6 h after PEG treatment and untreated samples were used as controls. Transcriptome analyses were carried out by Illumina RNA sequencing. A total of 5246 differentially expressed genes were identified. Gene ontology enrichment and biochemical pathway analyses of these 5246 genes showed that short-term osmotic stress affected osmotic adjustment, signal transduction, hormone responses, cell division, cell cycle and ribosome, and M08 may repress root growth to adapt osmotic stress. The results of this study describe the watermelon root transcriptome under osmotic stress and propose new insight into watermelon root responses to osmotic stress at the transcriptome level. Accordingly, these results allow us to better understand the molecular mechanisms of watermelon in response to drought stress and will facilitate watermelon breeding projects to improve drought tolerance.


Introduction
Water deficit is one of the main factors that affects agricultural productivity and limits the distribution of species worldwide [1]. Physiological and biochemical processes, such as transpiration, photosynthesis, respiration, carbohydrates and hormones, can be affected mostly by water deficits [2,3]. Plants have evolved a variety of strategies to deal with water deficits, including drought escape, drought avoidance and drought tolerance [4]. These strategies involve rich and complex metabolic networks, interactions and crosstalk between diverse molecular pathways.
Watermelon (Citrullus lanatus (Thunb.) Matsum. & Nakai) is an important cucurbit crop planted widely throughout the world [5]. Edible watermelon fruits provide essential nutritional compounds like sugars, lycopene and cardiovascular health-promoting amino acids [6,7]. However, insufficient water reduces watermelon fruit yield [8]. Few studies have examined the molecular mechanisms underlying watermelon responses to drought conditions, and these studies have focused on wild watermelon [4,[9][10][11]. Furthermore, root tissues in direct contact with dry soil are severely affected by drought, but there is little information about watermelon root responses to water deficits at the molecular level. To date, only two studies have examined drought stress in wild watermelons or closely related species at the molecular level. First, Si et al. used cDNA amplified fragment length polymorphism (cDNA-AFLP) to investigate the expression of diverse genes in Citrullus colocynthis (closely related to watermelon) roots under 20% PEG8000 induced drought stress [10]. In another report, protein profiling technology was utilized to identify drought-responsive proteins in wild watermelon root tissues [4]. In contrast with wild watermelon, domestic watermelon with relatively high resistance to water deficits has received less attention.
The low cost of high-throughput technology and the whole cultivated watermelon genome sequencing information offer unique opportunities for watermelon genomics and functional genomics research, and facilitate understanding regulatory networks involved in various biological processes under environment stress [12]. Transcriptome analyses have been performed to examine watermelon root responses to K + deficiency [13], but comparative transcriptome analyses of cultivated watermelon root responses to water deficit have not been reported and little is known about how watermelon roots respond to osmotic stress at the transcriptome level.
M08 is a cultivated inbred watermelon accession that exhibits relatively high resistance to water deficit [14]. The RNA-seq approach was used to investigate the root transcriptome profiles of hydroponically grown M08 after 6h PEG induced drought stress treatment. This study presents new insight into the molecular mechanisms of water deficit tolerance in the watermelon root system.

Plant material and osmotic stress treatment
Seeds of the watermelon inbred line M08 were germinated and planted in pots filled with a mixture of perlite and vermiculite (v/v = 1:1). Seedlings were cultivated in a greenhouse and 1/ 2 strength Hoagland's nutrient solution was applied daily to irrigate seedlings. When the second leaf of seedlings fully expanded, the perlite and vermiculite mixture attached to the seedling roots was carefully washed away and uniform seedlings were cultivated hydroponically in a chamber with a controlled environment. Opaque and black membrane-wrapped plastic boxes were used as hydroponic growth containers, and Hoagland's nutrient solution was used and replaced it every one week. Polystyrene foam plates with 15 holes of 25mm diameter were placed on plastic boxes. Each seedling that was individually wrapped with the sponge on its stem was planted in the hole and aerated continuously by air pumps. The hydroponic cultivation conditions were 28/20˚C, 16/8h of light/dark cycle, and relative humidity of 80%. When seedlings had 4-5 leaves, the PEG6000 was gradually added into Hoagland's solution and a final concentration of 20% (w/w, approximately -0.75MPa) was reached to simulate drought stress. The roots were harvested after 0, 3, 6, 12, and 24h treated by PEG6000. Phenotype of seedlings and transcript levels of Cla007307 (WRKY) and Cla006761 (MYB) were used to decide time point of mRNA sequencing. Samples at 6h (T-1, T-2 and T-3) were used to sequence mRNA and untreated samples (CK-1, CK-2 and CK-3) were used as controls. All samples were frozen in liquid nitrogen and stored at -80˚C until RNA extraction. Each treatment was comprised of five plant roots and three independent biological repeats were used for each treatment.

Quantitative RT-PCR (qRT-PCR) analysis
For qRT-PCR, specific primers were designed using Primer Premier 6.0, and ClUBCP (Cla010163) was used as an internal control [15]. qRT-PCR was carried out using the SYBR 1 Premix Ex Taq™ Kit (Takara, Dalian, China) on an StepOnePlus Real-Time PCR platform (Applied Biosystems, Foster City, CA, USA). The PCR reaction conditions were performed at 95˚C for 30 sec, followed by 40 cycles of 5 sec at 95˚C and 30 sec at 60˚C. Melting curves were used to verify PCR products. All reactions were performed in triplicate. The comparative ΔΔCt method was adopted to calculate the relative expression levels for each gene by normalizing the copies of target genes to the reference internal gene [16]. Primer pairs used in the qRT-PCR assays are listed in S1 Table. Library construction and RNA-seq analysis

Data processing and analysis
The raw RNA-seq reads were first filtered to eliminate adapter and low quality sequences. The high-quality clean reads obtained were mapped to watermelon reference sequences (http:// www.icugi.org/) using TopHat software with default parameters [17].
The mapped clean reads of each gene were counted and normalized into reads per million sequenced reads (RPKM) value using Cufflinks [18]. Differentially expressed genes (DEGs) between different treatment of samples were determined based on DESeq [19]. Genes were considered as DEGs if fold changes > = 2 were observed between different treatments and false discovery rates (FDRs) were < 0.01 [20].
Gene ontology (GO) analyses and GO enrichment were performed using the Gene Ontology database [21]. The Kyoto Encyclopedia of Genes and Genomics (KEGG) orthology database was adopted for pathway mapping [22].

Deciding appropriate time point for transcriptome sequencing
Cla017928 is a Δ 1 -pyrriline-5-carboxylate synthetase (P5CS) gene and Cla006761 is an orthologue of AtMYB77. The P5CS plays an important role in proline biosynthesis in higher plants [23] while AtMYB77 is involved in drought response [24]. For choosing the appropriate time point for transcriptome sequencing, the phenotype of seedlings and dynamic expression of Cla017928 (P5CS) and Cla006761 (MYB) were considered. Watermelon seedlings showed the symptoms of wilting at 6h and exhibited serious wilting at 12h and 24h ( Fig 1A). Cla006761 was up-regulated and reached peak level at 6h, while Cla017928 was down-regulated and reached the lowest at 6h (Fig 1B). Integrating phenotype of seedlings and dynamic expression of Cla017928 and Cla006761, the samples at 6h were used to RNA sequence.

Transcriptome sequencing and data analysis
Using the watermelon inbred line M08, RNA samples from 6h after 20% PEG6000 treated and untreated root tissues were used to create six independent libraries (CK-1, CK-2, CK-3, T-1, T-2 and T-3), which were sequenced using the Illumina HiSeq 2500 TM platform. After removing adaptors and low quality sequences, a total of 32.99 Gb clean and high quality data were obtained, and over 40 million (M) pair-end reads of 125 bp in length were obtained from each library ( Table 1). The ratio of unique mapped reads varied from 81.53% to 86.10%. The high quality data generated from six libraries provided a solid foundation for subsequent analyses.

Reliability of transcriptome sequencing data
To validate the reliability of transcriptome sequencing data, eighteen genes with various degrees of expression levels were evaluated by qRT-PCR. The data showed that the results of qRT-PCR were generally consistent with the transcriptome sequencing data (Fig 2A). Moreover, the linear regression equation y = 1.1155x -0.6465 with high correlation (R 2 = 0.9013) revealed that a positive correlation and significant similarity between the two analysis techniques ( Fig 2B). To further validate the reliability of RNA-seq data, correlations among biological replicates were assessed using the Pearson correlation coefficient (Fig 3). The libraries for the same treatment (i.e., biological replicates) were highly correlated. The weak correlation across treatments (CK and 20% PEG6000 treatments) suggests a large effect of water deficiency on the gene expression profile of watermelon root tissues.

Identification of differentially expressed genes (DEGs)
The total mapped reads were used to analyze differential expression in DESeq with FDR < 0.01 and fold change> = 2. In total, 5246 genes differentially expressed between PEGtreated and untreated watermelon root samples (S2 Table). Among these 5246 genes, 2753 were up-regulated and 2493 were down-regulated under water deficit stress. The DEGs comprised a large number of diverse genes. These results suggested that complex regulatory mechanisms were involved in watermelon root responses to water deficits.
BLAST2GO was used to retrieve DEGs in the GO database, and GO functional categorization was carried out using WEGO software. 4528 DEGs were classified in all three categories as follows: biological process (4148), molecular function (3670) and cellular component (4177) (Fig 4). It is not surprising that the terms "response to water deprivation (Corrected P-value 2.74E-3)," "hyperosmotic salinity response (Corrected P-value 3.83E-3)," "response to ethylene (Corrected P-value 9.31E-3)" appeared in the enriched biological process terms (S3 Table). Moreover, the highly enriched terms "regulation of plant-type hypersensitive response," "response to wounding," "defense response to fungus" and "respiratory burst involved in defense response" demonstrated crosstalk among diverse stress responses in watermelon roots, which were consistent with results in other plant species [25].
Plant hormones play important roles in plant responses to different stresses. In our data, most of the genes involved in "hormone-mediated signaling pathway," "response to ethylene" and "salicylic acid biosynthetic process" were up-regulated.
The protein folding is important for protein function. A variety of cues that disrupt protein folding in the endoplasmic reticulum lumen can activate the unfolded protein response, and even eventually cause programmed cell death (PCD) [26]. In our study, the highly enriched term of "endoplasmic reticulum unfolded protein response" indicated that osmotic stress affected protein folding. Accordingly, the roots initiated the mechanism to restore proper protein folding.
Autophagy PCD is activated in response to water deficits in the root apical meristem. In this way, the apical root dominance is deprived and the root system architecture is remodeled to adapt to water stress [27]. The GO terms "regulation of programmed cell death" and "root morphogenesis" were enriched in our data. These suggested that the PCD program may be activated in the root apical root meristem and root system architecture was altered under osmotic stress.
KEGG is a useful tool for the analysis of the roles of genes in various biological functions [28]. A total of 825 DEGs were assigned to 108 different biochemical pathways (S4 Table). The top five pathways comprising the most DEGs were "ribosome," "plant hormone signal transduction," "plant-pathogen interaction," "purine metabolism" and "starch and sucrose metabolism" (Fig 5). Of which, "ribosome" (Corrected P-value 5.80E-24) pathway was significantly changed. Most genes in the "ribosome" pathway were down-regulated, and these were consistent with the results of the GO analysis.

Protective genes in response to osmotic stress
Genes involved in osmotic adjustment. Under osmotic stress conditions, the accumulation of compatible solute molecules (such as proline, betaine and trehalose) plays a critical role   in maintaining cell membrane structure, cell turgor pressure and water status [23,29]. The P5CS is considered as a rate-limiting enzyme in proline biosynthesis in higher plants [23], and the ornithine aminotransferase (OAT) is a proline synthetase in the ornithine pathway of proline synthesis [30]. Surprisingly, P5CS (Cla017928) and OAT (Cla019569) were down-regulated by osmotic stress in roots. In contrast, Cla016474 which encodes proline dehydrogenase (ProDH) was up-regulated (S5 Table). Trehalose synthesis in plants mostly begins with the phosphorylated precursor, trehalose-6-phosphate (T6P) and then trehalose, with trehalose phosphate synthases (TPS) and trehalose phosphate phosphatases (TPP) [31]. TPS genes (Cla019181, Cla010101 and Cla009709) and TPP genes (Cla005675, Cl008123, Cla006270 and Cla014481) were up-regulated by osmotic stress. Genes involved in the metabolism of sucrose and galactinol showed significant up-regulation (S5 Table). Moreover, four amylase genes were up-regulated in response to osmotic stress.
Reactive oxygen species (ROS) scavenging. ROS act as important signal transduction molecules and toxic by-products of aerobic metabolism [32]. Abiotic stresses, such as drought, salinity, flooding, cold and heat can disorder the metabolic balance of cells and result in overproduction of ROS in plants. High concentrations of ROS damage proteins, lipids, carbohydrates and DNA owing to their reactivity and toxicity [33]. The respiratory burst oxidase homolog (Rboh) proteins are involved in the production of ROS [34]. To scavenge excess superoxide radicals, plants have developed highly efficient enzymes such as superoxide dismutase, catalase, peroxidase and ascorbate peroxidase, as well as non-enzymatic antioxidants such as ascorbic acid (ASH) and glutathione (GSH) [32]. In the present study, four genes encoding Rboh were notably up-regulated. DEGs involved in ascorbate peroxidase, glutathione reductase monodehydroascorbate reductase, dehydroascorbate reductase, glutathione peroxidase as well as most glutathione-S-transferase and peroxidase genes were up-regulated. However, the gene encoding superoxide dismutase (Cla011317) was down-regulated by water stress (S5 Table). These results suggested that osmotic stress might induce the complex antioxidant network. Moreover, the transcripts related to thioredoxins and glutaredoxins, which maintain the stable redox status of cells [35], were also abundant (S5 Table).
Other DEGs involved in osmotic protection. Multidrug and toxic compound extrusion proteins (MATE) are a group of newly characterized transporters that play crucial roles in stress responses [36,37]. In our experiment, Eighteen DEGs related to MATE transporters were detected (S5 Table). Aquaporins are membrane proteins that facilitate the transport of water across cell membranes and maintain cellular water homeostasis [38]. A total of eight genes encoding aquaporins showed differential expression under water stress in our study (S5 Table). The late embryogenesis abundant proteins (LEA) are important cellular dehydration protective proteins, and their expression levels are correlated with desiccation tolerance [39]. We found that two LEA (Cla015386 and Cla009416) were up-regulated by PEG treatment. Moreover, DEGs encoding molecular chaperones such as heat shock proteins (HSPs), chaperonins, and DnaJ-like were also differential expression.

DEGs involved in protein kinases, phosphatases and transcription factors
Osmotic stress elicits calcium signaling in plants [40]. Calmodulins (CaMs), CaM-like proteins (CMLs), calcineurin B-like proteins (CBLs) and Ca 2+ -dependent protein kinases (CDPKs) can sense calcium ion signals and then interact with their respective interacting partners under both biotic and abiotic stresses [41,42]. Some DGEs involved in Ca 2+ binding proteins were up-regulated under PEG treatment in the present study (S6 Table).
The mitogen-activated protein (MAP) kinase cascades, including mitogen-activated protein kinase kinase kinases (MAPKKKs), mitogen-activated protein kinase kinases (MAPKKs) and mitogen-activated protein kinases (MAPKs), are important in the response to osmotic stress [43]. A total of thirteen DEGs coding for MAP kinases were detected. Several DEGs related to SnF1-related protein kinases, serine/threonine protein kinases and phosphatase were up-regulated in our data (S6 Table).

Transcriptome of Watermelon Root under Osmotic Stress
DEGs referred to phospholipases, including phospholipase A (PLA), phospholipase C (PLC) and phospholipase D (PLD), were detected in the present study. In phospholipid signaling systems, the phospholipases catalyze the formation of messengers, such as inositol 1, 4, 5-trisphosphate (IP 3 ) and diacylglycerol (DAG), to modulate stress-responsive gene expression [40].
The expression of many functional genes is largely regulated by specific transcription factors. A large amount of transcription factors, such as MYB, bZIP, DREB/CBF, NAC and WRKY, are well characterized with respect to the regulation of biotic and abiotic stress responses [44]. The DEGs about transcription factor families are presented in S6 Table. Plant hormones Plant hormones mediate responses to both biotic and abiotic stresses and are important for plants to adapt to environment changes [45,46]. The DEGs related to phytohormone signaling including abscisic acid (ABA), auxin, cytokinin, gibberellic acid (GA), ethylene and jasmonic acid (JA) are listed in S7 Table. Cla009779 and Cla005404, which code for a 9-cis-epoxycarotenoid dioxygenase (NECD), a key enzyme for ABA biosynthesis [47], were up-regulated. The transcript levels of genes involved in ethylene synthesis, such as ERF transcription factors, 1-aminocyclopropane-1-carboxylate synthase (ACS) and 1-aminocyclopropane-1-carboxylate oxidase (ACO), were mostly up-regulated. Indole-3-acetic acid (IAA) is the main auxin in plants and is synthesized by tryptophan (Trp)-dependent or Trp-independent pathway [48]. The YUCCA (YUC) family of flavin monooxygenases and TRYPTOPHAN AMINOTRANSFERASE OF ARABIDOPSIS (TAA) family of amino transferases are essential for IAA biosynthesis in the Trp-dependent pathway [49,50]. Two YUCs and one TAA were down-regulated, as were the auxin transporter genes such as LAXs and PINs in our data. Bioactive GA has been recognized as an important phytohormone that regulates growth and development in plants [51]. In the present study, the gene (Cla021351) encoding ent-kaurenoic acid oxidase, the key enzyme of GA biosynthesis, was down-regulated, while the transcripts coding for GA 2-oxidase (GA2ox), which can inactivate endogenous bioactive GA [52], were up-regulated. This might result in decreases in the endogenous levels and bioactivity of GA. In addition, the transcripts about biosynthetic pathways of jasmonic acid and cytokinin were also affected by osmotic stress in the present study.

Other DEGs
Many DEGs involved in cell division and root growth, such as these encoding cyclin (CYC), cyclin-dependent kinase (CDK), microtubule-associated protein (MAP), Ras-like nuclear protein (Ran), E2 promoter-binding factor (E2F), DNA-binding-with-one-finger (DOF), ribosome, actin and tubulin were down-regulated (S8 Table). The ubiquitin 26S proteasome system (UPS) functions in removing misfolded or damaged proteins that may be produced after exposure to abiotic stress. In addition to playing an important role in UPS, the E3 ligases also regulate ABA-dependent stress signaling [53]. Numerous DEGs of E2-conjugating enzymes and E3-protein ligases emerged in our data (S9 Table), indicating that water stress could induce the UPS removal mechanism and E3 regulation.

Discussion
M08 is a cultivated watermelon inbred accession with relatively high resistance to water deficits. Drought exerts considerable influence on root tissue growth, prompting us to analyze the genes and pathways that regulate drought response and investigate the mechanism of adaptation under PEG-induced water deficit in M08 root tissues. RNA-seq was used to characterize drought responses in M08 root tissues at the transcriptome level and large-scale sequencing data were generated. In total, 5246 transcripts were differentially expressed in PEG-treated root samples as compared to untreated root samples. The DEGs were underwent GO and KEGG analyses, and the results provide novel insights into the mechanisms underlying water deficit stress in watermelon.

Osmotic adjustment in watermelon root tissue
Proline accumulation in plants plays an important role in response to environmental stresses owing to its osmoprotective function [54]. In plants, proline is synthesized from glutamate or ornithine, andP5CS, P5CR and OAT are important synthetases [55,56]. In our experiment, PEG treatment affected the two proline synthesis pathways. The genes encoding P5CS, P5CR and OAT were inhibited, while PDH that produces pyrroline-5-carboxylate (P5C) from proline was up-regulated under osmotic stress. Previous reports have shown that P5CS may be in or associated closely with the chloroplast [57,58] and the ability of proline synthesis is limited in the root apex under osmotic stress, but proline is transported from the leaf to the root [59]. Increased expression of PDH and low P5CS expression have been observed in some organs with high and increasing proline contents, e.g., in grape berries [60], maize and Arabidopsis root tip [59,61].
The expression levels of seven genes encoding TPP or TPS (Cla019181, Cla010101, Cla005675, Cla008123, Cla006270, Cla014481 and Cla009709) were up-regulated, and this may result in the accumulation of trehalose. Furthermore, the overexpression of some TPS genes improves the tolerance of rice seedlings to cold, drought and high salinity treatments [62]. Some OsTPS proteins interact with OsTPS1 and other OsTPS family members constituting protein complexes that potentially change trehalose-6-phospate (T6P) and trehalose levels, and then regulate stress responses [63].
Galactinol synthase (GolS) is involved in raffinose family oligosaccharides (RFO) metabolism and catalyzes an important step of the raffinose oligosaccharide biosynthetic pathway using galactose and myo-inositol as substrates [64]. RFO sugars play important roles in the drought tolerance of plants [29]. Overexpression of AtGolS2 in transgenic Arabidopsis enhances drought tolerance [65]. The gene (Cla010955) encoding GolS2 was up-regulated under osmotic stress in the present study.
Soluble carbohydrates are also produced by the hydrolysis of previously stored starch [66]. Water stress enhances the expression level of amylase, which catalyzes the breakdown of starch into soluble sugars [67,68]. Two α-amylases (Cla020676 and Cla010160) and two β-amylases (Cla021470 and Cla007635) were up-regulated in our study, suggesting that soluble carbohydrates derived from the hydrolysis of starch accumulated in watermelon roots under water stress.

ROS scavenging and detoxification
The overproduction of ROS in plants is often caused by various abiotic stresses and ultimately results in oxidative stress [69], and Rboh proteins are associated with ROS [70]. In our data, four genes encoding Rboh were up-regulated under water stress, while the transcripts of ROSscavenging enzymes and antioxidants were abundant. In addition, the thioredoxin system and glutaredoxin system were triggered by dehydration. Enhanced expression of the tomato glutaredoxin gene SlGRX1 in Arabidopsis increases abiotic tolerance against oxidative, drought, and salt stresses [71]. The thioredoxin system and glutaredoxin system are also induced in Populus euphratica by dehydration [68].

Transcription Factors
There are numerous transcription factors that regulate transcription of some genes in response to biotic and abiotic stresses. In our data, transcription factors involved in water stress were mostly enriched in AP2/ERF and WRKY, followed by MYB, bZIP, NAC and MYC family. A total of three DREBs, a sub-family of the AP2/ERF family involved in the ABA-independent pathway, were found in response to osmotic stress. Additionally, three bZIP family DEGs encoding abscisic acid-insensitive 5-like (ABI5-like) proteins which belong to the ABA-dependent pathway were up-regulated. These results showed that PEGinduced water stress activated the ABA-dependent pathway and ABA-independent pathway in watermelon roots.

Osmotic stress adaption
Environmental cues such as water, salinity and nutrients present many challenges for plant survival. Thus, plant roots have evolved to sense and integrate biotic cues in order to adjust their genetic program of post-embryonic root development [72]. Plant hormones such as ABA, auxins and gibberellins are involved in a complex signal system that plays important roles in both developmental regulation and environmental responses without mutual exclusive [73]. Auxin is a critical hormone for lateral root formation because it is involved in lateral root initiation, primordium development and lateral root emergence [74]. In contrast, ABA and ethylene play negative role in regulating Arabidopsis lateral root formation [75,76]. In our data, the transcripts encoding ACO, ACS, ETR1, EIN3, and NECD were up-regulated, while YUCCA, TAA, LAX, and PIN were repressed, which might lead to increases in ABA and ethylene accumulation, but a decrease in auxin accumulation.
In addition to plant hormones, root growth is also related to cell division, which progress involves the inactivation of cyclin (CYC) and cyclin-dependent kinase (CDK) [77]. CDKs bind to different CYCs to initiate the transition from post-mitotic interphase (G1) to the DNA synthetic phase (S) and the post-synthetic interphase (G2) to mitosis phase (M) [78]. In these progressions, the transcription factor E2 promoter-binding factor (E2F) induces genes that are required for cell cycle progression [79]. In the present experiment, all DEGs encoding cyclin, CDK and E2F were down-regulated under osmotic stress. On the other hand, the tubulin and actin are key components of the cytoskeleton and play regulatory roles in cell growth [80]. Most genes that referred to cytoskeleton components such as tubulin, actin and microtubuleassociated protein were down-regulated in our data. These results suggested that root cell division and growth were depressed. Ribosomes underlie the protein synthesis and this supports cell growth [81]. The KEGG pathway analysis suggested that the "ribosome" pathway was significantly changed in drought stress conditions. In addition, most DEGs involved in ribosomes were down-regulated (S8 Table). Moreover, the GO analysis showed that PCD was activated and this might result in depriving of the root apical dominance. Other biological processes in the GO analysis, such as "cell proliferation," "regulation of meristem growth," and "regulation of G2/M transition of mitotic cell cycle" were significantly enriched. Putting it all together, M08 may depress the root growth as a strategy to deal with short-term water stress. The same phenomenon has been observed in Arabidopsis roots under drought conditions [82]. According to our results and previous reports, we propose the hypothesis that plant roots may take defensive strategy and aggressive strategy to deal with osmotic stress. In term of the defensive strategy, plant may depress the root growth under the osmotic stress, just like in Arabidopsis roots [82] and our data. In contrast to defensive strategy, the aggressive strategy allows plants to induce root growth during drought stress [4,83].

Conclusions
To better understand the molecular mechanisms by which watermelon root tissues cope with short-term water stress, the watermelon inbred line M08 was subjected to 20% PEG6000 treatment and genome-wide differential gene expression was analyzed in comparison to untreated root tissues. RNA-seq technology was adopted to generate approximately 32.99 Gb of highquality data. This dataset provided elaborate gene expression profiles and enabled the identification of DEGs in response to water stress. The acclimation to short-term water stress of M08 root tissue involved osmotic adjustment, ROS scavenging, osmotic stress signal transduction and the inhibition of root growth (Fig 6). Our experiment provided novel insight into the molecular mechanisms of watermelon root tissues in response to water stress and improved our understanding of watermelon coping with water deficits. Supporting Information S1