Expression analysis of transcription factors in sugarcane during cold stress

Abstract Transcription factors (TF) are a wide class of genes in plants, and these can regulate the expression of other genes in response to various environmental stresses (biotic and abiotic). In the current study, transcription factor activity in sugarcane was examined during cold stress. Initially, RNA transcript reads of two sugarcane cultivars (ROC22 and GT08-1108) under cold stress were downloaded from SRA NCBI database. The reads were aligned into a reference genome and the differential expression analyses were performed with the R/Bioconductor edgeR package. Based on our analyses in the ROC22 cultivar, 963 TF genes were significantly upregulated under cold stress among a total of 5649 upregulated genes, while 293 TF genes were downregulated among a total of 3,289 downregulated genes. In the GT08-1108 cultivar, 974 TF genes were identified among 5,649 upregulated genes and 283 TF genes were found among 3,289 downregulated genes. Most transcription factors were annotated with GO categories related to protein binding, transcription factor binding, DNA-sequence-specific binding, transcription factor complex, transcription factor activity in RNA polymerase II, the activity of nucleic acid binding transcription factor, transcription corepressor activity, sequence-specific regulatory region, the activity of transcription factor of RNA polymerase II, transcription factor cofactor activity, transcription factor activity from plastid promoter, transcription factor activity from RNA polymerase I promoter, polymerase II and RNA polymerase III. The findings of above results will help to identify differentially expressed transcription factors during cold stress. It also provides a comprehensive analysis of the regulation of the transcription activity of many genes. Therefore, this study provides the molecular basis for improving cold tolerance in sugarcane and other economically important grasses.

cold damage, which is one of the vital mechanisms of cold acclimation (Thomashow, 1999).Presently, ICE1_CBF_COR transcriptional cascade is an obvious and well-known cold acclimation pathway, and in this pathway, C-repeat binding factors (CBF) also called dehydration responsive element binding factors (DREB) are induced rapidly during cold stress and bind to the cold-responsive genes promoter region and activate the transcription of downstream genes (Thomashow, 1999;Chinnusamy et al., 2006).Low temperature inducible, CBF/DREB elements are the transcription factors binding to C-repeats binding site that get activated under drought and cold stresses are direct evidence of ABA-dependent pathway.Ppy CBFs (Pyrus pyrifolia CBF) studies of the Asian pear suggested that all PpyCBFs showed expression during salinity, drought, abscisic acid stress, and cold or low temperature.Moreover, PpyCBFs showed expression both in ABA-mediated and independent signaling pathways (Ahmad et al., 2019).This CBF/DRE transcription factors family includes CBF1, CBF2, and CBF3 that contain (AP2) DNA-binding domain (Akhtar et al., 2012).The CBF transcription factors are recognized and bound to the C-repeat/dehydration response element motif (CRT/DRE) having conserved sequence (CCGAC) present in many cold-regulated genes promoters that mediate cold tolerance in plants (Shi et al., 2018).
Recently, exploration of transcription factors involved in abiotic and biotic stresses to understand functional divergence in sugarcane has become the main interest of various researchers (Pereira-Santana et al., 2017;Wang et al., 2020).Relatively, few cold-responsive genes have been reported and RNA-sequencing technology can provide useful transcriptional data for a deeper understanding of the molecular process in response to specific stress.Experiments performed on sugarcane (Saccharum spp) under cold stress have identified 34 EST and out of which 20 are novel COR genes and further one of them have shown similarity with CB (Miura and Furumoto, 2013).
The growth period of sugarcane crop in subtropics regions is restricted to eight to nine months due to winter's low temperature and it may usually occur from November to February every year (Dharshini et al., 2016).In sugarcane breeding programs, the tolerance to cold and development of cold-tolerant cultivars is important and requires comprehensive knowledge on molecular mechanisms naturally adapted by tolerant genotypes during cold stress (Dharshini et al., 2016).Very limited studies have been conducted in transcription factors activity during cold stress in sugarcane.Therefore, the present study aims to investigate the expression pattern of transcription factors in sugarcane during cold stress.Further, metabolic pathways were also identified and studied.

Introduction
Saccharum hybrid (Sugarcane) belongs to the Poaceae family that stores a high concentration of sugar in its stem and it is cultivated in both hemispheres about 35° N and 35° S from the equator (Cheavegatti-Gianotto et al., 2011;OECD, 2016).It has been known to mankind since prehistory due to the production of sugar and now it is an important economic crop cultivated in several countries including Pakistan (James, 2008;Figueroa-Rodriguez et al., 2019).Moreover, sugarcane is one of the most important commercial crops that gives essential by products such as sugar, ethanol, molasses, and bagasse (Sipaúba-Tavares et al., 2020).According to the UNFAO, the world's top 10 sugarcane producing countries includes Brazil, India followed by China, Thailand, Pakistan, Mexico, Colombia, Indonesia, Philippines, and United State (Sabiha-Hanim and Halim, 2018).Pakistan ranked fifth with average sugarcane production of 73.40 million tons in an area of 1.2 million hectares (Khan et al., 2019).
During the process of development and growth, plants may expose to various external environmental challenges such as cold, drought, high salinity, and pathogens, which can adversely affect the plant yield, growth, and development (Agarwal et al., 2006;Hussain et al., 2004, Srivastava andKumar, 2020).Cold stress or low temperature is key environmental stress that influences plant development and yield (Palva et al., 2002;Rehman et al., 2019;Rehman et al., 2020).The cold stress comprises of freezing stress (below 0°C) and chilling stress (ranges from 0°C to 15°C) and it can restrict the geographical distribution and cultivation of crops (Tsvetanov et al., 2000).Various plants in the region of tropical and subtropical are chilling sensitive and are injured by nonfreezing low temperatures, while temperate plants are chilling tolerant as these exhibit complex mechanisms to acclimatize themselves (Pareek et al., 2017;Wani and Herath, 2018).
Sugarcane is a cold-sensitive plant and cold stress causes severe damage to this plant (Li et al., 2015).According to Ebrahim et al.(1998), the chilling causes less growth and development of sugarcane, and thus sugarcane develops fewer shorter internodes and leaves, moreover, 37°C was considered the optimum temperature for better growth and development.According to field observation studies, sugarcane sensitivity to cold varies among varieties, moreover, subtropical hybrid plant species are reasonably cold tolerant in contrast to tropical species (Du et al., 1999).Recently, sugarcane has attained the interest of researchers for improving its resilience against cold, especially the molecular mechanism in response to cold stress has become an important issue for its improvement (Yang et al., 2017;Su et al., 2020;Rehman et al., 2020).
Expression analysis of Transcription factor in cold stress

Experimental materials and cold stress treatment
Two sugarcane cultivars namely ROC22 and GT08-1108 were selected (Tang et al., 2018).Cold stress experiment was performed at the Guangxi Sugarcane Research Institute, Guangxi Academy of Agricultural Sciences, China.RNA sequencing data (57.41GB) of the selected sugarcane cultivars were downloaded from the SRA, NCBI database, where the samples were arranged with three biological replicates with a cold stress treatment time point (4°C) and control (28°C) (Table 1).Total RNA was extracted from the selected sugarcane cultivars (cold treated as well as normal sugarcane plants) by using an RNA extraction kit (Tiangen).The extracted RNA was further confirmed by agarose gel electrophoresis and for RNA degradation and contamination detection, Nanodrop was used, while Qubit Fluorometer was used for RNA concentration.Furthermore, Agilent 210 was used to detect RNA integrity.For mRNA enrichment, magnetic beads with oligo (dT) were used, and fragmentation solution was added to interrupt mRNA.cDNA with a sixbase random primer was synthesized from the mRNA as a template.Buffer, dNTPS, and DNA polymerase were used to synthesize a double-stranded cDNA library.After cDNA library construction, the samples were sent for the Illumina High throughput sequencing to Beijing Nuohe Zhiyuan Biological Company.

Reading and splicing of transcriptome raw data
The transcriptomic raw data of the selected sugarcane cultivars were analyzed for quality control with FASTQC (https://www.bioinformatics.babraham.ac.uk/projects/ fastqc/) and low-quality sequences, as well as Illumina adaptors, were removed with Trimmomatic (Bolger et al., 2014).The sugarcane reference genome was used to align the reads by using STAR software and HTSeq-count was used for counting the reads aligned in each gene (Zhang et al., 2018).The number of reads from each library aligned in each gene was used as the measurement of gene expression levels.

Differential expression analysis and normalization of data and transcription factor
For differential expression analysis between treatments and normalization of RNA sequence data, R library edgeR was used (Robinson et al., 2010).Initially, the samples were filtered for normalization, and the expression of genes with 1 count per million (CPM) for at least two samples were excluded.TMM (Trimmed Mean) (Robinson and Oshlack, 2010) method was used for estimating normalization factors of each sample, and to equalize and minimize the volume and high expression of genes and differences in sample size.
During differential expression analysis, the package edgeR assumes the negative binomial (NBD) distribution of the sequences count digital discrete data related to each sample.Therefore, Conditional Maximum Likelihood (CML) estimation was used for sequential analysis of gene to gene (tag-wise) dispersions, and Bayes empirical strategy assumed the shrinkage in tag-wise (gene-wise) distribution of residual variation.Differential gene expression analysis (DGE) and over dispersed data were estimated by Fisher Exact test (Robinson and Smyth, 2008).The likelihood ratio test (LRT) was carried out for both low temperature and control experimental conditions.,Moreover, False Discovery Rate (FDR) (P-value of) gave us the significance of the differential expression with a threshold of an FDR<0.05.Log fold change ratio (LogFC) measured the magnitude of differential expression (Quackenbush, 2002).The upregulated and downregulated genes under cold stress were defined as the genes expressed with LogFC ≥ 2 and LogFC value ≤-2 value respectively.Further, the DEGs were grouped with cluster analysis, annotated, and visualized within biochemical pathways.

Gene ontology and gene ontology enrichment analysis
Agrigo (Tian et al., 2017) gene enrichment analysis was performed for the differentially expressed transcription factors in sugarcane cultivars during cold stress.Besides, Agrigo is web-based online database and kit of tools to provide deep support to the agriculture community for gene enrichment analysis (Zhou et al., 2010;Tian et al., 2017).The Gene Ontology and enrichment analysis were categorized into three that are: 1) Biological process, 2) Molecular function, 3) Cellular components.

RNA sequence data analysis
The RNA sequence data of both cultivars (ROC22 and GT081108) were downloaded from SRA, NCBI database.The adapter and unwanted reads were removed, and the reads were mapped to a recently published reference genome (Zhang et al., 2018).The edgeR package was used to perform expression analysis.Our targeted genes (Transcription factor activity genes) were selected, and their expression analysis was studied.In our data analysis, the transcription factor was found significantly upregulated in both cultivars (ROC22 and GT08-1108).In ROC22, 963 transcription factor genes were differentially upregulated in 5558 genes, while 283 genes were downregulated in 3252 genes (Table 2).While in GT08-1108, 973 transcription factor genes were upregulated in 5649 genes, and 293 genes were down-regulated in 3289 genes (Table 2).

Transcription factors upregulated in sugarcane cultivars during cold stress
The transcription factor activities were found upregulated in protein binding, transcription factor binding, the activity of DNA-sequence specific binding transcription factor, activity in transcription factor complex, transcription factor activity in RNA polymerase II, the activity of nucleic acid binding transcription factor, transcription corepressor activity, transcription factor binding of RNA polymerase II to DNA sequence-specific regulatory region, the activity of TF of RNA polymerase II, transcription factor cofactor activity, transcription factor activity at plastid promoter, transcription factor activity from RNA polymerase I promoter, transcription factor activity from the promoter of RNA polymerase II and RNA polymerase III, in both cultivars (Table 3 and 4).Mostly, the transcription factor activity was found upregulated in sequence-specific DNA binding up to 57% in total DGE (Figure 1).

Transcription factors downregulated in sugarcane cultivars during cold stress
In the present data analysis, the transcription factor activities were found highly upregulated.Although transcription factor activity in our analysis was also found downregulated in both cultivars.Genes (TF activity) are differentially downregulated in protein binding, transcription binding, the activity of sequence-specific DNA binding transcription factor, transcription factor complex, transcription factor activity of RNA polymerase II, the activity of nucleic acid binding transcription factor (Table 5 and 6).In downregulated transcription factor (TF) activity, the activity of transcription factor in protein binding was found 31% among all downregulated genes (Figure 2).

Metabolic pathways of the transcription factor in sugarcane during cold stress
The gene ontology and gene enrichment analysis were performed using Agrigo gene enrichment analysis.The transcription activity was enriched in molecular  function and cellular components while in biological processes no transcription factor activity was found.In both cultivars, the transcription factor activity was found highly upregulated in molecular function and cellular components (Figure 3 and 4).Further, the transcription factor activity was also found downregulated in molecular function and cellular components in both cultivars (Figure 5 and 6).

Discussion
Cold stress is one of the major abiotic stresses which reduces the quality, yield, and growth of many economically important crops worldwide (Zhou et al., 2011;Knight and Knight, 2012;Rehman et al., 2020).During cold stress, plants adjust or adjust their biochemical and physiological activities by reprogramming the expression of many coldresponsive genes (Bohnert et al., 1995;Chinnusamy et al., 2007).In subtropical regions, chilling stress is one of the major constraints for sugarcane cultivation, which hampers yield and sugar production (Mathur et al., 2021).In plants, mainly the physiological responses come from (COR) cold-responsive genes.The expression of these genes is controlled by transcription factors (Fu et al., 2018).
Currently, the most well-defined cold signaling pathway is the C-repeats binding transcriptional factor (CBF) regulatory cascade (Thomashow, 1999;Chinnusamy et al., 2003Chinnusamy et al., , 2007;;Fu et al., 2018;Knight and Knight, 2012).In our study, the transcription factor activity was found highly upregulated in sugarcane cultivar ROC22, 963 TF genes were significantly upregulated among a total of 5,649 upregulated genes, while 283 TF genes were found  downregulated among the total of 3,289 downregulated genes.Similarly, in cultivar GT08-1108, about 974 TF genes were identified among 5,649 upregulated genes, while 293 TF genes were found among 3,289 downregulated genes.Transcription factors play a significant role under diverse stress conditions in plants by either downregulation or up-regulation of genes at the transcriptional level (Javed et al., 2020).Almost in the vascular plant genome 7% coding capacity is attributed to TFs for regulating the gene expression (Rushton et al., 2008).In plants, thousands of TFs families have been identified, while the major TFs families including WRKY, NAC, MYB, AP2/ ERF, etc. mediated through signal transduction pathways have been used to cope with biotic and abiotic stresses in diverse crops (Javed et al., 2020).
Although genome-wide RNA expression analysis has become a routine tool in biological research, extracting   biological insight from such information remains a major challenge a powerful analytical method called Gene Enrichment Analysis (GEA) for interpreting gene expression data.The method derives its power by focusing on gene sets, that is, groups of genes that share common biological functions, chromosomal locations, or regulation (Subramanian et al., 2005).In the GO enrichment data analysis, the transcription factor activity was annotated with GO categories related to protein binding, transcription factor binding, the activity of sequence-specific DNA binding transcription factor, activity in transcription factor complex, transcription factor activity in RNA polymerase II, the activity of nucleic acid binding transcription factor, transcription corepressor activity, transcription factor binding of RNA polymerase II to the regulatory region of sequence-specific DNA, the activity of transcription factor of RNA polymerase II, transcription factor cofactor activity, transcription factor activity from plastid promoter, transcription factor activity from RNA polymerase I promoter, transcription factor activity from the promoter of RNA polymerase II and RNA polymerase III.Transcription factors (CBFs) are also known as DRE/CRT (Dehydration responsive elements) in the promoter region of cold-responsive genes, which regulate the expression as well as function (Stockinger et al., 1997;Liu et al., 1998;Vazquez-Hernandez et al., 2017).ICE1 (Inducer of CBF expression 1), an MYC-like bHLH transcriptional activator, acts as a positive regulator of transcription factors during cold signaling by binding to MYC recognition elements in the CBF promoter (Chinnusamy et al., 2003(Chinnusamy et al., , 2010;    Fursova et al., 2009).In our analysis, the transcription factor complex was highly upregulated by 4.9 logFc2.In GT08-1108, 630 upregulated genes for transcription factor activity, 43 for transcription factor activity in protein binding, 39 genes in RNAII polymerase transcription factor activity, 40 genes in transcription factor binding, 270 genes for transcription complex were recorded.the numbers of expressed upregulated transcription factors in ROC22 are 610 genes were upregulated for transcription factor activity, 43 for transcription factor activity in protein binding, 38 genes in RNAII polymerase transcription factor activity, 40 genes in transcription factor binding, 270 genes for transcription complex.

Conclusion
In Saccharum hybrid (Sugarcane), the lack of genome and transcriptome data has limited molecular analysis.In this research work, we provided the investigation about transcription factor activity during cold stress in two sugarcane cultivars namely GT08-1108 and ROC22.The results of this study showed that transcription factor genes play an important role in the expression of coldresponsive genes during cold stress in sugarcane.We further suggested to investigate the function of these TFS in drought, salinity, climate change, and extreme temperature in economically important grass crops like sugarcane, wheat, and rice.Finally, the development of cold and other environmental stress-tolerant sugarcane cultivars will be substantially beneficial for the sugar industry.

Figure 3 .
Figure 3. Schematic view of the metabolic pathway for transcription factors upregulated in sugarcane cultivar GT08-1108 during cold stress.

Figure 4 .
Figure 4. Schematic view of the metabolic pathway for transcription factors upregulated in sugarcane cultivar ROC22 during cold stress.

Figure 5 .
Figure 5. Schematic view of the metabolic pathway for transcription factors downregulated in sugarcane cultivar GT08-1108 during cold stress.

Figure 6 .
Figure 6.Schematic view of the metabolic pathway for transcription factors downregulated in sugarcane cultivar ROC22 during cold stress.

Table 1 .
Samples of sugarcane analyzed.

Table 2 .
Differentially expressed upregulated and downregulated transcription factors during cold stress in sugarcane cultivars ROC22 and GT08-1108.

Table 3 .
Transcription factor activity of upregulated genes during cold stress in ROC22 cultivar.
Brazilian Journal of Biology, 2023, vol.83, e242603 5/10 Expression analysis of Transcription factor in cold stress

Table 4 .
Transcription factor activity of upregulated genes during cold stress in GT08-1108 cultivar.

Table 6 .
Transcription factor activity of downregulated genes during cold stress in sugarcane cultivar GT08-1108.

Table 5 .
Transcription factor activity of downregulated genes during cold stress in sugarcane cultivar ROC22.