Uncovering the transcriptional responses of tobacco (Nicotiana tabacum L.) roots to Ralstonia solanacearum infection: a comparative study of resistant and susceptible cultivars

Background Tobacco bacterial wilt (TBW) caused by Ralstonia solanacearum is the most serious soil-borne disease of tobacco that significantly reduces crop yield. However, the limited availability of resistance in tobacco hinders breeding efforts for this disease. Results In this study, we conducted hydroponic experiments for the root expression profiles of D101 (resistant) and Honghuadajinyuan (susceptible) cultivars in response to BW infection at 0 h, 6 h, 1 d, 3 d, and 7d to explore the defense mechanisms of BW resistance in tobacco. As a result, 20,711 and 16,663 (total: 23,568) differentially expressed genes (DEGs) were identified in the resistant and susceptible cultivars, respectively. In brief, at 6 h, 1 d, 3 d, and 7 d, the resistant cultivar showed upregulation of 1553, 1124, 2583, and 7512 genes, while the susceptible cultivar showed downregulation of 1213, 1295, 813, and 7735 genes. Similarly, across these time points, the resistant cultivar had downregulation of 1034, 749, 1686, and 11,086 genes, whereas the susceptible cultivar had upregulation of 1953, 1790, 2334, and 6380 genes. The resistant cultivar had more up-regulated genes at 3 d and 7 d than the susceptible cultivar, indicating that the resistant cultivar has a more robust defense response against the pathogen. The GO and KEGG enrichment analysis showed that these genes are involved in responses to oxidative stress, plant–pathogen interactions, cell walls, glutathione and phenylalanine metabolism, and plant hormone signal transduction. Among the DEGs, 239 potential candidate genes were detected, including 49 phenylpropane/flavonoids pathway-associated, 45 glutathione metabolic pathway-associated, 47 WRKY, 48 ERFs, eight ARFs, 26 pathogenesis-related genes (PRs), and 14 short-chain dehydrogenase/reductase genes. In addition, two highly expressed novel genes (MSTRG.61386-R1B-17 and MSTRG.61568) encoding nucleotide-binding site leucine-rich repeat (NBS-LRR) proteins were identified in both cultivars at 7 d. Conclusions This study revealed significant enrichment of DEGs in GO and KEGG terms linked to glutathione, flavonoids, and phenylpropane pathways, indicating the potential role of glutathione and flavonoids in early BW resistance in tobacco roots. These findings offer fundamental insight for further exploration of the genetic architecture and molecular mechanisms of BW resistance in tobacco and solanaceous plants at the molecular level. Supplementary Information The online version contains supplementary material available at 10.1186/s12870-023-04633-w.

Moreover, transcriptomics, metabolomics, genomics, and proteomics are practical and diverse methodologies for elucidating and comprehensively understanding complex biological mechanisms under various stressful conditions in plants [17,[29][30][31].Of these, transcriptomics analysis with RNA-seq or microarrays has been used to unravel the molecular basis and genes associated with the specific biological processes in solanaceous crops, which could provide BW resistance for developing superior cultivars [9,30].In light of previous studies, transcriptomics and proteomics analysis showed that the methionine cycle (MTC), phenylpropanoid biosynthesis, glutathione metabolism, and gamma-aminobutyric acid biosynthesis pathways played essential roles in tomato BW [9,29,32].[32] used resistant and susceptible tomato genotypes to demonstrate that 140 up-regulated DEGs were related to hormones, lignin, and pathogenesis in the resistant genotype, while no changes occurred in the susceptible genotype.Similarly, Li et al. [8] identified 302 DEGs associated with potato BW, of which 81 were considered for BW resistance involved in signal transduction, terpenoids, pathogen recognition, hypersensitive response, and protection [8].
For TBW, Gao et al. [9] identified 158 and 835 DEGs in resistant and susceptible accessions, respectively, in a seedling root transcriptome study, while Shi et al. [29] reported that indole-3-acetic acid (IAA) and abscisic acid (ABA) plays significant roles in BW resistance in transcriptomics and metabolomics studies.A study reported the up-regulation of 6,233 genes in tobacco-resistant cultivars in response to BW infection.These genes were enriched in cell walls, ABC transporters, endocytosis, and glutathione metabolism [33].Similarly, Li et al. [17] and Pan et al. [33] used the tobacco stem transcriptome to detect genes, identifying the phenylpropanoid metabolic pathway for tobacco defense against BW.Furthermore, many studies have reported the root transcriptome in pepper [34], potato [35], eggplant [36], and tomato [37] as compared to tobacco.For example, root transcriptomes at 1, 3, and 5 d identified 115 genes in response to pepper BW infection [34].Moreover, transcription factors (TFs) [38], including WRKY, ERFs, NAC, MYB, bHLH, and bZip, and gene families like P450, MAPK, and DRGs identified in Solanaceae crops are involved in pathogen resistance [17,29,39].Further, recent studies have reported candidate genes associated with plant defense against BW, including StMKK1 in potato [40], CaNAC2c in pepper [41], SlNAP1 in tomato [42], and NbPDKs in tobacco [43].Limited transcriptomic studies have been conducted to investigate the genes involved in tobacco's resistance to R. solanacearum.Therefore, exploring and analyzing defense-related genes in tobacco through transcriptomics can provide valuable insights into the interaction between tobacco and R. solanacearum.
Controlling tobacco bacterial wilt in China is difficult due to its high severity and limited availability of resistant tobacco cultivars.Recently, we identified 52 candidate genes and 38 quantitative trait nucleotides (QTNs) using association mapping of 94 tobacco accessions with 126,602 SNPs [7].Findings from the studies mentioned above in other solanaceous species can serve as a reference for advancing tobacco breeding efforts.However, the genetic basis and molecular mechanism of BW resistance in tobacco may differ from those in other solanaceous species.Until now, few transcriptome profiling studies of BW resistance have been reported in tobacco.Therefore, in this study, we used RNA-seq to analyze and compare transcriptomes in the roots of D101 [resistant (R)] and Honghuadajinyuan [susceptible (S)] cultivars at early (6 h, 1 d) and late (3 d and 7 d) stage of seedlings after inoculation.The study aimed to understand the early and late response of gene expression to R. solanacearum infection in roots and to identify the genes conferring root resistance to R. solanacearum through gene ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG), and functional enrichment analysis.We identified significant divergence between the resistant and susceptible cultivars in their pathogen responses and also identified potential genes associated with BW resistance in tobacco plants for gene cloning, molecular studies, and breeding purposes.

Phenotypic divergence of tobacco cultivars under R. solanacearum infection
Firstly, we performed disease assays to determine the response of tobacco cultivars (D101 and Honghuadajinyuan) to bacterial wilt.We observed a significant difference between accessions in disease symptoms at 3 d and 7 d.Notably, the TBW symptoms (leaf wilting) started to appear in Honghuadajinyua after 3 d exposure to infection, with no symptoms in D101 (Fig. 1A).At 7 d, Honghuadajinyua stems looked black, with some plant death (Fig. 1B-C), whereas D101 seedlings had no similar or obvious symptoms (Fig. 1A).The roots of R and S cultivars showed clear disease symptoms at 7 d compared to control condition (Additional file 2: Fig. S1).These findings suggest that the D101 cultivar is highly resistant and Honghuadajinyua is susceptible to TBW infection (Fig. 1B).

Overview of RNA sequencing
A total of 1,458.74 million raw reads were collected with an average of 48.62; each library had 40.81-58.89 million raw reads (Additional file 1: Table S1).A total of 1,453.41(99.64%) million clean reads ranging from 40.66 (99.58%) to 58.67 (99.68%) were obtained after removing adopters, low-quality reads, and unknown bases (Additional file 1: Table S1).The clean reads data was high quality, with a Q30 base rate ranging from 92.78-94.94%,with an average GC content of 42.90%.Further, clean reads were aligned with the reference genome of tobacco (Nitab 4.5), with ~ 1,314.89(90.80%) million total reads mapped to the genome (Additional file 1: Table S1), of which 1,235.21 (85.29%) million reads were uniquely mapped to transcripts, with an average 41.17 (85.29%) and range 33.62-52.03(68.44-90.39%),indicating that the sequenced data were suitable for further analysis.Finally, 81,534 genes were identified, of which 69,500 were known, and 12,034 were new.The PC analysis showed 81.70% overall variation, of which PC1 and PC2 had 67.60% and 14.10%, respectively, variation in expression (Fig. 2A).The hierarchical clustering analysis showed similarities between the biological replicates of each R and S cultivar sample (Fig. 2B).Through PC and clustering analysis, the resistant and susceptible accessions clustered from each other, with their replicates located nearby (Fig. 2A-B).The R-7d and S-7d after infection were clustered separately from other time points.All the time points of R and S cultivars except R-7d and S-7d had significant positive correlations (r = 0.60-0.70)with other time points (Fig. 2C).The boxplots showed differences between different time points (Fig. 2D).

Mining differentially expressed genes (DEGs) under bacterial wilt
A pairwise comparison of gene expression levels between inoculation and control treatments identified 23,568 (19,130 known and 4,438 novels) DEGs in R and S at four time points (Fig. 3A and Additional file 1: Table S2).In the R cultivar, 20,711 (87.88%) genes were differentially expressed, of which 2,587 (1,553 up and 1,034 down), 1,873 (1,124 up and 749 down), 4,269 (2,583 up and 1,686 down), and 18,598 (7,512 up and 11,086 down) DEGs were identified at R-6 h, R-1d, R-3d, and R-7d, respectively, and 277 DEGs were common at all stages (Fig. 3B-C and Additional file 1: Table S2).Likewise, S samples at S-6 h, S-1d, S-3d, and S-7d had 3,166 (1,953 up and 1,213 down), 3,085 (1,790 up and 1,295 down), 3,147 (2,334 up and 813 down), and 14,115 (6,380 up and 7,735 down) DEGs, respectively (Fig. 3B), with 554 commonly identified at all time-points (Fig. 3D and Additional file 1: Table S2).In addition, 157 DEGs were commonly induced at 6 h, 1 d, 3 d, and 7 d in R and S cultivars (Fig. 3B-C).The R cultivar had 1.10 and 1.17 times more up-regulated genes than the S cultivar at 3 d and 7 d, respectively, while the S cultivar had 1.26 and 1.59 times more up-regulated genes than the R cultivar at 6 h and 1 d (Additional file 1: Table S2).These results indicate that TBW infection had a more substantial regulatory effect on gene regulation at 3 d and 7 d in R cultivars and 6 h and 1 d in S cultivars.

KEGG enrichment analysis of DEGs
Compared with the control, the DEGs at 6 h, 1 d, 3 d, and 7 d were significantly assigned to 11, 13, 17, and 35 pathways in the R cultivar (Fig. 5A-B and Additional file 1: Table S4  respectively.In the R cultivar, the highest number of DEGs were enriched in phenylpropanoid biosynthesis, pyruvate metabolism, plant-pathogen interaction, phenylalanine metabolism, carbon metabolism, plant hormone signal transduction, glutathione metabolism, MAPK signaling pathway, steroid biosynthesis, peroxisome, and circadian rhythm, involved in response to R. solanacearum at all time-points (Fig. 5A-B and Additional file 1: Table S4).Similarly, biosynthesis of various antibiotics, ribosomes, biosynthesis of amino acids, glycolysis/gluconeogenesis, plant hormone signal transduction, glutathione metabolism, phenylpropanoid biosynthesis, photosynthesis, plant-pathogen interaction, and phenylalanine metabolism were involved in the response of the S cultivar to BW infection at all time-points (Additional file 2: Fig. S8-S9 and Additional file 1: Table S4).KEGG analysis revealed the enrichment of several pathways, including MAPK signaling, phenylpropanoid biosynthesis, plant hormone signal transduction, and phenylalanine and glutathione metabolism in both R and S cultivars (Fig. 5).Additionally, a substantial number of highly up-regulated genes in these cultivars were associated with phenylalanine metabolism, which serves as a precursor for various secondary metabolites in plants.

Identification of potential candidate genes associated with bacterial wilt resistance
Based on GO and KEGG enrichment and functional annotation analysis, we evaluated potential candidate genes related to R. solanacearum resistance (Figs. 6 and 7 and Additional file 1: Tables S5-S7), including those involved in phenylalanine metabolism (stilbenoid, diarylheptanoid, and gingerol biosynthesis//ko00945 and flavonoid biosynthesis//ko00941) and glutathione metabolism pathways, and TFs (WRKYs, ARFs, and ERFs) and gene families (SDR, NBS-LR, PRs, and CKX).Hence, DEGs within these pathways or TFs were selected as potential candidate genes for subsequent investigations.

Gene families involve in TBW resistance
NBS-LRR genes encode nucleotide-binding site leucinerich repeat (NBS-LRR) proteins involved in pathogen recognition and signaling.We identified two novel genes (MSTRG.61386-R1B-17and MSTRG.61568)up-regulated in both cultivars at 7 d (Additional file 1: Table S7).Further, 26 genes that encode pathogenesis-related (PRs) proteins were identified in both cultivars, which may be involved in various aspects of the plant defense response (Fig. 7D and Additional file 1: Table S7).Among these, Nitab4.5_0014031g0010(PRB1), Nitab4.5_0004861g0030(PRB1), Nitab4.5_0008835g0020(STH-2), and Nitab4.5_0000360g0100 (TLP1) were differentially expressed in the R cultivar than the S cultivar at 1 d, 3 d, and 7 d (Fig. 7D).Additionally, PRB5 genes were down-regulated in both cultivars.In this study, 14 DEGs encode shortchain dehydrogenase/reductase-detected at different time points in the R and S cultivars after bacterial infection-involved in the biosynthesis of phytohormones that play a key function in defense response of plants against pathogens and BW (Fig. 7E and Additional file 1: Table S7).Nitab4.5_0011992g0020 and Nitab4.5_0000001g0060were significantly up-regulated in both cultivars at all time points.One gene, Nitab4.5_0001703g0050(SDR2a), was highly expressed in the R cultivar but low in the S cultivar (Fig. 7E).The remaining SDR genes showed an initial upregulation during the early stage of infection, followed by a down-regulation after post-infection.

Validation of transcriptomic data through quantitative real-time PCR
To measure the expression levels of potential resistance-related genes, 13 DEGs (11 up-regulated and 2 down-regulated) were randomly selected and analyzed by qRT-PCR at 6 h, 1 d, 3 d, and 7 d after infection.Of these, two genes (Nitab4.5_0000356g0010 and Nitab4.5_0012932g0010)associated with GST, one gene (Nitab4.5_0005779g0060)related to CKX5, and one gene (Nitab4.5_0008666g0010)related to short-chain dehydrogenase TIC32.As a result, the RT-qPCR produced similar results as the RNA-seq data at 6 h, 1 d, 3 d, and 7 d (Fig. 8A-H).The correlation analysis revealed a highly significant positive correlation between RT-qPCR and RNA-seq data at R-6 h (r = 0.95 and P = 8.5 × 10 -5 ) and S-6 h (r = 0.96 and P = 9.8 × 10 -5 ), indicating methodological reliability and supports the consistency of the RNAsequencing data (Fig. 8A-H).

Discussion
Tobacco bacterial wilt caused by R. solanacearum is a major soil-borne disease that adversely affects the quality and yield of tobacco crops.The survival of R. solanacearum in non-host or soil, and its ability to cause severe damage once the host appears makes TBW a significant economic threat [7,9,17,44,45].The plant's transcriptional response to pathogen infection is a dynamic process involving the modulation of signals and compounds that change over time [29,33,46].This study identified 239 potential candidate genes associated with the phenylpropane pathway, glutathione metabolism, WRKY, ERFs, ARFs, PRs, and NBS-LRR using root transcriptomic sequencing data of R and S cultivars at 6 h, 1 d, 3 d, and 7 d (Fig. 6-7 and Additional file 1: Tables S5-S7).
The Illumina sequencing data was reliable and accurately verified using RT-qPCR expression (Fig. 8) with a high Q30 base rate ranging from 92.78-94.94%.A total of 81,534 genes were identified, with 69,500 known and 12,034 novel genes, in line with published studies [33,46].The results indicate that BW infection significantly changes gene expression in R and S cultivars, with 20,711 DEGs in the R cultivar and 16,663 DEGs in the S cultivar (Fig. 3A-B and Additional file 1: Table S2).Among these, 554 were common between all samples (Fig. 3B-C), suggesting that the response to BW infection is cultivarspecific and occurs at different time-points in the two cultivars, similar to studies on tobacco [9], tomato [37], pepper [34], and Arabidopsis [47].For example, Ishihara et al. [32] identified 140 DEGs in a resistant tomato variety and none in a susceptible variety.
The identified DEGs were involved in 1,236 and 1,121 GO terms in the R and S cultivars, respectively (Fig. 4 and Additional file 1: Table S3).The R cultivar had a stronger defense response to R. solanacearum infection than the S cultivar, as evidenced by the up-regulation of genes involved in transferase activity, cell wall organization, hormone signal transduction, abiotic stress response, and oxidative stress response.In contrast, the S cultivar had a weaker defense response, as evidenced by the up-regulation of genes involved in tetrapyrrole binding, oxidoreductase activity, and carbohydrate metabolic processes (Fig. 4 and Additional file 1: Table S3).Our results are consistent with previous studies on BW resistance in solanaceous crops, linking resistance to increased expression of genes related to defense mechanisms, cell wall organization, and oxidative stress response [29,[32][33][34]48].For example, Pan et al. [33] found that 3,967 DEGs were associated with endocytosis, cell wall, signal transduction, and starch/sucrose metabolism in response to BW infection.
The KEGG analysis showed that the DEGs at 6 h, 1 d, 3 d, and 7 d were significantly assigned to different pathways in the R and S cultivars (Fig. 5 and Additional file 1: Table S4).The highest up-regulated genes in the R cultivar (Additional file 2: Fig. S6-S7 and Additional file 1: Table S4) were enriched in plant hormone signal transduction, glutathione metabolism, Fig. 8 Validation of 13 DEGs in response to bacterial wilt resistance via RT-qPCR.The RT-qPCR and RNA-seq log 2 FC values were compared using correlation analysis at R-6 h (A), R-1d (B), R-3d (C), R-7d (D), S-6 h (E), S-1d (F), S-3d (G), and S-7d (H).R and P indicate the correlation coefficient and corresponding P-value.The correlation showed a significant positive relationship between the RT-qPCR and RNA-seq at each time point in R and S cultivars phenylpropanoid biosynthesis, plant-pathogen interactions, and MAPK signaling pathway, consistent with studies investigating the molecular mechanisms of BW resistance in tomato [37], eggplant [36], tobacco [17,30,33], and potato [8,35].Furthermore, other studies reported the involvement of the enriched pathways in tobacco defense against R. solanacearum [26,29,46].These results also align with studies that reported the activation of multiple defense mechanisms in response to TBW infection in resistant tomato and eggplant cultivars, including the activation of glutathione metabolism and phenylpropanoid biosynthesis [32,34,49].
The present study identified 49 potential candidate genes in the phenylpropanoid pathway (Fig. 6A).Activation of the phenylpropanoid pathway has been associated with increased resistance to BW in tobacco [9,35,46].Six PAL genes (Nitab4.5_0002130g0010,Nitab4.5_0000754g0260,Nitab4.5_0012382g0010,Nitab4.5_0018680g0010,Nitab4.5_0000582g0180, and Nitab4.5_0000071g0170)were up-regulated at 3 d and 7 d post-infection (Fig. 6A).Several studies have investigated the role of the phenylpropanoid pathway in defending against BW [17,35,37].For example, Park et al. [50] reported that the up-regulation of PAL genes plays a significant role in the plant defense response to BW, while the OsPAL2 mutant gene in rice contributed to the control of pathogen disease resistance [51].This study revealed two genes (CYP98A2 and CYP98A3) of cytochrome CYP450 (Fig. 6A) that play an essential role in BW resistance in many crops [52][53][54][55].Moreover, we identified 46 candidate genes involved in glutathione metabolism, including APX, GST, ODC, OXP, and PRP (Fig. 6B), which play a crucial role in the detoxification of ROS and the regulation of cellular signaling pathways that lead to the activation of defense mechanisms [6,30,37,56,57].Our results showed that the R cultivar up-regulated more genes involved in glutathione metabolism than the S cultivar at 3 d and 7 d after TBW inoculation (Fig. 6B), consistent with other studies [7,36,58].Overexpression of GST and GPX genes increased expression activity twofold in wild-type tobacco seedlings, contributing to chilling or salt stress tolerance [59].Li et al. [17] reported up-regulated DEGs involved in the phenylpropanoid pathway and glutathione metabolism for TBW resistance.Therefore, the glutathione metabolism and phenylpropanoid pathway are essential defense mechanisms in tobacco plants against BW caused by R. solanacearum.Further research is needed to understand how these pathways contribute to disease resistance in tobacco.
TFs are crucial plant regulators [38,60], with WRKY TFs playing essential roles in plant immunity against various stresses.WRKY40, WRKY6, WRKY27, and WRKY22 have been identified as positive regulators of Solanaceae crop resistance to BW in several [38,[60][61][62].In this study, WRKY6 (Nitab4.5_0004702g0020,Nitab4.5_0005781g0050, and Nitab4.5_0000005g0190)and WRKY11 (Nitab4.5_0000286g0080)were found to exhibit a resistance response to TBW in the R cultivar (Fig. 7A), particularly in the later stage, consistent with the response pattern observed in pepper to BW [62].These findings indicate that these genes may be important regulators of tobacco resistance against BW [63,64].Moreover, ERFs are a class of plant TFs involved in biological and abiotic stress responses [17,65,66].In this study, 48 ERFs genes were significantly expressed in R and S cultivars at all time points (Fig. 7B).The expression pattern of ERF5 (Nitab4.5_0000734g0070) in the R cultivar was distinct from that in the S cultivar (Fig. 7B), especially at 3 d.The eight ARF genes identified in the present study were up-regulated at 3 d and 7 d after infection (Fig. C); these genes play a critical role in regulating plant growth and development to improve plant defense mechanisms against pathogens [17,33,39].Pathogenesis-related proteins (PRs) are essential for BW resistance [67,68], with the present study identifying 26 PR genes (PR-1, PR-4, and PR-5) in both cultivars at all time points (Fig. 7D).Two PR-1 genes (Nitab4.5_0004861g0040 and Nitab4.5_0005400g0020)were up-regulated at 3 d and 7 d in both cultivars (Fig. 7D), similar to previous studies [68,69].Rivière et al. [69] found that PR-1b was up-regulated in tobacco plants resistant to BW and that silencing PR-1b increased disease susceptibility.Similarly, Alamillo et al. [70] and Pruss et al. [71]reported up-regulated expression of PR-2 in tobacco infected with R. solanacearum, while silencing PR-2 increased disease susceptibility [72].Further, we identified 14 SDR genes involved in various metabolic pathways, including detoxification and hormone biosynthesis (Fig. 7E).In tobacco, SDR genes play a role in the defense response against bacterial pathogens: Yu et al. [73] reported that SDR genes are involved in the biosynthesis of jasmonic acid, which is critical for defense response against pathogens.Overall, the present study identified 239 potential candidate genes associated with pathways, TFs, and gene families (Figs. 6 and 7 and Additional file 1: Tables S5-S7), which will help understand the defense response mechanism and could be potential targets for improving resistance to BW disease.In summary, our results revealed that DEGs associated with phenylpropane/flavonoids pathway, glutathione metabolic pathway, WRKY, ERFs, ARFs, pathogenesis-related genes (PRs), and short-chain dehydrogenase/reductase genes essential to control the resistance and programmed cell death.Of these, most of the defense-related genes were up-regulated in the R cultivar compared to the S cultivar at an early stage of infection.Therefore, our hypothesized model indicates the involvement of defense-related and susceptible-related genes networking of R. solanacearum in the tobacco system (Fig. 9).

Conclusion
This study used high throughput sequencing to identify DEGs and related pathways in resistant and susceptible cultivars under R. solanacearum infection at early and late seedling stages.In early plant response after 6 h, 2 ARFs, 10 WRKY, and 11 ERFs genes were differentially expressed, which might be the plants' primary response.The R and S cultivars had similar response levels at the early stage, but the R cultivar-related DEGs had specific expression patterns in response to infection after 7 d.Finally, 239 potential candidate genes were identified using GO, KEGG, functional enrichment analysis, and literature search.These DEGs in the resistant cultivar are related to glutathione metabolism and phenylpropanoid biosynthesis, TFs, and PRs, and are likely the main substances conferring root resistance against BW infection in tobacco.Additionally, two novel genes (MSTRG.61386-R1B-17and MSTRG.61568) were upregulated in both cultivars at 7 d, which can be used for functional studies to find the exact role in BW resistance.The finding of this study shows the molecular mechanisms involved in tobacco resistance to bacterial wilt and provides an important source for controlling and breeding tobacco plants resistant to TBW.

Plant material and growth conditions
Two tobacco cultivars, D101 (R) and Honghuadajinyuan (S), were received from Nanxiong Scientific Research Institute of Guangdong Research Institute of Tobacco Science, Shaoguan, China.The experiments were conducted in the greenhouse of Guangzhou University, China.Seeds of both cultivars were surface-sterilized with distilled water and spread on a black plastic tray with white and pit seedling sponges.The tray was covered with a layer of fresh-keeping film and placed in the growth chamber to germinate at 28℃, 60% relative humidity, and 16 h/8 h light/dark photoperiod.After germination, seedlings in the crosswise stage were transferred into hydroponic boxes (41 × 24 × 14 cm; length × width × height).Hoagland's nutrient solution was added to the aerated hydroponic containers, comprising 945 mg/L Ca(NO 3 ) 2

Inoculation treatment and sample collection
The highly pathogenic R. solanacearum strain HBLC5 was isolated from tobacco plants obtained from Hukou experimental station in Nanxiong City as described in our previous studies [7,15].For inoculation preparation, Fig. 9 A model pathway involved in resistance against tobacco bacterial wilt.The tobacco reaction to R. solanacearum infection and the systematic pathway involved in resistance.The signals from the pathogens were recognized through related receptors (R genes) and transcription factors, which triggered the plant defense mechanism.After that, defense-related molecules/genes, including PR-proteins, metabolites, glutathione and phenylpropane genes, were activated to cause resistance to pathogens through high expression levels the R. solanacearum strain HBLC5 was cultured on BG medium and reproduced in Liquid LB medium by shaking at 30℃ and 180 rpm for two days.Five-leaf stage plants were infected with 300 mL cell suspension (adjusted to 1 × 10 8 CFU/mL and OD 600nm = 0.1) by irrigating roots in one pot.The control was treated with sterile water instead of the bacterial solution.Each tobacco variety was grown in three pots, each with 60 seedlings.Tobacco roots from three individual seedlings taken from one pot were considered one biological replicate.Therefore, three biological replicates of each variety were harvested from three pots at every time point.After inoculation, the roots from R and S cultivars were collected at R/S-0 h, R/S-6 h, R/S-1d, R/S-3d, and R/S-7d, and the samples were stored at -80 °C after freezing in liquid nitrogen until RNA extraction.

RNA extraction, library preparation, and Illumina sequencing
The total RNA of tobacco seedlings roots at five-six leaf stage (~ 30 days) was extracted from three biological replicates at each time-point using the TRIzol Reagent Kit (Invitrogen, Carlsbad, CA, USA) with the manufacturer's instructions.In short, RNA was extracted from the frozen samples, and DNase I (TakaraBio, Japan) was used to remove DNA contamination.The RNA quality and integrity were assessed through 1% agarose gel electrophoresis and a microplate spectrophotometer (BioTek Company, USA).To develop the cDNA libraries, mRNA with polyA tails was purified using Oligo(dT) magnetic beads and fragmented.The first and second strands of cDNA were synthesized with random hexamer primers and DNA polymerase I, respectively.Adapters were ligated to double-stranded cDNAs, which were then exposed to PCR amplification and purification.The quality of the cDNA libraries was measured using an Agilent 2100 Bioanalyzer, and sequencing was performed on the Illumina Novaseq 6000 platform by Genedenovo Biotechnology Co., Ltd in Guangzhou, China.

Sequencing data analysis and genomic annotation
The sequenced raw reads were processed using fastp v0.18.0 software [74] to determine the reliability and quality of sequenced data and remove the reads that contained adaptors, ≤ 10% unknown bases, and low-quality reads at Q ≤ 20%.Clean reads were used to determine Q20, Q30, and GC contents and for all further downstream analyses.The Nitab4.5 reference genome of N. tabacum [75] and gene annotation files were downloaded from the solgenomics database (https:// solge nomics.net/ organ ism/ Nicot iana_ tabac um/ genome).The first step was to align the clean reads (i.e., reads filtered for quality and adapter sequences) to a reference genome.In this case, the clean reads were aligned to Nitab4.5 using HISAT2 v2.1.0[76] and Bowtie2 v2.2.5 [77].After alignment, the sequences of all assembled genes were aligned to different databases for functional annotation, including the Nitab4.5 reference genome, UniProt, NT, KEGG, GO, and COG using Trinotate (http:// trino tate.github.io/) software with default parameters except E ≤ 1e-5.Then, transcripts were reconstructed using StringTie v2.2.0 [78] based on the HISAT2 alignment.Transcript abundance was measured as FPKM values using RSEM software [79].Finally, the principal component analysis (PCA) and Pearson's correlation coefficients were calculated using the R4.1.0environment (http:// www.r-proje ct.org/) based on the FPKM expression data.

Identification of differentially expressed genes
DESeq2 R package [80] was used to perform differential expression analysis by pairwise comparison of the control and treatments at 6 h, 1 d, 3 d, and 7 d.The P-value was adjusted using the Benjamini-Hochberg method with an FDR threshold of 0.05, and genes were identified as upor down-regulated based on the |log2FC|> 1 criteria.The hierarchical K-mean clustering method was used to analyze the DEGs using the hclust function in R version 4.1.0(http:// www.r-proje ct.org/), and DEGs were then visualized using a heatmap as log 2 FC input.

GO and KEGG enrichment analysis
GO enrichment analysis of DEGs was performed using the GO database (http:// www.geneo ntolo gy.org/) with a significance threshold of corrected P ≤ 0.05.KEGG enrichment analysis was conducted using the online KEGG database (https:// www.genome.jp/ kegg/ kegg1.html) at P ≤ 0.05 as a significance threshold to identify significant biological functions and metabolic pathways [81].The clusterProfiler v3.4.4 software was used to visualize the top GO terms and KEGG pathways [82].

Gene expression analysis through RT-qPCR
Thirteen genes were randomly selected with different expression profiles to validate sequencing data by RT-qPCR using specific primers presented in Additional file 1: Table S8.Reverse transcription was completed using the HiScript ® III RT SuperMix for PCR (+ gDNA wiper) Kit (Vazyme Biotech Co., Ltd).The qRT-PCR was performed using 1 µL cDNA template mixed with 10 µL 2 × ChamQ SYBR qPCR Master Mix (Vazyme Biotech Co., Ltd) and 0.4 µL each of the forward and reverse primers in a final volume of 20 µL.Amplification consisted of 40 cycles of 95℃ for 10 s, 60℃ for 30 s, and 95℃ for 15 s, followed by a denaturation step to generate the melting curves.The Ntubc2 (accession number: AB026056) was selected as an internal reference gene, and gene expression profiles were measured using the 2-∆∆CT method.Three technical and biological replicates were used for each sample of cultivars at all time points.
• fast, convenient online submission • thorough peer review by experienced researchers in your field • rapid publication on acceptance • support for research data, including large and complex data types • gold Open Access which fosters wider collaboration and increased citations maximum visibility for your research: over 100M website views per year

•
At BMC, research is always in progress.

Learn more biomedcentral.com/submissions
Ready to submit your research Ready to submit your research ?Choose BMC and benefit from: ? Choose BMC and benefit from:

Fig. 1
Fig. 1 Disease symptoms in D101 (resistant; R) and Honghuadajinyuan (susceptible; S) tobacco cultivars after inoculation with R. solanacearum.Seedlings of R (A) and S (B) cultivars at 0 d, 1 d, 3 d, and 7 d after infection.The leaves show clear disease symptoms at 3 and 7 days in both cultivars as compared to 1 and 3 days, but the leaves of the S cultivar at 7 d show black color and plant death.C Bacterial wilt symptoms in leaves and stem of tobacco seedlings ) and 10, 14, 19, and 40 pathways in the S cultivar (Fig. 5C-D and Additional file 2: Fig. S8-S9),

Fig. 2
Fig. 2 Overview of tobacco root transcriptomic data based on bacterial wilt resistant (R) and susceptible (S) tobacco cultivars at five time points.A Principal component analysis represents 81.70% overall variation in expression; B Cluster analysis using FPKM values of each sample with replicates indicates the similarities between the biological replicates of both cultivars; C Correlation analysis between resistant and susceptible tobacco cultivars at four time points between the samples; D log 10 transformation of FPKM values at four time-points show the differences between the samples at each time points

Fig. 3
Fig. 3 Analysis of differentially expressed genes in resistant (R) and susceptible (S) tobacco cultivars in response to bacterial wilt infection.A Volcano plots represent significant genes in R and S cultivars at 6 h, 1 d, 3 d, and 7 d relative to CK-0 h, with log 2 FC values transformed into -log 10 (FDR).Blue and red dots indicate up-regulated and down-regulated DEGs (log2FC ≥ ± 1, FDR ≤ 0.05); black dots represent non-significant.B The number of genes differentially expressed at different time points.The highest number of DEGs were identified at 7 d in R and S cultivars.C-D The Venn diagram shows the number of unique and common DEGs in R and S tobacco cultivars at each time point

Fig. 4
Fig. 4 Gene ontology (GO) enrichment analysis of up-and down-regulated genes in bacterial wilt resistant (R) and susceptible (S) tobacco cultivars at 6 h and 1 d.The top significant GO terms of the three categories are shown at R-6 h (A), R-1d (B), S-6 h (C), and S-1d (D) at P ≤ 0.05.The circular diagram represents the enrichment, GO terms with the number of up-and down-regulated genes.The outer circular line indicates the GO ID of each term, the middle circular line indicates the q-value of enrichment with the total number of genes, and the inner circular lines indicate the up-regulated (purple) and down-regulated (blue) number of genes

Fig. 5
Fig. 5 KEGG analysis of differentially expressed genes in bacterial wilt resistant (R) and susceptible (S) tobacco cultivars at 6 h and 1 d.The top KEGG pathways are listed at R-6 h (A), R-1d (B), S-6 h (C), and S-1d (D) at P ≤ 0.05.The circular diagram represents the enrichment KEGG pathways with the number of up-and down-regulated genes.The outer circular line indicates the pathway ID of each pathway, middle circular line indicates the q-value of enrichment with total number of genes, and inner circular lines indicate the up-regulated (purple) and down-regulated (blue) number of genes

Fig. 6
Fig. 6 Heatmaps of differentially expressed genes enriched in glutathione and phenylpropane pathways involved in bacterial wilt resistance in tobacco.A Number of DEGs involved in the phenylpropane pathway; B Number of DEGs involved in glutathione pathway and showed different expression levels at each time point

Fig. 7
Fig. 7 Heatmaps of transcription factors and gene families involved in bacterial wilt resistance in tobacco.The small Fig represent the significant DEGs of (A) WRKYs, (B) ERFs, (C) ARFs, (D) PRs, and (E) SDRs in both cultivars at different time points after the inoculation of bacterial wilt

Table 1
List of up-regulated and down-regulated DEGs, GO, and their expression spectrum under bacterial wilt infection