Skip to main content

Advertisement

Comparative physiology and transcriptome analysis allows for identification of lncRNAs imparting tolerance to drought stress in autotetraploid cassava

Article metrics

Abstract

Background

Polyploidization, pervasive among higher plant species, enhances adaptation to water deficit, but the physiological and molecular advantages need to be investigated widely. Long non-coding RNAs (lncRNAs) are involved in drought tolerance in various crops.

Results

Herein, we demonstrate that tetraploidy potentiates tolerance to drought stress in cassava (Manihot esculenta Crantz). Autotetraploidy reduces transpiration by lesser extent increasing of stomatal density, smaller stomatal aperture size, or greater stomatal closure, and reducing accumulation of H2O2 under drought stress. Transcriptome analysis of autotetraploid samples revealed down-regulation of genes involved in photosynthesis under drought stress, and less down-regulation of subtilisin-like proteases involved in increasing stomatal density. UDP-glucosyltransferases were increased more or reduced less in dehydrated leaves of autotetraploids compared with controls. Strand-specific RNA-seq data (validated by quantitative real time PCR) identified 2372 lncRNAs, and 86 autotetraploid-specific lncRNAs were differentially expressed in stressed leaves. The co-expressed network analysis indicated that LNC_001148 and LNC_000160 in autotetraploid dehydrated leaves regulated six genes encoding subtilisin-like protease above mentioned, thereby result in increasing the stomatal density to a lesser extent in autotetraploid cassava. Trans-regulatory network analysis suggested that autotetraploid-specific differentially expressed lncRNAs were associated with galactose metabolism, pentose phosphate pathway and brassinosteroid biosynthesis, etc.

Conclusion

Tetraploidy potentiates tolerance to drought stress in cassava, and LNC_001148 and LNC_000160 mediate drought tolerance by regulating stomatal density in autotetraploid cassava.

Background

Cassava (Manihot esculenta Grantz) is a diploid plant (2n = 2× = 36) and an important cash and energy crop cultivated in Asia, Africa, and Latin America for its storage roots, making it critical for food security and economic development [1]. Water scarcity harms production when cassava is cultivated in severely water-deficit regions, and although cassava can tolerant a wide range of adverse environmental conditions including drought, this can limit growth and survival [2].

Drought stress increases oxidative damage in plants [3] and reduces photosynthesis [4, 5]. To cope with this, plants have evolved complex mechanisms such as deeper and thicker root systems [6], stomatal modulation systems including reduced aperture size and/or density to reduce water loss from transpiration [7, 8], and accumulation of osmotic adjustment compounds [9]. At low to moderate concentrations, reactive oxygen species (ROS) such H2O2 may act as second messengers in stress signalling. However, excessive H2O2 production can trigger progressive oxidative damage to plant cells. Antioxidant enzymes scavenge excess H2O2 and other ROS to protect plant cells from damage [10, 11].

Doubling of the whole genome to generate polyploidy is ubiquitous among higher plant species, and the change to a polyploidy state increases abiotic stress tolerance in crop species. For example, tetraploidy in cenchrus, Arabidopsis and paulownia improves drought tolerance by lowering H2O2 accumulation and enhancing ROS clearance [12,13,14]. Lower stomatal conductance and the abscisic acid (ABA) signalling pathway are involved in drought tolerance in the leaves of autotetraploid Rangpur lime (Citrus limonia) [15]. Autopolyploidization increases the potassium content and promotes tolerance to salinity stress in Arabidopsis [16]. Altered anatomy induced by polyploidy may confer stronger drought tolerance upon autotetraploid cassava [17]. However, the physiological and molecular advantages underlying these adaptations have not been widely investigated.

Transcriptome analyses have confirmed subtle changes due to autopolyploidy in Arabidopsis and Paulownia tomentosa × Paulownia fortunei under drought stress [13, 14]. Studies of long non-coding RNAs (lncRNAs), which are longer than 200 nucleotides in length and lack a coding sequence, have expanded our understanding of eukaryote transcriptome [18]. LncRNAs play important roles in many different biological processes in plants [19]. Thousands of lncRNAs associated with drought responses have been identified in several crop species due to the rapid development of omics sequencing technologies [20,21,22,23,24]. However, to date, only a limited number of lncRNAs, including COOLAIR, COLDAIR, npc536, IPS1, LDMAR, PMS1T, DRIR, ELENA1 and TL, have been cloned and characterised [19, 25,26,27,28,29,30,31]. Overexpressing lncRNA DRIR can enhance drought tolerance in Arabidopsis [29].

In the current study, we characterised physiological differences between diploid progenitor (2×) and autotetraploid (4×) cassava under standard and drought conditions. Comparative transcriptome analysis was used to identify differentially expressed genes (DEGs) and lncRNAs in the leaves of 2× and 4× cassava under well-watered and dehydration conditions. Co-expression analysis revealed that two differentially expressed (DE) lncRNAs regulated six DEGs that improve drought tolerance in cassava by modulating stomatal density.

Results

Autotetraploid cassava displays stronger drought tolerance than diploid plants

We exposed ‘Xinxuan 048’ 2× and 4× plants to soil with the same relative water content to compare their responses to drought stress. The RMSC was 30% for controls (Fig. 1a) and 4.5% for the drought stress treatment (Fig. 1b). After 15 days of withholding watering, at which time the RSMC was 4.5%, leaves of the 2× plants displayed moderate drooping, while leaves of 4× plants displayed only slight drooping (Fig. 1b). After withdrawing water for 30 days, followed by 2 days of recovery, all 2× plants were dead, while the growing points of 4× plants survived, and the upper leaves of autotetraploids remained green (Fig. 1c). Furthermore, the RWC of detached 4× leaves was higher than that of 2× plants (Fig. 2a). Analysis of water loss rate showed that 2× detached leaves lost water much more rapidly than those of 4× plants (Fig. 2b). Moreover, the transpiration rate of 4× seedlings was significantly lower than that of 2× seedlings (Fig. 2c). Thus, overall, 4× cassava plants were significantly more drought-tolerant than 2× plants.

Fig. 1
figure1

Phenotype changes in cassava plants in response to drought stress. On the left are 2× cassava ‘Xinxuan 048’ plants, and 4× plants are shown on the right. a Cassava plants (100-day-old) under normal conditions with 30% RSMC, b Plants after 15 days of water deprivation, with 4.5% RSMC, c Plants after water deprivation for 30 days (1% RSMC) and recovery for 2 days

Fig. 2
figure2

4× cassava plants enhanced drought tolerance by reducing transpiration. a Comparison of the RWC, b water losses percentage in detached leaves, c transpiration rate of seedlings in 3-month-old 2× and 4× cassava plants. Values are means ± SD (*p < 0.05, **p < 0.01, Student’s t-tests)

Autotetraploidy alters drought-mediated stomatal function and photosynthetic capacity

Since water loss occurs mainly through stomatal movement, the reduced water loss in 4× plants prompted us to investigate the stomatal function of the two genotypes under control and 4.5% RSMC conditions. The stomatal density of 4× plants (7.24) showed a significant reduction compared with that of 2× plants (10.64) under control conditions, and the density of both genotypes increased with increasing drought stress, but the density of 4× plants (7.88) remained significantly lower than that of 2× plants (12.24) in dehydrated leaves (Fig. 3a). The average stomatal aperture in 2× plants under control conditions was 0.12 μm, compared with 0.09 μm in dehydrated 2× plants. The average stomatal aperture in control 4× plants was ~ 0.16 μm, but this dropped 0.06 μm under drought condition (Fig. 3b). Thus, 4× leaves showed enhanced stomatal closure in response to drought, and stomatal conductance was markedly decreased in 4× leaves compared with 2× leaves under 4.5% RSMC (Fig. 3c). These results indicate that stress tolerance in 4× plants may be due to reduced transpiration rate via lesser extent increasing of stomatal density, smaller stomatal aperture size, and/or greater stomatal closure.

Fig. 3
figure3

Changes in stomatal function and photosynthetic parameters induced by drought stress. a Stomatal density, b Stomatal aperture, c Stomatal conductance, d Net photosynthetic rate, e SPAD value. Values are means ± SD. Different letters represent significant differences at p < 0.05 (Duncan’s tests)

Next, we investigated photosynthesis parameters of the two genotypes. Photosynthesis is expected to be hindered due to lack of CO2 under drought stress. Although the stomatal conductance of 2× plants was greater than that of 4× plants under drought stress, the net photosynthetic rate of 4× plants was significantly higher (Fig. 3d). This result is consistent with the higher SPAD value of 4× plants, although the difference was not significant (Fig. 3e).

Physiological effects of autotetraploidy

Six physiological traits were used to investigate dynamic changes in response to drought stress in 2× and 4× plants. An increase of more than 150% in H2O2 production was observed in 2× leaves under drought stress, compared with less than 50% in drought-treated 4× leaves (Fig. 4a). By contrast, the relative increase in T-AOC and CAT activity following drought treatment was much less in 2× plants than 4× plants (Fig. 4b and c). The MDA content is an indicator of the degree of lipid peroxidation, which can reflect damage to plant cell membranes [32]. Compared with controls, the MDA content was dramatically elevated in 2× dehydrated leaves but decreased in 4× dehydrated leaves, although the reduced value is not significant (Fig. 4d). The proline and soluble sugar content were also increased significantly in 2× plants with increasing drought stress, but were relatively unchanged in 4× plants (Fig. 4e and f). These results indicate that osmolyte biosynthesis is not influenced by tetraploidy. Overall, tetraploidy resulted in reducing accumulation of H2O2, and hence the ability to alleviate cell membrane injury, thereby increasing drought tolerance in 4× plants.

Fig. 4
figure4

Changes in the physiological traits of 2× and 4× plants in response to drought stress. a H2O2 content, b T-AOC, c CAT activity, d MDA content, e Proline content, f Soluble protein content. Values are means ± SD. Different letters represent significant differences at p < 0.05 (Duncan’s tests)

Identification of the lncRNAs and hierarchical clustering

In total, 1,232,088,086 bp raw reads were generated of 12 libraries by paired end sequencing. After moving the adapters and low-quality reads, 1,182,713,428 bp clean reads were mapped to the cassava genome. The mapped reads for each sample were assembled into transcripts using a reference-based approach. The correlation coefficient of the expression quantity of all transcriptome was more than 0.824 between each other among the 12 samples (Additional file 1: Figure S1). Expression levels of all coding genes and lncRNA transcripts were systematically estimated using FPKM values. In total, 2372 lncRNAs, including 821 antisense_lncRNAs and 1551 lincRNAs, were identified in this study (Additional file 2: Table S1). To calculate the degree of differential expression (DE) among lncRNAs, hierarchical clustering was performed using FPKM values of lncRNAs under control and drought stress conditions in 2× and 4× leaves. The results suggest that tetraploidization may have a limited effect on the mRNA transcriptome and lncRNA expression, since 2XCK and 4XCK clustered together, and 2XDR and 4XDR formed another cluster (Fig. 5a). Remarkably, significant DE was observed following drought stress in both diploid and autotetraploid cassava, suggesting that DE lncRNAs may play an important role in responses to drought stress (Fig. 5b).

Fig. 5
figure5

Hierarchical clustering of a mRNAs and b lncRNAs under control and drought stress conditions in 2× and 4× plants

LncRNAs as potential miRNA precursors

By aligning miRNA precursors to the 2372 lncRNAs, we identified 18 lncRNAs as 21 known cassava miRNA precursors, including miR162, miR166c, miR408 and miR477c (Additional file 3: Table S2). A single lncRNA could serve as a precursor of several miRNAs, and different miRNAs could be targeted by the same lncRNA. Lnc_001314-miR171a was DE in response to drought.

Comparison of mRNA transcripts under drought stress

In order to investigate transcriptional changes under water deficit conditions, we performed two pairwise transcriptomes comparisons; (2XDR_vs._2XCK) vs. (4XDR_ vs._ 4XCK), and (4XCK_vs._2XCK) vs. (4XDR_vs._2XDR).

Using a q-value < 0.05 and a fold change > 1 or < − 1 as thresholds, 1562 DEGs were found to be specifically responsive to drought in 4× sample (Fig. 6a), of which 687 genes were up-regulated and 875 were down-regulated (Additional file 4: Table S3). A total of 5484 genes were commonly expressed (Fig. 6a; Additional file 5: Table S4). We identified 2412 DEGs in 2× samples, including 1032 up-regulated and 1380 down-regulated genes (Fig. 6a; Additional file 6: Table S5). In 4× versus 2× plants subjected to dehydration, 814 DEGs were identified (Fig. 6b), including 288 down-regulated and 526 up-regulated DEGs (Additional file 7: Table S6). The reliability of the deep sequencing data was validated by quantitative real time PCR (qPCR) with gene-specific primers (Additional file 8: Table S7) for six randomly selected mRNAs (Additional file 9: Figure S2). The results showed that all the six genes displayed similar expression patterns in both RNA-seq and qPCR data. Notably, the number of drought-responsive mRNAs in 2× plants (2412) was higher than in 4× plants (1562), suggesting that 2× leaves might be more sensitive to drought than those of 4× plants, consistent with the more pronounced phenotype differences in 2× plants under water deficit stress.

Fig. 6
figure6

Venn diagrams showing the number of specific and common a DEGs and c DE lncRNAs in 2× and 4× plants in response to drought. The number of b DEGs and d DE lncRNAs in the comparison between 4× stressed leaves vs. 2× stressed leaves and 4× vs. 2× controls

A closer inspection of the 5484 common genes identified a number of down-regulated genes. In particular, almost all the genes encoding subtilisin-like proteases, except MANES_05G206800 and MANES_06G139800, were more strongly down-regulated in 2× leaves than in 4× leaves following dehydration (Table 1). We used the qPCR to detect the expression of seven genes in the four samples, the results showed that MANES_05G206800 and MANES_06G139800 turned out to be opposite trend with the RNA-seq data, while the other five were consistent (Additional file 10: Figure S3). The discrepancy between qPCR and RNA-seq could be attributed to different statistical method. Consistently, a subtilisin-like protease gene mutant, sdd1–1, exhibited a two- to four-fold increase in stomatal density in Arabidopsis [33]. IiSDD1 was found to be autopolyploidy responsive and down-regulated by drought stress in autopolyploidy Isatis indigotica [34]. In the present work, the stomatal density was elevated in stressed 4× leaves (8.84%), but to a lesser extent than in 2× leaves (15.04%; Fig. 3a). Interestingly, the four subtilisin-like protease genes including, MANES_18G044300, MANES_18G039200, MANES_17G062600 and MANES_06G013200 were also found to be up-regulated expressed in 4× versus 2× plant subjected to dehydration (Additional file 7: Table S6). These results indicate a correlation between the regulation of stomatal density-related genes, resulting in a lesser extent of increasing stomatal density under dehydration stress, and consequently enhanced drought tolerance in 4× plants.

Table 1 Fold change in expression levels of subtilisin-like protease genes in 2× and 4× cassava plants under drought stress conditions

Among the 5484 common genes that were down-regulated, genes involved in photosynthesis were particularly pronounced, including ferredoxin, ferredoxin-NADP reductase, beta-carbonic anhydrase, chlorophyll a-b binding protein, fructose-bisphosphate aldolase, photosynthetic NADH dehydrogenase subunit, and PsbP domain-containing protein 3. Expression of these genes was reduced in dehydrated 4× leaves more than in 2× leaves (Table 2), suggesting that photosynthesis was affected more in 4× plants than in 2× plants following drought treatment, consistent with the more pronounced reduction in stomatal aperture in 4× leaves under drought stress (Fig. 3b).

Table 2 Changes in differentially expressed genes involved in photosynthesis in 2× and 4× cassava leaves under drought stress conditions

Regulating the expression of transcription factors (TFs) is an efficient strategy for amplifying transcriptional responses. As shown in Table 3, members belonging to the WRKY, MYB, and ERF TF families were represented. In particular, all seven ERF family members were down-regulated in 4× plants under drought stress, and all 14 WRKY members except MANES_10G127100 and MANES_11G066500 were up-regulated under drought stress (Table 3). Consistently, these three families of TFs appear to play important roles in drought stress signalling [35,36,37].

Table 3 Members of three transcription factor families differentially expressed in 2× and 4× cassava plants under drought stress conditions

Comparison of lncRNA transcripts under drought stress

Similarly, (2XDR_vs._2XCK) vs. (4XDR_ vs._ 4XCK) and (4XCK_vs._2XCK) vs. (4XDR_vs._2XDR) comparisons identified DE lncRNAs. Among the lncRNAs with a q-value <0.05 and a fold change > 1 or < − 1, 69 DE lncRNAs were specific to drought-stressed 4× leaves (Fig. 6c), including 45 up-regulated and 24 down-regulated DE lncRNAs (Additional file 11: Table S8). A total of 138 DE lncRNAs were identified in both genotypes (Fig. 6c; Additional file 12: Table S9), including 104 drought-responsive lncRNAs specific to 2× plants (Fig. 6c), of which 72 were up-regulated while 32 were down-regulated (Additional file 13: Table S10). A Venn diagram (Fig. 6d) showed among 17 DE lncRNAs, 11 were down-regulated and six were up-regulated in 4× stressed leaves vs. 2× stressed leaves (Additional file 14: Table S11). The sequences of DE lncRNAs are listed in Additional files 15 and 16: Tables S12 and S13.

In order to explore the functions of the 86 lncRNAs specific to 4× plants, we constructed a co-expression network and identified target genes associated with drought tolerance in 4× vs. 2× plants under drought stress (Additional files 17 and 18: Tables S14 and S15). We used the qPCR to detect the expression of four lncRNAs with the corresponding trans target genes in the two samples, Me4XDR and Me4XCK. The results indicated that the expression is consistent with the RNA-seq result (Fig. 7). Remarkably, we found that two of the DE lncRNAs might regulate six subtilisin-like protease DEGs; LNC_001148 (lincRNA) may regulate MANES_18G039200 (SBT1.7), MANES_06G139800 (SBT1.7), MANES_05G175300 (SBT1.6) and MANES_06G020400 (SBT2.5), while LNC_000160 (antisense lncRNA) may target MANES_16G094900 (SBT1.7) and MANES_18G044300 (SBT1.7). To validate the putative relationship between the two DE lncRNAs and six DEGs, expression levels were examined by qPCR. The results revealed similar lncRNA expression patterns to those obtained in the RNA-seq (Fig. 8), which suggests that the two lncRNAs identified using deep sequencing may target genes encoding subtilisin-like proteases via co-expression.

Fig. 7
figure7

The qPCR validation of the expression of four lncRNAs with the co-expressed target genes between Me4XDR and Me4XCK. Data are means ± SD of three biological experiments. Cassava β-actin was used as an internal control

Fig. 8
figure8

Comparison of expression results from RNA-Seq and qPCR methods for the two lncRNAs and co-expressed subtilisin-like protease family members. Data are means ± SD of three biological experiments. Cassava β-actin was used as an internal control

The functions of most of the 86 DE lncRNAs remain largely unknown, therefore, the genes co-expressed with 86 DE lncRNAs overlapping with DEGs were subjected to functional enrichment analysis using KEGG (Additional file 19: Table S16). The results indicated that 20 pathways were found responsive to drought stress. They included galactose metabolism, pentose phosphate pathway, plant hormone signal transduction, glycolysis/gluconeogenesis, biosynthesis of secondary metabolites and brassinosteroid biosynthesis, etc. (Additional file 20: Table S17).

Discussion

Autopolyploidy confers advantages over diploid counterparts, often manifested in enhanced adaptation to adverse environmental conditions. However, the molecular advantages remain largely unknown. Herein, we demonstrated the importance of autotetraploidy in response to drought in cassava. The 4× cassava plants potentiated stronger antioxidant and photosynthesis capacities than 2× plants, helping them cope better with drought stress. Autotetraploidy elevated the expression levels of two lncRNAs co-expressed with genes encoding subtilisin-like proteases regulating stomatal density under drought stress. Additionally, autotetraploidy decreased the stomatal density on the leaf epidermis, which in turn improves the RWC and reduces water loss, enhancing drought tolerance in cassava plants.

ROS participate in signal transduction, but are also toxic to cell membranes when present in excessive quantities, and can affect plant drought susceptibility [38]. We observed that H2O2 in 4× stressed leaves were significantly lower than in 2× leaves (Fig. 4a). To detoxify excess drought-induced ROS, plants have developed a complex antioxidant system [39]. We found that 4× plants possessed a more efficient enzymatic antioxidant system involving CAT and T-AOC for reducing accumulation of H2O2 under drought stress (Fig. 4b and c). MANES_05G1307001, homolog of catalase isozyme 1 in Ricinus communis, is up regulated in both of 2× and 4× stressed leaves compared with normal condition. However, the fold change of MANES_05G1307001 in 4× is more than that in 2×. Therefore, autotetraploidy may be associated with the regulation of antioxidation ability. Some antioxidant enzymes were found to be up-regulated in stressed leaves, while others were down-regulated, suggesting cellular redox status may be complex and dynamic. We found that levels of UDP-glucosyltransferase were increased more or reduced less in leaves of 4× plants in response to drought compared with 2× leaves (Additional file 5: Table S4). Consistently, overexpression of Arabidopsis UGT79B2/B3 significantly enhances plant tolerance to drought stress by modulating anthocyanin accumulation, which enhances ROS scavenging [40]. Our results suggest that DEGs encoding UDP-glucosyltransferase may be tightly associated with improved drought tolerance in 4× by modulating ROS levels.

Stomata are surrounded by two guard cells in the epidermis that regulate the shape and size of the pore aperture [41]. One of the earliest adaptive responses to drought in cassava leaves is stomatal closure and/or decreased stomatal density to reduce water loss [42]. Stomata control the balance between the uptake of CO2 for photosynthesis and the release of water by modulating transpiration, thereby governing water use and abiotic stress tolerance [43]. It has been demonstrated that reducing the number of stomata does not affect carbon fixation due to increased CO2 concentration [44]. Therefore, the decreased photosynthetic capacity of 4× plants is consistent with a lower transpiration rate, which results, at least in part, from enhanced stomatal closure during drought treatment (Fig. 3b-d). It is known that photosynthesis-related genes are down-regulated after drought treatment in many plants [45]. However, we found that the net photosynthesis rate of 4× plants remained higher than that of 2× plants under dehydration conditions, which might be attributed to higher SPAD values in 4× plants (Fig. 3e). Thus, 4× plants may maintain a higher photosynthetic rate, facilitating better adaptation to drought stress conditions by meeting increased energy demand.

Our transcriptomic data indicates that autotetraploidy influences the expression of genes encoding TFs involved in drought stress in cassava, and general response mechanisms integrating hormone signalling in response to external stimuli. TFs are efficient amplifiers of transcriptomic responses that regulate differential gene expression attributed to autotetraploidy. ERFs are key regulatory hub proteins in hormone and regulatory ROS-responsive gene expression that confer abiotic stress tolerance [37]. Some WRKY and MYB TFs are components of ABA-mediated stomatal movement involved in drought responses [36, 46]. Three DE TF families were found to be in autotetraploid Arabidopsis under drought stress conditions [13].

A complex regulatory system controls drought tolerance in cassava. Most previous research in this area has focused on coding genes. LncRNAs are an important class of regulators in diverse biological processes involving complex mechanisms [47], and numerous lncRNAs have been identified in a few crop species including wheat (Triticum aestivum) [48], Medicago truncatula [49], Brassica napus [50], maize (Zea mays) [51] and cotton (Gossypium spp.) [52]. In our current study, we discovered 2372 lncRNAs, including 821 antisense_lncRNAs and 1551 lincRNAs, from 12 libraries. Li et al. [53] identified 682 lncRNAs from nine cassava samples. Differences between our current results and these previous results might be attributed to (i) differences in experimental design and transcriptome analysis, since the work of Li et al. [53] was based on 15-day-old seedling tissues under polyethylene glycol-simulated drought stress, and samples were collected from shoot apices and leaves; (ii) differences in bioinformatics strategies, since CPAT [54] and Pfam Scan [55, 56] were employed to filter transcripts with coding potential and screen candidate lncRNAs. LncRNAs could execute their functions to respond to stress in either cis-acting or trans-acting in the genome via diverse mechanisms in plant [57, 58]. Based on comparative transcriptome analysis, lncRNA16397 was found to be involved in Phytophthora infestans resistance by co-expression glutaredoxin in tomato (Lycopersicon esculentum Mill.) [59]. In the present study, LNC_001148 and LNC_000160, 898 bp and 2688 bp, respectively, were down-regulated by drought treatment in both 2× and 4× cassava plants. Co-expression analysis revealed that they appeared to regulate six genes encoding subtilisin-type proteinases. SDD1, belong to subtilisin serine proteinase family, which are known to negatively regulate stomatal density and distribution in Arabidopsis [33]. IiSDD1 participates not only in the drought responsive pathways, but also involves in autopolyploidy I. indigotica evolution [34]. Thus, relative to 2× dehydrated leaves, the higher expression of LNC_001148 and/or lower expression of LNC_000160 in 4× leaves may result in less down-regulation of target genes encoding subtilisin-type proteinases, and hence lesser extent of increasing stomatal density, thereby enabling 4× plants to better adapt to drought stress. Exactly how the two lncRNAs confer drought tolerance in cassava will be studied in future work.

In Table S14, LNC_000211 was down regulated with fold change more than 8 times following drought treatment in 4× plant specifically. LNC_000211 could trans regulate MANES_16G111000 (encoding Dehydration-responsive element-binding protein), MANES_11G139300 (encoding ERF), MANES_18G039200 (encoding Subtilisin-like protease SBT1.7), MANES_08G148200 (encoding thioredoxin-like 1–1), respectively. Previous studies indicated that these four target genes were involved in drought tolerance in plants [33, 37, 60, 61]. Therefore, LNC_000211 may play an important role in mediating the tolerance to drought in 4× cassava plant. Based on the KEGG enrichment analyses, we found that the trans target genes of the 86 DE lncRNAs were involved in galactose metabolism, pentose phosphate pathway, brassinosteroid biosynthesis, etc. The three pathways were reported to be associated with drought responsive in sugarcane (Saccharum officinarum L.), purging nut (Jatropha curcas) and Arabidopsis, respectively [62,63,64]. Taken together, our results suggest that 4× cassava implements divergent mechanisms to modulate the response to drought stress.

The current understanding of lncRNA regulation in response to drought stress is in its infancy in 4× cassava. These findings provide a comprehensive view of 4× DE lncRNAs, which will enable in-depth functional analysis.

Conclusion

Our study demonstrated that tetraploidy potentiates tolerance to drought stress in cassava. The co-expressed network analysis indicated that LNC_001148 and LNC_000160 in autotetraploid dehydrated leaves regulated six genes encoding subtilisin-like protease, thereby result in increasing the stomatal density to a lesser extent in autotetraploid cassava. This study helps to explain the role of autotetraploidy in conferring drought tolerance, and indicates the evolutionary potential of polyploidy in abiotic stress tolerance.

Methods

Water deficit treatment and water recovery

Cassava variety ‘Xinxuan 048’ (2× and 4×) used in this study were original from our lab [65]. Stem-propagated plants were grown in plastic pots (30 cm in height × 40 cm in diameter) containing well-mixed soil in March 2017. Each pot contained one cutting, and were placed in a greenhouse under a 16 h light/8 h dark photoperiod at the Guangxi Academy of Agricultural Sciences (GXAAS). Two-month-old cassava plants of each genotype were well-watered before drought stress treatment. Before dehydration treatment, each potted plant was watered until it was saturated to ensure consistency of water content. The moisture content of soil was measured by a AZS-100 soil moisture sensor (TRIME-PICO32, Germany). Control plants were well watered every 4 days. A relative soil moisture content (RSMC) of 30, 4.5% (after 15 days of withholding watering) and 1% (after 30 days of withholding watering until the top buds of diploid cassava displayed obvious wilting) was used for controls and drought stress treatments for two time points, respectively. After recovery for 2 days, the survival rate was determined and plants were photographed. All experiments were repeated in triplicate.

Measurement of water loss rate

Seventy-day-old plants were used to calculate the water loss rate, for which the fourth, fifth and sixth leaves (counting from the top of the plant) and petioles were excised from each plant. Three plants were included for each genotype. Detached leaves were placed on filter paper in a culture room under a 16 h light/8 h dark photoperiod. The abaxial surface was placed facing up for dehydration, and leaves were weighed every hour. Water loss was estimated from the percentage of fresh weight lost relative to the initial fresh weight. All experiments were repeated in triplicate.

Relative water content (RWC) assay

Four plants each of each genotype were used to measure the RWC. Seventy-day-old plants were detached at the fifth leaf (counting from the top of the plant) and weighed immediately (M1), then placed in water for 18 h under dark conditions, dried using a filter paper, then weighed again (M2). Saturated leaves were finally oven-dried for 24 h at 65 °C to a constant weight (M3), and the RWC was measured using eq. 1:

$$ \left(\mathrm{M}1\hbox{-} \mathrm{M}3\right)/\left(\mathrm{M}2\hbox{-} \mathrm{M}3\right)\times 100\% $$
(1)

Drought stress treatment and plant sampling

One hundred-day-old cassava plants were used for drought stress analysis. Before dehydration treatment, each potted plant was watered until it was saturated to ensure consistency of water content. The moisture content of soil was measured by a AZS-100 soil moisture sensor (TRIME-PICO32, Germany). Control plants were well watered every 4 days. RSMC of 30 and 4.5% (after 15 days of withholding watering) was used for controls and drought stress treatments, respectively, and all treatments included three plants per genotype. The fifth leaf of each plant was used to measure stomatal aperture size, photosynthetic capacity, and physiological indices. The lower epidermis of leaves was used for examination as described previously by Zhou et al. [65]. The length and width of stomata were measured using 200 guard cells, and the number of stomatal openings was counted in 10 microscopic fields using three independent replicates. A LI-6400 portable photosynthesis system was used to measure the net photosynthetic rate, transpiration rate, and stomatal conductance. The relative chlorophyll content in the fifth leaves of each plant was determined by a soil plant analysis development (SPAD) value measured by an SPAD-502 Plus chlorophyll meter (Konica Minolta, Japan).

Physiological indices, namely catalase (CAT) activity, total antioxidant capacity (T-AOC), H2O2 content, and malondialdehyde (MDA) content, were estimated using reagent kit cat. # A007–2, A015, A064 and A003 (Institute of Nan Jing Jian Cheng Bioengineering), respectively. Proline content and total soluble protein content were measured using the method of Bates et al. [66] and Guy et al. [67], respectively.

Plants used for RNA-seq were treated as described above, and parallel leaves were frozen in liquid nitrogen and stored at − 80 °C.

Library construction and deep sequencing

A total of 3 μg of RNA per sample was used as input material for RNA sample preparation. Whole-transcriptome library preparation and high-throughput sequencing were performed by Novogene Bioinformatics Technology Co., Ltd. (Beijing, PR China). Three replicates were generated for each control (2XCK and 4XCK) and treatment (2XDR and 4XDR). A total of 12 libraries were constructed, which included Me2XCK-1, Me2XCK-2, Me2XCK-3, Me2XDR-1, Me2XDR-2, Me2XDR-3, Me4XCK-1, Me4XCK-2, Me4XCK-3, Me4XDR-1, Me4XDR-2 and Me4XDR-3. All the libraries were sequenced on an Illumina Hiseq 2500 platform, and 125 bp paired-end reads were generated. Clean reads were obtained by removing reads containing adapters and poly-N sequences, and low-quality reads were also removed from raw data. Clean reads were aligned to the reference genome using TopHat v2 [68]. Both scripture (beta2) [69] and Cufflinks (v2.1.1) [70] were employed to assemble mapped reads for each sample into transcripts using a reference-based approach.

Coding potential was predicted using one or more of CNCI [71], CPC [72], Pfam-scan [73, 74] and PhyloCSF [75], and sequences without coding potential were considered candidate lncRNAs.

Identification of DEGs and DE lncRNAs

Cuffdiff (v2.1.1) was used to calculate fragments per kilobase of transcript per million mapped reads (FPKM) values of both lncRNAs and coding genes in each sample [70]. Gene FPKMs were computed by summing the FPKMs of transcripts in each gene group. To assess the three biological replicates, the log2 (fold_change)-transformed FPKM values were performed. Cuffdiff provides statistical routines for determining DE in digital transcript using a model based on the negative binomial distribution [70]. A default q-value < 0.05 was set as the threshold for DE. Genes with log2 (fold_change) values > 0 were deemed up-regulated, while genes with log2 (fold_change) values < 0 were considered down-regulated.

Trans target genes prediction of lncRNAs

Trans targets genes of lncRNAs were identified from expression correlations between lncRNAs and coding genes using custom scripts. The pearson correlation coefficient method was used to calculate the correlation between lncRNA and mRNA among samples. Correlations corresponding to a coefficient > 0.95 or < − 0.95 were selected for the functional enrichment analysis.

KEGG enrichment analysis

We used KOBAS software to test the statistical enrichment of the genes co-expressed with DE lncRNAs overlapping with DEGs in KEGG pathways (http://www.genome.jp/kegg/) [76].

Quantitative real time PCR validation

Total RNA was used to synthesise cDNA using a PrimeScript RT Reagent Kit (TaKaRa). Three technical replicates and three biological replicates were included for each experiment, and qPCR was performed using SYBR Premix Ex Taq [77].

Availability of data and materials

Not applicable.

Abbreviations

ABA:

abscisic acid

CAT:

catalase

DE:

differentially expressed / differential expression

DEG:

differential expressed gene

LncRNA:

long non-coding RNA

MDA:

malondialdehyde

ROS:

reactive oxygen species

RSMC:

relative soil moisture content

RWC:

relative water content

SD:

standard deviation

SPAD:

soil plant analysis development

TF:

transcriptional factor

References

  1. 1.

    FAO. Why cassava? http://www.fao.org/ag/agp/agpc/gcds/index_en.html, last accessed 6 July 2016. 2014.

  2. 2.

    Burns A, Gleadow R, Cliff J, Zacarias A, Cavagnaro T. Cassava: the drought, war and famine crop in a changing world. Sustainability. 2010;2(11):3572–607.

  3. 3.

    Munné-Bosch S, Jubany-Mari T, Alegre L. Drought-induced senescence is characterized by a loss of antioxidant defences in chloroplasts. Plant Cell Environ. 2001;24(12):1319–27.

  4. 4.

    Havaux M, Tardy F. Loss of chlorophyll with limited reduction of photosynthesis as an adaptive response of Syrian barley landraces to high-light and heat stress. Funct Plant Biol. 1999;26(6):569–78.

  5. 5.

    Kyparissis A, Petropoulou Y, Manetas Y. Summer survival of leaves in a soft-leaved shrub (Phlomis fruticosa L., Labiatae) under Mediterranean field conditions: avoidance of photoinhibitory damage through decreased chlorophyll contents. J Exp Bot. 2000;46(12):1825–31.

  6. 6.

    Hu H, Xiong L. Genetic engineering and breeding of drought-resistant crops. Annu Rev Plant Biol. 2014;65:715–41.

  7. 7.

    Berry JA, Beerling DJ, Franks PJ. Stomata: key players in the earth system, past and present. Curr Opin Plant Biol. 2010;13(3):233–40.

  8. 8.

    Casson SA, Hetherington AM. Environmental regulation of stomatal development. Curr Opin Plant Biol. 2010;13(1):90–5.

  9. 9.

    Chen TH, Murata N. Enhancement of tolerance of abiotic stress by metabolic engineering of betaines and other compatible solutes. Curr Opin Plant Biol. 2002;5(3):250–7.

  10. 10.

    Hou X, Xie K, Yao J, Qi Z, Xiong L. A homolog of human ski-interacting protein in rice positively regulates cell viability and stress tolerance. Proc Natl Acad Sci. 2009;106(15):6410–5.

  11. 11.

    Fang YJ, Xiong LZ. General mechanisms of drought response and their application in drought resistance improvement in plants. Cell Mol Life Sci. 2015;72(4):673–89.

  12. 12.

    Chandra A, Dubey A. Effect of ploidy levels on the activities of D1-pyrroline-5-carboxylate synthetase, superoxide dismutase and peroxidase in Cenchrus species grown under water stress. Plant Physiol Biochem. 2010;48(1):27–34.

  13. 13.

    del Pozo JC, Ramirez-Parra E. Deciphering the molecular bases for drought tolerance in Arabidopsis autotetraploids. Plant Cell Environ. 2014;37(12):2722–37.

  14. 14.

    Xu E, Fan GQ, Niu SY, Zhao ZL, Deng MJ, Dong YP. Transcriptome-wide profiling and expression analysis of diploid and autotetraploid Paulownia tomentosa × Paulownia fortunei under drought stress. PLoS One. 2014;9(11):e113313.

  15. 15.

    Allario T, Brumos J, Colmenero-Flores JM, Iglesias DJ, Pina JA, Navarro L, Talon M, Ollitrault P, Morillon R. Tetraploid Rangpur lime rootstock increases drought tolerance via enhanced constitutive root abscisic acid production. Plant Cell Environ. 2013;36(4):856–68.

  16. 16.

    Chao DY, Dilkes B, Luo HB, Douglas A, Yakubova E, Lahner B, Salt DE. Polyploids exhibit higher potassium uptake and salinity tolerance in Arabidopsis. Science. 2013;341(6146):658–9.

  17. 17.

    Nassar NM, Graciano-Ribeiro D, Fernandes SD, Araujo PC. Anatomical alterations due to polyploidy in cassava, Manihot esculenta Crantz. Genet Mol Res. 2008;7(2):276–83.

  18. 18.

    Rinn JL, Chang HY. Genome regulation by long noncoding RNAs. Annu Rev Biochem. 2012;81:145–66.

  19. 19.

    Chekanova JA. Long non-coding RNAs and their functions in plants. Curr Opin Plant Biol. 2015;27:207–16.

  20. 20.

    Qi X, Xie S, Liu Y, Yi F, Yu J. Genome-wide annotation of genes and noncoding RNAs of foxtail millet in response to simulated drought stress by deep sequencing. Plant Mol Biol. 2013;83(4–5):459–73.

  21. 21.

    Di C, Yuan J, Wu Y, Li J, Lin H, Hu L, Zhang T, Qi Y, Gerstein MB, Guo Y, Lu ZJ. Characterization of stress-responsive lncRNAs in Arabidopsis thaliana by integrating expression, epigenetic and structural features. Plant J. 2014;80(5):848–61.

  22. 22.

    Shuai P, Liang D, Tang S, Zhang Z, Ye CY, Su Y, Xia X, Yin W. Genome-wide identification and functional prediction of novel and drought-responsive lincRNAs in Populus trichocarpa. J Exp Bot. 2014;65(17):4975–83.

  23. 23.

    Zhang W, Han Z, Guo Q, Liu Y, Zheng Y, Wu F, Jin W. Identification of maize long non-coding RNAs responsive to drought stress. PLoS One. 2014;9(6):e98958.

  24. 24.

    Li WB, Li CQ, Li SX, Peng M. Long noncoding RNAs that respond to Fusarium oxysporum infection in ‘Cavendish’ banana (Musa acuminata). Sci Rep. 2017;7:16939.

  25. 25.

    Franco-Zorrilla JM, Valli A, Todesco M, Mateos I, Puga MI, Rubio-Somoza I, Leyva A, Weigel D, García JA, Paz-Ares J. Target mimicry provides a new mechanism for regulation of microRNA activity. Nat Genet. 2007;39(8):1033–7.

  26. 26.

    Amor BB, Wirth S, Merchan F, Laporte P, d'Aubenton-Carafa Y, Hirsch J, Maizel A, Mallory A, Lucas A, Deragon JM, Vaucheret H, Thermes C, Crespi M. Novel long non-protein coding RNAs involved in Arabidopsis differentiation and stress responses. Genome Res. 2009;19(1):57–69.

  27. 27.

    Ding JH, Lu Q, Ouyang YD, Mao HL, Zhang PB, Yao JL, Xu CG, Li XH, Xiao JH, Zhang QF. A long noncoding RNA regulates photoperiod-sensitive male sterility, an essential component of hybrid rice. Proc Natl Acad Sci. 2012;109(7):2654–9.

  28. 28.

    Fan YR, Yang JY, Mathioni SM, Yu JS, Shen J, Yang XF, Wang L, Zhang QH, Cai Z, Xu CG, Li XH, Xiao JH, Meyers BC, Zhang QF. PMS1T, producing phased small-interfering RNAs, regulates photoperiod-sensitive male sterility in rice. Proc Natl Acad Sci. 2016;113(52):15144–9.

  29. 29.

    Qin T, Zhao HY, Cui P, Albesher N, Xiong LM. A nucleus-localized long non-coding RNA enhances drought and salt stress tolerance. Plant Physiol. 2017;175(3):1321–36.

  30. 30.

    Seo JS, Sun HX, Park BS, Huang CH, Yeh SD, Jung C, Chua NH. ELF18-INDUCED LONG-NONCODING RNA associates with mediator to enhance expression of innate immune response genes in Arabidopsis. Plant Cell. 2017;29:1024–38.

  31. 31.

    Liu X, Li DY, Zhang DL, Yin DD, Zhao Y, Ji CJ, Zhao XF, Li XB, He Q, Chen RS, Hu SN, Zhu LH. A novel antisense long noncoding RNA, TWISTED LEAF, maintains LEAF blade flattening by regulating its associated sense R2R3-MYB gene in rice. New Phytol. 2018;218(2):774–88.

  32. 32.

    Møller IM, Jensen PE, Hansson A. Oxidative modifications to cellular components in plants. Annu Rev Plant Biol. 2007;58(1):459–81.

  33. 33.

    Berger D, Altmann T. A subtilisin-like serine protease involved in the regulation of stomatal density and distribution in Arabidopsis thaliana. Genes Dev. 2000;14(9):1119–31.

  34. 34.

    Xiao Y, Yu XJ, Chen JF, Di P, Chen WS, Zhang L. IiSDD1, agene responsive to autopolyploidy and environmental factors in Isatis indigotica. Mol Biol Rep. 2010;37:987–94.

  35. 35.

    Chen L, Song Y, Li S, Zhang L, Zou C, Yu D. The role of WRKY transcription factors in plant abiotic stresses. Biochim Biophys Acta. 2012;1819(2):120–8.

  36. 36.

    Baldoni E, Genga A, Cominelli E. Plant MYB transcription factors: their role in drought response mechanisms. Int J Mol Sci. 2015;16(7):15811–51.

  37. 37.

    Müller M, Munné-Bosch S. Ethylene response factors: a key regulatory hub in hormone and stress signaling. Plant Physiol. 2015;169(1):32–41.

  38. 38.

    Miller G, Suzuki N, Ciftci-Yilmaz S, Mittler R. Reactive oxygen species homeostasis and signalling during drought and salinity stresses. Plant Cell Environ. 2010;33(4):453–67.

  39. 39.

    Noctor G, Foyer CH. ASCORBATE AND GLUTATHIONE: Keeping active oxygen under control. Annu Rev Plant Physiol Plant Mol Biol. 1998;49:249–79.

  40. 40.

    Li P, Li YJ, Zhang FJ, Zhang GZ, Jiang XY, Yu HM, Hou BK. The Arabidopsis UDP-glycosyltransferases UGT79B2 and UGT79B3, contribute to cold, salt and drought stress tolerance via modulating anthocyanin accumulation. Plant J. 2017;89(1):85–103.

  41. 41.

    Bergmann DC, Sack FD. Stomatal development. Annu Rev Plant Biol. 2007;58:163–81.

  42. 42.

    Turyagyenda LF, Kizito EB, Ferguson M, Baguma Y, Agaba M, Harvey JJ, Osiru DS. Physiological and molecular characterization of drought responses and identification of candidate tolerance genes in cassava. AoB Plants. 2013;5:plt007.

  43. 43.

    Hetherington AM, Woodward FI. The role of stomata in sensing and driving environmental change. Nature. 2003;424(6951):901–8.

  44. 44.

    Woodward FI. Stomatal numbers are sensitive to increases in CO2 from pre-industrial levels. Nature. 1987;327(6123):617–8.

  45. 45.

    Mohammadi PP, Nouri M-Z, Komatsu S. Proteome analysis of drought-stressed plants. Current Proteomics. 2012;9(4):232–44.

  46. 46.

    Rushton DL, Tripathi P, Rabara RC, Lin J, Ringler P, Boken AK, Langum TJ, Smidt L, Boomsma DD, Emme NJ, Chen XF, Finer JJ, Shen QX, Rushton PJ. WRKY transcription factors: key components in abscisic acid signalling. Plant Biotechnol J. 2012;10(1):2–11.

  47. 47.

    Liu X, Hao LL, Li DY, Zhu LH, Hu SN. Long non-coding RNAs and their biological roles in plants. Genomics Proteomics Bioinformatics. 2015;13(3):137–47.

  48. 48.

    Xin MM, Wang Y, Yao YY, Song N. Identification and characterization of wheat long non-protein coding RNAs responsive to powdery mildew infection and heat stress by using microarray analysis and SBS sequencing. BMC Plant Biol. 2011;11:61.

  49. 49.

    Wang TZ, Liu M, Zhao MG, Chen RJ, Zhang WH. Identification and characterization of long non-coding RNAs involved in osmotic and salt stress in Medicago truncatula using genome-wide high-throughput sequencing. BMC Plant Biol. 2015;15:131.

  50. 50.

    Joshi RK, Megha S, Basu U, Rahman MH, Kav NN. Genome wide identification and functional prediction of long non-coding RNAs responsive to Sclerotinia sclerotiorum infection in Brassica napus. PLoS One. 2016;11(7):e0158784.

  51. 51.

    Lv YD, Liang ZK, Ge M, Qi WC, Zhang TF, Lin F, Peng ZH, Zhao H. Genome-wide identification and functional prediction of nitrogen-responsive intergenic and intronic long non-coding RNAs in maize (Zea mays L.). BMC Genomics. 2016;17:350.

  52. 52.

    Wang MJ, Yuan DJ, Tu LL, Gao WH, He YH, Hu HY, Wang PC, Liu N, Lindsey K, Zhang XL. Long noncoding RNAs and their proposed functions in fibre development of cotton (Gossypium spp.). New Phytol. 2015;207(4):1181–97.

  53. 53.

    Li SX, Yu X, Lei N, Cheng ZH, Zhao PJ, He YK, Wang WQ, Peng M. Genome-wide identification and functional prediction of cold and/or drought-responsive lncRNAs in cassava. Sci Rep. 2017;7:45981.

  54. 54.

    Wang L, Park HJ, Dasari S, Wang S, Kocher JP, Li W. CPAT: coding-potential assessment tool using an alignment-free logistic regression model. Nucleic Acids Res. 2013;41(6):e74.

  55. 55.

    Eddy SR. A new generation of homology search tools based on probabilistic inference. Genome informatics. 2009;23(1):205–11.

  56. 56.

    Finn RD, Bateman A, Clements J, Coggill P, Eberhardt RY, Eddy SR, Heger A, Hetherington K, Holm L, Mistry J, Sonnhammer E, Tate J, Punta M. Pfam: the protein families database. Nucleic Acids Res. 2014;42:D222–30.

  57. 57.

    Zhu QH, Wang MB. Molecular functions of long non-coding RNAs in plants. Genes. 2012;3(1):176–90.

  58. 58.

    Wang J, Meng X, Dobrovolskaya OB, Orlov YL, Chen M. Non-coding RNAs and their roles in stress response in plants. Genome Proteomics Bioinformatics. 2017;15(5):301–12.

  59. 59.

    Cui J, Luan YS, Jiang N, Bao H, Meng J. Comparative transcriptome analysis between resistant and susceptible tomato allows the identification of lncRNA16397 conferring resistance to Phytophthora infestans by co-expressing glutaredoxin. Plant J. 2017;89(3):577–89.

  60. 60.

    Broin's M, Cuiné S, Peltier G, Rey P. Involvement of CDSP 32, a drought-induced thioredoxin, in the response to oxidative stress in potato plants. FEBS Lett. 2000;467(2–3):245–8.

  61. 61.

    Yang YL, Liao WB, Yu XL, Wang B, Peng M, Ruan MB. Overexpression of MeDREB1D confers tolerance to both drought and cold stresses in transgenic Arabidopsis. Acta Physiol Plant. 2016;38:243.

  62. 62.

    Sapeta H, Lourenço T, Lorenz S, Grumaz C, Philipp Kirstahler P, Barros P, Costa JM, Sohn K, Oliveira M. Transcriptomics and physiological analyses reveal co-ordinated alteration of metabolic pathways in Jatropha curcas drought tolerance. J Exp Bot. 2016;67(3):845–60.

  63. 63.

    Yang YT, Fu ZW, Su YC, Zhang X, Li GY, Guo JL, Que YX, Xu LP. A cytosolic glucose-6-phosphate dehydrogenase gene, ScG6PDH, plays a positive role in response to various abiotic stresses in sugarcane. Sci Report. 2014;4:7090.

  64. 64.

    Northey J, Liang SY, Jamshed M, Deb S, Foo E, Reid J, McCourt P, Samuel M. Farnesylation mediates brassinosteroid biosynthesis to regulate abscisic acid responses. Nature plants. 2016;2:16114.

  65. 65.

    Zhou HW, Zeng WD, Yan HB. In vitro induction of tetraploids in cassava variety ‘Xinxuan 048’ using colchicine. Plant Cell Tissue Organ Cult. 2017;128(3):723–9.

  66. 66.

    Bates L, Waldren R, Teare I. Rapid determination of free proline for water-stress studies. Plant Soil. 1973;39(1):205–7.

  67. 67.

    Guy C, Haskell D, Neven L, Klein P, Smelser C. Hydration-state-responsive proteins link cold and drought stress in spinach. Planta. 1992;188(2):265–70.

  68. 68.

    Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL. TopHat2: parallel mapping of transcriptomes to detect InDels, gene fusions. Genome Biol. 2013;14(4):R36.

  69. 69.

    Guttman M, Garber M, Levin JZ, Donaghey J, Robinson J, Adiconis X, Fan L, Koziol MJ, Gnirke A, Nusbaum C, Rinn J, Lander E, Regev A. Ab initio reconstruction of cell type-specific transcriptomes in mouse reveals the conserved multi-exonic structure of lincRNAs. Nat Biotechnol. 2010;28(5):503–10.

  70. 70.

    Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, Salzberg SL, Wold BJ, Pachter L. Transcript assembly and quantification by RNA-seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28(5):511–5.

  71. 71.

    Sun L, Luo HT, Bu DC, Zhao GG, Yu KT, Zhang CH, Liu YN, Chen RS, Zhao Y. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Res. 2013;41(17):e166.

  72. 72.

    Kong L, Zhang Y, Ye ZQ, Liu XQ, Zhao SQ, Wei L, Gao G. CPC: assess the protein-coding potential of transcripts using sequence features and support vector machine. Nucleic Acids Res. 2007;36:W345–9.

  73. 73.

    Bateman A, Birney E, Cerruti L, Durbin R, Etwiller L, Eddy SR, Griffiths-Jones S, Howe KL, Marshall M, Sonnhammer EL. The Pfam protein families database. Nucleic Acids Res. 2002;30:276–80.

  74. 74.

    Punta M, Coggill PC, Eberhardt RY, Mistry J, Tate J, Boursnell C, Pang N, Forslund K, Ceric G, Clements J, Heger A, Holm L, Sonnhammer EL, Eddy SR, Bateman A, Finn RD. The Pfam protein families database. Nucleic Acids Research, Database Issue. 2012;40:D290–301.

  75. 75.

    Lin MF, Jungreis I, Kellis M. PhyloCSF: a comparative genomics method to distinguish protein coding and non-coding regions. Bioinformatics. 2011;27(13):i275–82.

  76. 76.

    Kanehisa M, Araki M, Goto S, Hattori M, Hirakawa M, Itoh M, Katayama T, Kawashima S, Okuda S, Toshiaki T, Yamanishi Y. KEGG for linking genomes to life and the environment. Nucleic Acid Research. 2008;36:D480–4.

  77. 77.

    Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2–ΔΔCT method. Methods. 2001;25(4):402–8.

Download references

Acknowledgments

We thank Dr. Mengbing Ruan at Institute of Tropical Bioscience and Biotechnology, Chinese Academy of Tropical Agricultural Sciences for proof-reading.

Funding

This work was supported by the Fund for Guangxi Crop Genetic Improvement and Biotechnology Key Lab (17–259-55), Guangxi Natural Science Foundation Youth Science Foundation Project (2017GXNSFBA198223) and Guangxi Natural Science Foundation Project (2016GXNSFBA38010). These funding bodies did not directly participate to the design of the study and collection, analysis, interpretation of data and writing the manuscript.

Author information

LX and HBY conceived and designed the work. LX conducted the experiments. LX analyzed the data and wrote the paper. SBCH revised the paper. XHS, WDZ and LYL participated in the preparation of the materials. XYX and SC participated in the cultivation and management of the plants. All authors read and approved the manuscript.

Correspondence to Hua-Bing Yan.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1:

Figure S1. The pearson correlation coefficient of lncRNA between the 12 samples in this study. (TIF 19628 kb)

Additional file 2:

Table S1. Expression levels of 2372 lncRNAs detected in this study. (XLSX 649 kb)

Additional file 3:

Table S2. LncRNAs as miRNA precursors. (XLSX 10 kb)

Additional file 4:

Table S3. The 1562 DEGs specific to 4× plants in response to drought. (XLSX 150 kb)

Additional file 5:

Table S4. The 5484 DEGs common to both 2× and 4× plants in response to drought. (XLSX 467 kb)

Additional file 6:

Table S5. The 2412 drought-responsive DEGs specific to 2× plants. (XLSX 244 kb)

Additional file 7:

Table S6. The 814 DEGs in 4× stressed leaves vs. 2× stressed leaves. (XLSX 91 kb)

Additional file 8:

Table S7. Sequences of qPCR primers used in this study. (XLSX 12 kb)

Additional file 9:

Figure S2. Validation of the expression of six mRNAs selected randomly identified by RNA-seq using qPCR. Six mRNAs were selected randomly from 4XDR and 4XCK libraries. Bars represent means ± SD of three biological replicates. Cassava β-actin was used as an internal control. (TIF 11225 kb)

Additional file 10:

Figure S3. Comparison of the expression of seven genes encoding subtilisin-like protease between Me2XDR_vs._Me2XCK and Me4XDR_vs._Me4XCK detected by qPCR. Bars represent means ± SD of three biological replicates. Cassava β-actin was used as an internal control. (TIF 16466 kb)

Additional file 11:

Table S8. The 69 DE lncRNAs specific to 4× plants in response to drought stress. (XLSX 13 kb)

Additional file 12:

Table S9. The 138 DE lncRNAs common to both 2× and 4× leaves in drought-stressed plants. (XLSX 17 kb)

Additional file 13:

Table S10. The 104 DE lncRNAs specific to 2× plants in response to drought stress. (XLSX 15 kb)

Additional file 14:

Table S11. The 17 DE lncRNAs in 4× vs. 2× plants subjected to drought stress. (XLSX 9 kb)

Additional file 15:

Table S12. Sequences of the 69 DE lncRNAs specific to 4× plants in response to drought stress. (XLSX 51 kb)

Additional file 16:

Table S13. Sequences of the 17 DE lncRNAs in 4× vs. 2× stressed leaves. (XLSX 21 kb)

Additional file 17:

Table S14. A list of trans target genes of the 69 DE lncRNAs specific to 4× plants in response to drought stress, and Pearson correlations between lncRNAs and their corresponding target genes. (XLSX 411 kb)

Additional file 18:

Table S15. A list of trans target genes of the 17 DE lncRNAs in 4× vs. 2× stressed leaves, and Pearson correlations between lncRNAs and their corresponding target genes. (XLSX 193 kb)

Additional file 19:

Table S16. The transcript IDs co-expressed with 86 DE lncRNAs overlapping with DEGs in 4× leaves compared with 2× leaves under drought stress. (XLSX 47 kb)

Additional file 20:

Table S17. The top 20 KEGG enrichment pathways of the genes co-expressed with 86 DE lncRNAs overlapping with DEGs in 4× leaves compared with 2× leaves under drought stress. (XLSX 50 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Xiao, L., Shang, X., Cao, S. et al. Comparative physiology and transcriptome analysis allows for identification of lncRNAs imparting tolerance to drought stress in autotetraploid cassava. BMC Genomics 20, 514 (2019) doi:10.1186/s12864-019-5895-7

Download citation

Keywords

  • Autotetraploid
  • Cassava
  • Comparative transcriptomics
  • Drought stress
  • lncRNAs
  • Stomatal density