Morphological and physiological responses to drought
The overall morphological changes of oil palm seedlings to drought stress with a period of 14 days were observed in both leaves and roots. The effects of drought stress were first observed in leaf morphology, showing initial edge and tip necrosis and then wilting and yellowing for the drought treated samples (Figure 1a). In terms of roots, in comparison to the control, the drought treatments showed significant decrease in the number of roots, root volume and overall biomass (Figure 1b). Trypan blue staining further showed that the drought treatment roots experienced not only obvious cell deformation but also more cell membrane injury than that of the well-watered controls (Figure 1c). Overall, these observations are consistent with those of previous studies in plant species under drought stresses and undergoing water deprivation [18]. These results indicate substantial physiological responses of the oil palm seedlings under drought stress [19], and provided useful starting materials to study the genes, pathways and networks involved in drought responses of the important crop species using RNA-seq [20] and bioinformatics analysis.
DEGs in drought responses of oil palm
An average of 70.2 million (M) cleaned reads were obtained across all six samples. The controls had higher reads coverage than the drought treatments (84.8 M vs 55.7 M) (Supplementary Table S1). Nevertheless, the sequence coverage of the treatments (> 200x of transcriptomes) was sufficient to construct transcripts and identify differentially expressed genes. Approximately 70% of cleaned reads were uniquely mapped to the reference genome of 31,640 annotated protein coding genes [21]. The drought stressed seedlings showed slightly higher uniquely mapping rates than the controls (72.1% vs 66.3%), indicating the duplicated genes also play important roles in response to drought stress [22] in oil palm, a species of palaeotetraploid origin [21] and future studies should also focus on paralogous genes and their potential functions in stress responses [23].
A total of 2,084 and 1,358 DEGs were identified using two different approaches: DESeq2 and EdgeR, respectively, within which 1,293 were shared by the two data sets (Figure 2). In detail, DESeq2 identified 944 down-regulated and 1,140 up-regulated DEGs, while EdgeR screened 624 down-regulated and 734 up-regulated DEGs. The number of common down- and up-regulated DEGs were 614 and 679, respectively, between the two approaches (Figure 2 & Supplementary Table S2). To obtain confident results, only the common DEGS between these two approaches were kept for further analysis. Based on the relative expression of DEGs across samples, the drought treatment and control samples were clearly differentiated by both PCA and hierarchical clustering analyses and showed substantial differences in expression profiles (Figure 3). We further assessed the accuracy of the RNAseq data set by comparing to qPCR of randomly selected nine genes (Supplementary Table S3). We observed an overall high consistency of the expression patterns of these genes between RNA-seq and qPCR (Supplementary Figure 1a), with a correlation coefficient of 0.978 (P < 0.0001), as examined using Pearson’s correlation test (Supplementary Figure 1b). In total, these data indicate that the RNA-seq data is solid.
Interestingly, we observed most of the DEGs in subcategories phenylalanine metabolism and tryptophan metabolism to be down-regulated (Table 1). As studied in plant species, phenylalanine and tryptophan metabolisms are more involved in pathogen related immune responses [38, 39]. The down regulation of most DEGs within these categories likely implies the effects of metabolic compensation to drought stress responses by sacrificing the less important biological functions. Except for these enriched pathways, we still found two genes: two-component response regulator ORR9 and two-component response regulator ORR24, that were down- and up-regulated, respectively, were enriched into the subcategory: zeatin biosynthesis (Table 1), which is also revealed to play important roles in drought stress response in Populus simonii [40]. Further studies are required to know how the expression patterns of the two genes are related to responses to drought stress in oil palm.
Ontology enrichment analysis of genes responding to drought stress
To understand the transcriptomic responses to the drought stress, we first carried out gene ontology enrichment analysis using the 1,293 DEGs consistently identified by both EdgeR and DESeq2. A total of 89 GO terms were significantly enriched, involving a number of categories of diverse functions (Supplementary Table S4). The most significant enrichment entities included GO terms related to cell wall biogenesis and functions (e.g., GO: 0009834, GO: 0044036, GO: GO:0009664, GO: 0016998 and GO: 2000652), which was consistent with the above observation that the cell wall of the drought treatments has likely been damaged by severe drought stress and thus has triggered the mechanism of damage and repair (Figure 4). Moreover, we also observed significant enrichments related to phenylpropanoid biosynthesis and metabolisms (e.g., ath00940 and GO: 0009698 and GO: 0046271). It has been widely studied that phenylpropanoid pathway is activated by stress conditions, such as drought, salinity and extreme temperature, and leads to accumulation of phenolic compounds, which play critical physiological roles in regulation under abiotic stress to cope with environmental challenges [24]. Moreover, we found some GO terms classified into the groups related to ion transport and homeostasis (e.g., GO: 0006811, GO: 0030004 GO: 0030007, GO: 0015698 and GO: 0034220) and response to osmotic stress and water homeostasis (e.g., GO:0006970 and GO: 0030104). Differential expression of genes in these functions likely results from the responses of plants to water deprivation by direct regulation of osmotic pressure [25, 26]. In addition, a number of genes were enriched into the biological categories related to regulation of cellular ketone metabolic process (GO: 0010565), suggesting that genes involved in ketone metabolic process played important roles in drought stress in oil palm [27, 28]. Hormone regulations are also indispensable to stress responses of plant species. Here, we identified two enriched GO terms related to hormone regulation and metabolism (e.g., GO:0010817 and GO: 0042447). Previous studies have shown that production of numerous secondary metabolites is essential for physiological processes to respond to abiotic stress [29, 30]. Consistent with these results, we found several significant enrichments related to these terms: small molecule biosynthetic process (GO: 0044283), amino sugar and nucleotide sugar metabolism (ath00520), galactose metabolism (ath00052), benzene-containing compound metabolic process (GO: 0042537), linoleic acid metabolism (ath00591) and xyloglucan metabolic process (GO: 0010411). Interestingly, we also identified significant enriched GO terms, like response to jasmonic acid (GO:0009753) and ABC transporters (ath02010), playing crucial roles in abiotic stress responses (Figure 4).
The interactions of these enriched GO terms were further investigated using network analysis and eight enriched GO networks were identified, with each consisting of no less than 3 genes (Figure 5). The major GO networks involved those related to cell wall related biogenesis and metabolism (GO: 0009834 and GO: 0044036), small molecule related biosynthetic and metabolic processes (ath00940, GO: 0009698, GO: 0044283 and GO: 0010565) and ion transport and homeostasis related processes (GO: 0006811, GO: 0034220 and GO: 0030004). These data imply that genes in these networks are more extensively induced to differentially express to respond to drought stress [25, 30, 31]. We further investigated the enriched KEGG pathways and found that the functions of the enriched pathways were generally consistent with that of the enriched GO terms as shown above (Supplementary Table S5). Above all, these enrichment analyses suggest that many genes, pathways and networks respond to the drought stress in the roots of oil palm seedlings. The DEGs, pathways and networks identified in this study provide valuable resources for future studies on their functions to improve drought tolerance of oil palm.
Plant hormone signal transduction in drought stress responses
Plant hormones not only play crucial roles in controlling growth and development, but also are indispensable in regulation of stress responses [32]. Herein, we first focused on the DEGs involved in plant hormone signal transduction pathway and found significant enrichments of DEGs within subcategories of KEGG pathways including a-Linolenic acid metabolism, carotenoid biosynthesis, phenylalanine metabolism, tryptophan metabolism and zeatin biosynthesis (Table 1). Genes related to a-Linolenic acid metabolism have often been revealed to play important roles in drought stress responses [33, 34]. We found four DEGS were involved in a-Linolenic acid metabolism and three out of them were up-regulated. Interestingly, the down-regulated DEG, jasmonic acid-amido synthetase JAR1 (LOC105048226), was a duplicated copy of the up-regulated one, jasmonic acid-amido synthetase JAR1 (LOC105046997), suggesting functional divergence of paralogous genes since genome duplication events. Nevertheless, consistent with the expression patterns of most DEGs in this subcategory, the up-regulated jasmonic acid-amido synthetase JAR1 might be more important in regulation of drought stress response in oil palm. Moreover, we found all of the DEGs, including probable protein phosphatase 2C 24 and two duplicated copies of probable protein phosphatase 2C 75, in the subcategory carotenoid biosynthesis, were up-regulated. These three DEGs were also enriched into the subcategory abscisic acid pathway (ABA) within MAPK signalling pathway (Table 1). Carotenoid biosynthesis signalling pathway is specifically induced by root and contributes to induce ABA production to regulate ion homeostasis, as studied in Arabidoposis [35]. ABA-independent signalling pathways have been widely studied to be involved in the regulation of drought stress response in many kinds of plant species [36]. These results suggest that ABA related genes also play important roles in drought stress responses of oil palm.
ABC transporters in drought responses
Membrane transporters play vital roles in regulation of water and ion homeostasis of organisms, among which ATP-binding cassette (ABC) transporters constitute one of the largest protein families and act as both exporters and importers, driven by ATP hydrolysis [37]. An increasing number of studies have shown that ABC transporters play irreplaceable roles in transmembrane allocations of various molecules to adapt to rapidly changing environments, such as water scarcity, heavy metal stress and pathogen stress [38]. In order to survive in these changing abiotic conditions, it is necessary for cells to absorb nutritious chemical substances and discharge endogenous toxins, as well as exchange signalling molecules [38]. Thus, the ABC transporters occupy a diverse range of functions and hence the regulations upon stress responses are also complicated. Here, we found six DEGs were enriched into the pathway of ABC transporters (Table 1). Five of them were ABCB subfamily members, among which three (ABC transporter B family member 11, ABC transporter B family member 19 and putative multidrug resistance protein (LOC105038824)) were down-regulated and two (ABC transporter B family member 9 and another putative multidrug resistance protein (LOC105060251)) were up-regulated. Such differential expression patterns of these ABCB subfamily transporters indicate the complicated functions in controlling of influx and efflux of chemical molecules [39, 40]. Other than these ABC transporters, we also identified an ABCC subfamily member, ABC transporter C family member 5, which was up-regulated. Interestingly, two putative multidrug resistance protein genes (LOC105038824 and LOC105060251) were differentially expressed against drought stress. As shown in previous studies, multidrug resistance‐associated proteins are widely involved in regulation of stress responses, such as salt stress, water deprivation, oxidative stress and fungal stress [41]. Taken together, these different types of ABC transporters likely play important roles in responses to drought stress of oil palm.
Protein-protein interaction networks in response to drought responses
Other than significantly enriched GO terms and KEGG pathways, we also identified three protein-protein interaction networks, focused on ion transport, reactive nitrogen species metabolic process and nitrate assimilation (Figure 6 & Table 2). Eight DEGs were involved in the ion transport network, among which five (ammonium transporter 2 member 1 (AMT2-1), amino acid transporter ANT1 (ANT1), cation/H(+) antiporter 20 (CHX20), plasma membrane ATPase 4 (PMA4) and potassium channel AKT1 (AKT1)) and three (receptor-like protein kinase HSL1 (HSL1), plasma membrane ATPase (PMA) and ABC transporter G family member 42 (ABCG42)) were up- and down-regulated, respectively (Table 2). Interestingly, most of the cation channel and transporter genes were up-regulated, including ammonium transporter 2 member 1 (AMT1), amino acid transporter ANT1 (ANT1), cation/H(+) antiporter 20 (CHX20) and potassium channel AKT1 (KT1), indicating their positive effects in regulating ion homeostasis in oil palm [42]. Nevertheless, we also observed three DEGs were down-regulated in the same network, implying both positive and negative feedback regulations are acting on this network [43]. Reactive nitrogen species metabolic process is also suggested to have critical roles in stress responses, such as drought and salinity [44]. Consistently, we identified three up-regulated genes: magnesium transporter MRS2-1 (MGT2), putative chloride channel-like protein CLC-g (AT5G33280) and serine/threonine protein kinase OSK1 (KIN10), in this protein-protein interaction network. Nitrate assimilation is another biological process affecting salt and water stress tolerance in plants [45]. Here, we found four DEGs involved in this network: two were up-regulated (cationic amino acid transporter 6, chloroplastic (CAT6) and sodium/hydrogen exchanger 4 (NHX4)), while the other two were down-regulated (amino acid permease 8 (AAP8) and vacuolar cation/proton exchanger 1a (CAX1)). As these protein-protein interaction networks play crucial roles in drought stress response, the DEGs involved in these networks provide excellent candidates to improve drought tolerance of oil palm by genetic engineering and/or selective breeding.
Transcription factors in drought responses
To date, more and more studies have paid attention to the biological functions of transcription factors as regulatory elements binding proteins [10]. Transcription factors are vital for development, response to intercellular and environmental signals and pathogenesis [10]. The expression changes are often associated with important cellular processes [7]. In this study, we identified 96 differentially expressed transcription factors that were classified into 28 families (Supplementary Table S6 & Table 3). Previous studies have shown that transcription factors are broadly involved in drought/abiotic stress responses, such as the members of family MYB, WRKY, DREB, NAC and AP2/EREBP [16, 46-48]. Here we also observed genes in these transcription factor families were differentially expressed against drought stress in oil palm, further supporting their important roles in drought tolerance in plant species. Interestingly, we found several families of transcription factors that were rarely studied and involved in abiotic stress responses, such as the C2H2, LFY and TALE transcription families. Therefore, it is important to understand the mechanisms of regulatory functions of these genes, which might be useful to help improve drought tolerance of related plant species.