Nội dung

Qin et al. BMC Plant Biology (2020) 20:30 https://doi.org/10.1186/s12870-020-2241-9 RESEARCH ARTICLE Open Access Transcriptome analysis reveals key genes involved in the regulation of nicotine biosynthesis at early time points after topping in tobacco (Nicotiana tabacum L.) Yan Qin1†, Shenglong Bai1†, Wenzheng Li2†, Ting Sun1, David W. Galbraith1,3, Zefeng Yang4, Yun Zhou1, Guiling Sun1* and Bingwu Wang2* Abstract Background: Nicotiana tabacum is an important economic crop. Topping, a common agricultural practice employed with flue-cured tobacco, is designed to increase leaf nicotine contents by increasing nicotine biosynthesis in roots. Many genes are found to be differentially expressed in response to topping, particularly genes involved in nicotine biosynthesis, but comprehensive analyses of early transcriptional responses induced by topping are not yet available. To develop a detailed understanding of the mechanisms regulating nicotine biosynthesis after topping, we have sequenced the transcriptomes of Nicotiana tabacum roots at seven time points following topping. Results: Differential expression analysis revealed that 4830 genes responded to topping across all time points. Amongst these, nine gene families involved in nicotine biosynthesis and two gene families involved in nicotine transport showed significant changes during the immediate 24 h period following topping. No obvious preference to the parental species was detected in the differentially expressed genes (DEGs). Significant changes in transcript levels of nine genes involved in nicotine biosynthesis and phytohormone signal transduction were validated by qRT-PCR assays. 549 genes encoding transcription factors (TFs), found to exhibit significant changes in gene expression after topping, formed 15 clusters based on similarities of their transcript level time-course profiles. 336 DEGs involved in phytohormone signal transduction, including genes functionally related to the phytohormones jasmonic acid, abscisic acid, auxin, ethylene, and gibberellin, were identified at the earliest time point after topping. Conclusions: Our research provides the first detailed analysis of the early transcriptional responses to topping in N. tabacum, and identifies excellent candidates for further detailed studies concerning the regulation of nicotine biosynthesis in tobacco roots. Keywords: Tobacco, Topping, Transcriptome analysis, Differentially expressed genes, Nicotine biosynthesis and regulation * Correspondence: sung@vip.henu.edu.cn; bwwang76@hotmail.com † Yan Qin, Shenglong Bai and Wenzheng Li contributed equally to this work. 1 Key Laboratory of Plant Stress Biology, State Key Laboratory of Crop Stress Adaptation and Improvement, State Key Laboratory of Cotton Biology, School of Life Sciences, Henan University, Kaifeng 475004, China 2 Tobacco Breeding Center, Yunnan Academy of Tobacco Agricultural Sciences, Kunming 650021, Yunnan, China Full list of author information is available at the end of the article © The Author(s). 2020 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. Qin et al. BMC Plant Biology (2020) 20:30 Background For tobacco (Nicotiana tabacum L.) plants, topping (defined as the removal of the flowering head and young leaves) is an essential cultivation practice. Topping switches the plant from a seed reproductive to a leaf vegetative phase, and this significantly increases leaf nicotine contents [1]. A number of studies, employing a variety of experimental techniques, have investigated tobacco responses to topping [2–4]. An up-regulation of nicotine biosynthesis, found to occur exclusively in roots, and particularly in growing root tips, is one of the typical responses of tobacco plants to topping [5]. Nicotine plays pivotal roles both in establishing the commercial quality of tobacco, and in defending plants against herbivores. Biosynthesis of nicotine, a secondary metabolite associated with the tobacco stress response, is reproducibly promoted by topping in tobacco roots [6]. Nicotine comprises two main nitrogen-containing rings, the pyrrolidine ring and the pyridine ring [7]. Biosynthesis of the pyrrolidine ring involves arginine decarboxylase (ADC) [8], ornithine decarboxylase (ODC) [9], S-adenosylmethionine decarboxylase (SAMDC), S-adenosyl-L-methionine synthetase (SAMS), putrescine N-methyltransferase (PMT) [10, 11], and Nmethylputrescine oxidase (MPO) [12, 13]. Biosynthesis of the pyridine ring starts with the nicotinic acid dinucleotide (NAD) biosynthetic pathway. Enzymes participating in the early metabolic conversion steps of this pathway include aspartate oxidase (AO), quinolinate synthase (QS), and quinolinic acid phosphoribosyl transferase (QPT) [14–16]. The A622 gene (encoding an isoflavone reductase-like protein) is responsible for nicotine ring coupling, and the BBL genes (encoding berberine bridge enzyme-like proteins) are involved in the subsequent oxidation step that leads to nicotine [17, 18]. A recent report [19] employed suppression subtractive hybridization (SSH) techniques to further examine the transcriptional responses of tobacco roots during the first 24 h after topping. Of the 129 high quality expressed sequence tags identified as representing DEGs, most were involved in stress/defense, in secondary metabolism, and in signaling/transcription [19]. The regulation of nicotine biosynthesis has long been considered a complex physiological response, and many TFs are directly or indirectly involved in its regulation [20, 21]. Further insights into the transcriptional regulation of the nicotine biosynthetic pathway have come from the analysis of two subtractive cDNA libraries of jasmonate-treated Nicotiana benthamiana roots, and through examination of the effects of virus-induced gene silencing (VIGS) technologies. Of the sixty-nine TFs, six (from three TF families) affect nicotine metabolism, with NbbHLH1 and NbbHLH2 (basic helixloop-helix) genes positively regulating the jasmonate activation of nicotine biosynthesis, as evidenced by overexpression [22]. Page 2 of 15 Although specific genes regulating nicotine synthesis after tobacco topping have been identified, a detailed description of the transcriptional regulatory network that responds to topping is not available. The situation is further complicated by the allotetrapoid status of N. tabacum, formed through the hybridization of N. sylvestris (S-subgenome) and N. tomentosiformis (T-subgenome), and how these two subgenomes respond to topping is unclear. In this study, we have sequenced tobacco root transcriptomes at seven different time points (0, 0.5, 1, 3, 5, 8 and 24 h) after topping. These time points were chosen to identify candidate genes associated with the regulation of nicotine biosynthesis at the earliest stages, as well as to allow discovery of upstream regulators of nicotine synthesis through clustering of the time-course profiles of TF gene expression, and to compare the responses of the two subgenomes to topping. This comprehensive approach to characterization of the transcriptional responses of tobacco, especially focusing on the early regulation of nicotine biosynthesis, should serve to advance genetic improvement in this crop. Results Transcriptome sequencing and quality assessment Total RNA of tobacco roots, isolated separately from 18 individual plants, was employed for RNA sequencing (RNA-Seq) library construction. The 18 RNA-Seq libraries were sequenced using the Illumina platform. After filtering out low quality sequences (quality scores < 25), 105 Gb of cleaned data was obtained, representing approximately 6 Gb per sample. The cleaned sequence GC content varied from 42.1 to 42.7% (Additional file 6: Table S1). The mapping rates for the cleaned sample reads aligned against the reference genome sequence ranged from 91.6 to 97.8% (Additional file 6: Table S1). The sequencing quality and gene expression levels were generally consistent across the sequenced samples (Additional file 1: Figure S1). Identification and verification of differentially expressed genes (DEGs) The expression levels of the genes from the tobacco transcriptomes were calculated and normalized to FPKM values (Fragments Per Kilobase of transcript per Million fragments mapped). The values of Pearson’s Correlation Coefficient across biology replicates exceeded 0.82. In terms of the correlation between samples from different time points, some samples displayed higher values with those from other time points. For example, BWR3-2A showed a correlation coefficient of 0.96 as compared to BWR24-1A, and 0.95 with BWR5-2A (Additional file 2: Figure S2). Further experiments will be required to elucidate this unexpected observation. Through comparison of the samples at each time point to the t = 0 sample, Qin et al. BMC Plant Biology (2020) 20:30 Page 3 of 15 and using a fold-change (FC) > 2, and a false discovery rate (FDR) < 0.05 as selection criteria, 4830 DEGs were identified after topping. An almost identical number (2082 and 2075 genes) came from the N. tomentosiformis and N. sylvestris genomes, respectively (Additional file 7: Table S2). Notably, the number of DEGs at 0.5 h (2,562) was much greater than those at any other time point, indicating that more genes respond to topping at earlier times. The number of DEGs dropped to its lowest level (815) at 1 h after topping (Fig. 1a). However, a second burst of differential gene expression was observed at t = 8 h (1,756), followed by a decrease at t = 24 h (Fig. 1a). The results imply that the N. tabacum root produces two discrete peaks of transcriptional activity, at 0.5 h and 8 h after topping. This result is consistent with the gene numbers identified as specifically induced at each of the six time points after topping, the largest number being 1186 at t = 0.5 h and the second highest number being 585 at t = 8 h after topping (Fig. 1b). To validate the transcriptional results obtained by RNA-Seq, we selected nine genes related to nicotine biosynthesis and phytohormone signal transduction, and examined their transcriptional responses by qRT-PCR. The expression trends of these genes analyzed by qRTPCR were consistent with the RNA-Seq analysis performed at the corresponding time points (Fig. 2). The changes of selected DEGs obtained by RNA-Seq analysis had good correlations with those obtained by qRT-PCR a (R2 = 0.674). These results confirm that the alterations in gene expression detected by RNA-Seq accurately reflect transcript differences at the different time points after topping. Functional classification and enrichment analysis of DEGs 4830 DEGs showing significant variation at the different time points after topping were selected for further analysis. Based on their relative expression levels, the DEGs were divided into different categories using hierarchical clustering, being distinguishable in terms of the temporal patterns of the transcriptional responses of the roots at the various time points after topping (Additional file 3: Figure S3). The predicted functions of the DEGs were then obtained from their GO (Gene Ontology) annotations, and using KEGG (Kyoto Encyclopedia of Genes and Genomes) pathway analysis. According to GO term annotation, the DEGs were distributed across 42 functional terms, as follows: 19 terms for biological process, 12 terms for molecular functions, and 11 terms for cellular component (Additional file 4: Figure S4). GO enrichment analyses were performed to classify the putative functions of DEGs in the comparisons of libraries prepared from the different time points (Fig. 3). The DEGs in GO enriched categories of biological process were mainly involved in response to oxidative stress (GO:0006979), the phenylpropanoid metabolic process (GO:0009698), the lignin metabolic process b up down 1hvs0h 92 0.5hvs0h 2000 232 46 1300 Gene Number 170 62 25 11 6 10 10 4 14 31 7 3hvs0h 8 16 15 1186 60 39 14 18 5 9 56 699 27 12 35 602 1000 352 23 15753 17 493 380 34 1057 435 408 24hvs0h 0 101 11 10 33 45 8 24 4 25 27 2 37 116 8 8 17 11 245 162 5hvs0h h hv s0 8hvs0h 24 8h vs 0h 0h vs 5h vs 0h 3h h vs 0 1h vs 0h 580 5h 0. 417 8 5 19 767 758 55 9 44 34 1262 251 10 Fig. 1 DEGs statistics at different stages after tobacco topping. a In the stacked bar-graphs, the up-regulated DEGs are located in the red regions, and the down-regulated DEGs in the blue regions. b Venn diagram of DEGs at different time points after topping in tobacco Qin et al. BMC Plant Biology (2020) 20:30 Page 4 of 15 Fig. 2 Validation of RNA-Seq data by qRT-PCR. a Expression levels of 9 randomly-selected DEGs of the nicotine anabolic pathway as measured by qRT-PCR (the columns) and the corresponding expression trends measured by RNA-Seq (the lines). The error bars represent SDs (n = 3). Asterisks represent significantly different transcript levels between the topping treatment and control plants at the indicated times. (t-test; *, P < 0.05; **, P < 0.01; ***, P < 0.001). b Correlation analysis of fold-change data between qRT-PCR and RNA-Seq. Scatterplots are generated from the log2 expression ratios of qRT-PCR analyses (x-axis) and from RNA-Seq analyses (y-axis). Each scatter point depicts a time point at which significant differences in gene expression levels were found. The equation of the linear regression relationship and the associated correlation coefficient (R2) are provided Qin et al. BMC Plant Biology (2020) 20:30 Page 5 of 15 Fig. 3 Gene Ontology (GO) term enrichment analysis. Significantly enriched GO terms were selected based on a FDR < 0.05. GO terms of the categories of Biological Processes, Cellular Components, and Molecular Functions are depicted in red, green, and blue, respectively (GO:0009808), and response to abiotic stimulus (GO: 0009628). The DEGs of GO enriched categories of cellular component were mainly involved in the apoplast (GO:0048046), the extracellular region (GO:0005576), the external encapsulating structure (GO:0030312), and the cell wall (GO:0005618). The DEGs of GO enriched categories of molecular function were mainly associated with peroxidase activity (GO:0004601), antioxidant activity (GO:0016209), and a series of transporter activities (GO: GO:0006857, GO:0008272, GO:0008509, and GO: 0008271) (Fig. 3). To further investigate the functions of differentially expressed transcripts in response to topping, we performed enrichment analyses by mapping the sequences to the KEGG database categories. The DEGs with KEGG annotation were assigned to 28 classes, mainly related to signal transduction (221), carbohydrate metabolism (212), biosynthesis of other secondary metabolites (166), and metabolism of terpenoids and polyketides (69) (Additional file 5: Figure S5). KEGG enrichment analyses also indicated the DEGs were significantly enriched in the main pathways of phenylpropanoid biosynthesis (ko00940), of starch and sucrose metabolism (ko00500), and in the plant MAPK signaling pathway (ko04016, their responses to wounding and their roles in the biosynthesis of secondary metabolism have been illustrated previously [23–25] (Fig. 4). DEGs involved in nicotine synthesis and transport We further investigated whether the genes activated by topping were involved in nicotine biosynthesis and transport. As expected, nine gene families involved in nicotine biosynthesis (AO, QS, ODC, ADC, SAMS, PMT, A622, MPO, and BBL) (Fig. 5), and two gene families involved in nicotine transport (MATE, NUP), as identified by showing at least 93% identities with the primary sequences of previously reported enzymes, were found in the DEG dataset (Additional file 8: Table S3). All genes showed transcriptional up-regulation, with most being up-regulated at 8 h and 24 h after topping; our qPCR assay also verified the expression changes of four genes (PMT1a, PMT1b, MPO, ODC) at the corresponding time points (Fig. 2). Both gene families encoding MATE and NUP in nicotine transport were found to be upregulated (Additional file 8: Table S3). Similar to the situation across all DEGs, most of those involved in nicotine synthesis and transport were found in both subgenomes. One DEG encoding AO was derived from the T-subgenome and all the DEGs encoding MPO and ODC were from the S-subgenome. Transcription factors (TFs) of DEGs, and gene clustering by expression patterns To investigate the upstream regulatory mechanisms of nicotine biosynthesis after topping, we next focused on Qin et al. BMC Plant Biology (2020) 20:30 Page 6 of 15 Fig. 4 KEGG enrichment analysis. Each circle in the figure represents a KEGG metabolic pathway, and the number of genes enriched in a pathway corresponds to the size of the circle. The degree of significance of the enrichment of DEGs in a pathway is represented by -log10 (qvalue). The abscissa indicates the ratio of the number of DEGs annotated to a particular pathway to the number of the DEGs annotated to all pathways Fig. 5 The DEGs involved in nicotine biosynthesis. Solid and dashed lines indicate defined and undefined reactions, respectively, with the DEGs shown in red. Abbreviations: ADC, arginine decarboxylase; ODC, ornithine decarboxylase; SAMS, S-adenosyl-L-methionine synthetase; SAMDC, Sadenosylmethionine decarboxylase; PMT, putrescine N-methyltransferase; MPO, N-methylputrescine oxidase; AO, aspartate oxidase; QS, quinolinate, synthase; QPT, quinolinic acid phosphoribosyl transferase; A622, isoflavone reductase-like protein; BBLs, berberine bridge enzyme-like proteins. The genes with significantly up-regulated transcription levels are shown in red Qin et al. BMC Plant Biology (2020) 20:30 Page 7 of 15 the types of TFs represented in the DEGs from the tobacco root transcriptome. In our study, a total of 549 DEGs encoding TFs were identified (Additional file 9: Table S4), being divided into 49 TF families. Amongst these, the number of TFs was highest at t = 0.5 h (355), accounting for 65% of all TFs, with 240 being upregulated and 115 down-regulated. This was followed by the t = 8 h timepoint (253), accounting for 46% of all TFs, with 94 being up-regulated and 159 downregulated. This suggests that many TFs genes participate in instant early gene activation. Notably, and representing most of these TFs, 18 families were found to contain more than 10 gene members: AP2-EREBP (75), MYB (69), bHLH (44), NAC (30), bZIP (30), Orphans (26), HB (24), WRKY (23), HSF (19), C2C2-Dof (17), GRAS (16), LOB (15), MADS (12), GNAT (11),AUX / IAA (11), G2like (11), C3H (10) and C2H2 (10) (Fig. 6). To further examine the contributions of specific TFs to the regulatory network of nicotine biosynthesis, we performed clustering using the 549 TFs and the upregulated structural genes associated with the nicotine 0 20 40 60 AP2−EREBP MYB bHLH NAC all bZIP Orphans HB 0.5hvs0h 1hvs0h 3hvs0h 5hvs0h WRKY HSF 8hvs0h 24hvs0h C2C2−Dof GRAS LOB MADS GNAT AUX/IAA G2−like C3H C2H2 Fig. 6 TFs classification of DEGs in tobacco. DEGs at different time points after topping are represented by different colors, the abscissa representing the number of transcription factors at each time point biosynthesis pathway. Fifteen clusters showing similar expression profiles were obtained (Fig. 7 and Additional file 10: Table S5). It can be observed that several clusters are similar but with minor differences. For instance, the TFs of Clusters 2 and 12 were up-regulated at t = 0.5 h, and the TFs in Cluster 9 and 10 were upregulated at 0.5–1 h. They then returned to the expression levels found before topping (Fig. 7). Notably, most of the up-regulated DEGs in nicotine biosynthesis were in Cluster 11 (20 DEGs), which displayed the greatest upregulation at t = 8 h and at t = 24 h after topping (Fig. 7 and Additional file 10: Table S5). Seventeen genes from the bHLH family and the AP2-EREBP family were found in Cluster 11, including ERF189 (Nitab4.5_0003090g0030 and Nitab4.5_0015055g0010), and ERF91 (Nitab4.5_ 0004620g0030) (Additional file 10: Table S5). DEGs connected with phytohormone signal transduction As phytohormones are known to quickly respond to tobacco topping and to also influence nicotine biosynthesis in tobacco roots, we examined the role of phytohormone signal transduction in the transcriptional responses induced by topping. We identified 336 DEGs, including those related to the biosynthesis, metabolism and action of auxin (IAA), abscisic acid (ABA), ethylene, gibberellin (GA), and jasmonic acid (JA) (Additional file 11: Table S6). The 53 DEGs involved in IAA signal transduction included the ARF (auxin response factor) family (4), the AUX/IAA (auxin responsive protein) family (11), the AUX1 (amino acid transporter protein) family (22), the GH3 (GH3 auxinresponsive promoter) family (7), and the SAUR (auxin responsive SAUR protein) family (9). Most of the DEGs associated with the IAA signaling pathway showed significant up-regulated expression changes, 21 of 36 genes being upregulated at t = 0.5 h, and 10 of 18 genes being up-regulated at t = 8 h. For the ABA signal transduction pathway, six gene families were identified, including the PYL/PYR (abscisic acid receptor) family (5), the SAPK (Serine threonine protein kinase) family (2), the PP2C (protein phosphatase 2C) family (28), the CIPK (CBL-interacting protein kinase) family (19), the CDPK (Calcium-dependent protein kinase) family (9), and the Calmodulin (Calmodulin-like protein) family (8). 43 expression changes were detected at t = 0.5 h, and 19 at t = 8 h. 83 DEGs were implicated in ethylene signaling, including the AP2-EREBP (ethylene responsive transcription factor) family (75), and the ETR (ethylene receptor) family (8), with most DEGs identified at t = 0.5 h (55). The GA and the JA signaling pathways (four and three gene families, respectively) also showed significant transcriptional changes after topping. Quantification of phytohormones and nicotine Phytohormones play a vital role in regulating plant defense and development. To gain insights into the Qin et al. BMC Plant Biology (2020) 20:30 Page 8 of 15 Fig. 7 The clustering of gene expression pattern of DEGs on the TFs and the genes involved in nicotine biosynthesis and transport across different time points after topping in tobacco. The x-axis represents treatment conditions and the y-axis represents centralized and normalized expression values. The red lines indicate the mean expression trends of the TFs (dotted lines) belonging to each cluster. The gene number is marked following the cluster ID mechanisms whereby phytohormones affect the responses of tobacco to topping, we measured the levels of IAA, JA, JA-Ile, and ABA in the root samples at the various timepoints after topping. Both JA and auxin signaling pathways were induced by topping at t = 3 h. JA levels at t = 3 h were significantly increased by almost 34% (P = 0.035, paired t test), and reduced by 23.5 and 18.9% at t = 8 h and t = 24 h (Fig. 8). The dynamics of JA-Ile levels elicited by topping closely followed those of JA, the levels of JA-Ile significantly increasing to about 3-fold at t = 24 h (P = 0.014, paired t test) as compared with untreated plants. The levels of IAA significantly increased at t = 3 h (P = 0.024, paired t test), whilst declining to initial levels at t = 24 h (Fig. 8). The levels of ABA gradually increased to 2.3-fold at t = 8 h (P = 0.0003, paired t test), and to 1.6-fold at t = 24 h (P = 0.009, paired t test), as compared to the untreated plants. We also measured nicotine levels after topping. Our analyses indicated the nicotine levels significantly increased to 1.5-fold at t = 24 h (P = 0.01, paired t test) after topping (Fig. 8). Discussion Transcriptome sequencing and DEGs responses to topping Nicotine is a characteristic secondary product of tobacco. In most Nicotiana species, it is synthesized in the roots, being then transported to the leaves where it accumulates [26]. Nicotine synthesis and accumulation is controlled not only by various environmental cues but also by managerial practices including topping [2, 6]. The factors controlling the topping-induced increase in alkaloid biosynthesis are not well understood, but involve a complex physiological response in the plant as a result of altered phytohormone induced signaling [11]. In order to better understand the mechanism of the tobacco response to topping, we have sequenced and analyzed the transcriptomes of N. tabacum roots at early time points after topping. We identified a total of 4830 topping-responsive DEGs, with representatives being distributed across a number of different molecular functional categories, including secondary metabolism, plant Qin et al. BMC Plant Biology (2020) 20:30 Page 9 of 15 Fig. 8 Mean (+SE) concentrations of phytohormones and nicotine from three replicates of roots harvested after topping treatment within the indicated times. a Mean (+SE) concentrations of JA, JA-Ile, IAA, and ABA from N. tabacum roots were measured using HPLC-MS/MS within 24 h after topping. b The level of nicotine from roots within 24 h after topping treatment, untreated plants served as controls. Asterisks represent significantly different hormone and nicotine levels between control and treatment plants after the indicated times. (t-test; *, P < 0.05; **, P < 0.01) hormone signaling transduction, stress defense, and other metabolism. DEGs involved in nicotine biosynthesis and subgenome transcriptional preference We detected 1.5-fold changes in nicotine levels at t = 24 h after topping with each biological replicates containing 4 individual plants, which is consistent with previous report [27]. It was worth to mention that no significant changes of nicotine content were detected at this time point with each biological replicates containing a single plant (data not shown), indicating that the individual plants response differently in the short time of decapitation. We then inspected more closely the transcriptional changes of genes known to be involved in nicotine biosynthesis and transport. Nine gene families in the nicotine biosynthesis pathway were identified within the DEGs. Their functions included pyridine ring synthesis (AO and QS), pyrrolidine ring synthesis (ODC, ADC, PMT, SAMS, MPO), and the coupling of the two nicotine rings (A622 and BBL). Ornithine decarboxylase (ODC) catalyzes the first and rate-limiting step of polyamine biosynthesis that converts ornithine into putrescine. Down-regulation of ODC transcript levels using RNAi led to lower leaf levels of nicotine in N. tabacum [28, 29]. Correspondingly, in our study, one ODC gene was significantly up-regulated at t = 8 h and t = 12 h. A second example involves PMT, which converts putrescine into N-methylputrescine [30]. In that previous study, five PMT genes were investigated (NtPMT1a, NtPMT1b, NtPMT2, NtPMT3 and NtPMT4). Transcripts derived from NtPMT2 and NtPMT1b showed the greatest increase in abundance (about 3-fold) during the first 24 h after topping [31]. In line with these findings, Qin et al. BMC Plant Biology (2020) 20:30 the expression levels of all the five PMT genes in our study were significantly up-regulated at different time points after topping (Additional file 7: Table S2). A further example is provided by SAMS, which indirectly contributes to nicotine biosynthesis by supplying the S-adenosylmethionine cofactor for the PMT reaction [32]. In our study, five SAMS genes were significantly up-regulated at early time points after topping (Additional file 8: Table S3). The enzyme QPT plays a critical role in the synthesis of the pyridine moiety of nicotine in Nicotiana, in addition to its ubiquitous role in NAD(P)(H) synthesis [33]. The tobacco genome contains two duplicated QPT genes (designated QPT1 and QPT2). QPT1 is expressed at a constitutive basal level in all plant tissues, with somewhat higher levels of expression within the apical meristem. In comparison, QPT2 is expressed exclusively in the tobacco root and is regulated coordinately with other structural genes for nicotine biosynthesis [33, 34]. Although QPT1 and QPT2 were not present in our list of DEGs, the gene of QPT2 (Nitab4.5_ 0000742g0010) showed significant up-regulation at 3 h and 24 h in our qPCR assay (data not shown). Finally, we consider A622, which is expressed in the root, and may be involved in the final condensation reaction of nicotine biosynthesis [35]. The capacity of N. glauca to produce anabasine was markedly reduced when an RNAi approach was used to down-regulate gene expression, thereby decreasing A622 protein levels. This resulted in plants having almost undetectable levels of pyridine alkaloids in their leaves, even after suffering damage to apical tissues [18]. In our work, the expression of A622 was significantly up-regulated at t = 8 h and t = 24 h after topping (Additional file 7: Table S2), implying that A622 positively regulates the biosynthesis of nicotine in the final coupling ring step. Consistent with the results for A622, BBL genes were also found to be induced after topping, as previously reported [17]. N. tabacum is hypothesized to be a consequence of hybridization of two parental genomes (N. sylvestris and N. tomentosiformis). N. tomentosiformis displayed much lower nicotine levels than N. sylvestris in both roots and leaves [36]. All DEGs, especially those involved in nicotine biosynthesis and transport, showed no obvious preference to any of the two subgenomes. More experimental work will be required to elucidate the molecular basis of heterosis and the dramatic domestication selection following the hybridization. Differential expressed transcription factors related to regulation of nicotine synthesis The pattern of expression changes induced by topping for the structural genes of nicotine biosynthesis was initially derived from bioinformatics analyses of sequence data, but was confirmed by the qRT-PCR experiments. Page 10 of 15 It was therefore reasonable to speculate that the sequence data could also be used to identify the TF(s) that actively regulate the production of nicotine at the early time points after topping. As reported previously, many TFs play important roles in regulating nicotine biosynthesis, including members of the AP2/EREBP, bHLH, ARF, and WRKY families [37, 38]. The AP2/EREBP family is the largest TF family in the tobacco genome [39], and ERF-type TFs of the group IX subfamily, including ERF1, ERF189, and ERF32, have been recently identified as direct regulators of the structural genes of nicotine biosynthesis [22, 33, 40]. The second largest class of TFs shown to induce alkaloid biosynthesis in Nicotiana is the MYC2-like bHLH family. MYC2, belonging to the bHLHfamily of TFs, is a key component in conserved jasmonate signaling [41]. It positively regulates nicotine production either directly, through G-box-mediated binding and activation of nicotine structural genes, or indirectly, through the activation of the ERF genes [42, 43]. Wang et al. [44] found that overexpression of NtMYC2a led to great enhancement, under field testing, of nicotine levels in the transgenic lines. Although the mechanism by which ARF1 regulates nicotine synthesis remains unclear, VIGS (virus-induced gene silencing) of bHLH3 and ARF1 results in a significant increase in nicotine content as compared to control plants [22]. In addition, WRKYR1, the Group II member of the WRKY family, was specifically and highly expressed in tobacco roots. This suggests it regulates the expression of genes related to nicotine synthesis, like PMT [37]. Screening the DEGs in our study led to the identification of 549 DEGs annotated as TFs, including members of the AP2/EREBP (75), bHLH (44), WRKY (23), and ARF (4) families. (Fig. 6 and Additional file 9: Table S4). To elucidate patterns of co-regulation of TFs, we clustered all DEG TFs along with the structural genes involved in nicotine biosynthesis. We found that Cluster 11 contained 20 of 28 structure and transporter genes in the DEGs, as well as a total of 17 genes from the bHLH and AP2-EREBP families. We selected one of seven genes in the bHLH family (Nitab4.5_0000093g0110) for study using RNAi, finding the resultant plant showed an altered nicotine level (unpublished data), and implying an important regulatory role in nicotine biosynthesis. Other newly identified AP2/ EREBP, bHLH, WRKY, ARF, MYB, and NAC TF genes might also be involved in nicotine biosynthesis, since all of these TF families have been described as functioning in the regulation of plant secondary metabolism [45–47]. These are therefore good targets for further experiments. Involvement of DEGs responsive to topping in phytohormone signal transduction Plant hormones play pivotal roles in regulating numerous aspects of plant growth and development, including