Genome-wide identification and characterization of SnRK family genes in Brassica napus

DOI: https://doi.org/10.21203/rs.3.rs-17109/v1

Abstract

Background: Sucrose non-fermenting 1 related protein kinases (SnRK) play crucial roles in responding to biotic and abiotic stresses through activating protein phosphorylation pathways. However, little information of SnRK genes was available in Brassica napus, one of important oil crops. Recently, the released sequences of the reference genome of B.napus provide a good chance to perform genome-wide identification and characterization of BnSnRK gene family in the rapeseed.

Results: Totally 114 SnRK genes distributed on 19 chromosomes were identified in the genome of B.napus and classified into three subfamilies on the basis of phylogenetic analysis and the domain types. According to gene structure and motif composition analysis, the BnSnRK sequences showed obvious divergence among three subfamilies. Gene duplication and synteny between the genomes of the rapeseed and Arabidopsis were also analyzed to provide insights into the evolutionary characteristics of BnSnRK family genes. Cis-element analysis revealed that BnSnRKs may response to diverse environmental stresses. Moreover, the expression patterns of BnSnRKs in various tissues and under diverse abiotic stresses were distinct difference. Besides, Single Nucleotide Polymorphisms (SNP) distribution analysis suggests the function disparity of BnSnRK family genes in different genotypes of the rapeseed.

Conclusion: We examined genomic structures, evolution features, expression patterns and SNP distribution of 114 BnSnRKs. The results provide valuable information for functional characterization of BnSnRK genes in future studies.


Introduction

Plants develop various molecular defense mechanisms to cope with abiotic stresses, including salinity, drought and cold stresses. Gene expression regulation and protein modification are two important ways for plants to deal with these stresses (Bohnert et al., 2006). The processes of phosphorylation and dephosphorylation mediated by protein kinase play crucial roles in protein modification (Hunter et al., 1995). Among the reported protein kinase genes, sucrose non-fermenting 1 (SNF1)-related protein kinases (SnRKs) are involved in different physiological processes (Hrabak et al., 2003).

In plants, SnRKs could be divided into three subfamilies: SnRK1, SnRK2 and SnRK3 based on their sequence similarity and gene structures (Hrabak et al., 2003; Halford and Hardie, 1998). In detail, the SnRK1 subfamily contain a highly conserved N-terminal protein kinase (Pkinase) domain, which are the homologous genes of SNF1 in yeasts and AMPKs in mammals (Celenza et al., 1986). The other subfamilies SnRK2/3 are unique in plants, both are more diverse in plants than the SnRK1 subfamily members. The member of SnRK2 harbors a conserved Pkinase domain and a C-terminal variable adjusting conserved domain (Kulik et al., 2011). In addition, SnRK3 known as CIPK (CBL-interacting protein kinases), also contain conserved N-terminal protein kinase domains and NAF domains, PPI domains in C-terminal (Albrecht et al., 2011; Masaru et al., 2003).

The SnRK1 family genes involve in the response of plant cells to starvation and metabolic stress. In Solanum tuberosum, SnRK1 was proved to be involved in induction of sucrose synthase expression and played an important role in carbohydrate metabolism regulation (Purcell et al. 1998). Besides,ow energy stress (e.g. darkness and hypoxia) could trigger SnRK1α nuclear translocation and further induce SnRK1 target genes to refresh cellular energy for plant growth (Ramon et al., 2019). Furthermore, considerable evidences indicated that SnRK1 genes were hubs in a network of various signaling pathways including cell cycle regulation, pathogen responses and meristem development (Halford et al., 2009).

On the other hand, the SnRK2 genes play important roles in responding to abiotic stresses in plants, especially for osmotic and salt stress. For instance, SnRK2.10 phosphorylate several target genes including dehydrins ERD10 and ERD14 to deal with osmotic stress in Arabidopsis thaliana (Maszkowska et al., 2019). SnRK2.1 positivly regulate salt tolerance in Nicotiana tabacum (Diédhiou et al., 2008; Zhang et al., 2014). SnRK2 genes also respond to drought stress through being involved in regulation of ABA pathway. Besides, OST1 (AtSnRK2.6) plays a core role in the ABA signaling pathway in stomatal guard cells, and OST1 protein stability can be regulated via E3-ubiquitin ligase HOS15 to reduce ABA signal sensitivity in Arabidopsis (Anna-Chiara et al., 2002; Ali et al., 2019).

SnRK3 kinases known as CIPKs (CBL(calcium sensor calcineurin B-like proteins)-interacting protein kinases), perform vital functions in resistance to various stresses in plants (Hrabak et al., 2003; Kim et al., 2003). For example, SOS (salt overly sensitive) system was the first discoverd CBL-CIPK pathway in A.thaliana. In detail, the calcium signal produced by salt stress was sensed by SOS3 (AtCBL4) on cell membrane, and then SOS3 combined with SOS2 (AtCIPK24) forming complex to phosphorylate SOS1 (Na+/H+ antiporter) to remove excess Na+ out of root cells (Guo et al., 2002; Liu et al., 2000). Besides, MdCIPK13 and MdCIPK22 enhanced salt and drought tolerance through targeting sucrose transporter MdSUT2.2 for phosphorylation in apple (Ma et al., 2019a; Ma et al., 2019b). Overexpression of BnCBL1-BnCIPK6 in B.napus could enhance high salinity tolerance and low potassium tolerance (Chen et al., 2012). In conclusion, increasing evidences emphasized the importance of SnRKs function in nutrition utilization and stress response, and finally researchers could improve plants resistance to stresses by genetic manipulation using these genes.

Brassica napus is an important oil crop in the world. Recently, the genomes of Darmor-bzh (winter ecotype), Zhongshuang 11 and NY7 (seni-winter ecotype) were successfully sequenced and assembled (Chalhoub et al., 2014; Sun et al., 2017; Zou et al., 2019); however, the systematic analysis of BnSnRK gene family has been not well reported. In this study, 114 SnRK gene members were identified in the B.napus genome. We systematically analyzed their phylogenetic relationships, protein motifs, gene structures, chromosome distributions and cis-elements in promoter regions. Moreover, the expression profiles of BnSnRKs in diverse tissues as well as under abiotic stresses were determined. In addition, SNPs of each BnSnRK gene were systematically identified in a worldwide collection with 300 core germplasm accessions. These results will provide useful information for further investigation of molecular mechanisms of BnSnRK genes for abiotic stress tolerance and molecular breeding in B.napus.

Methods

Identification of BnSnRK family genes in the B.napus genome

The amino acid sequences of SnRKs gene family in A.thaliana were downloaded from the NCBI (https://www.ncbi.nlm.nih.gov/). The homologous genes of AtSnRKs in B.napus genome were blasted in the reference genome of rapeseed cultivar Ningyou 7 (http://ibi.zju.edu.cn/bnpedigome/download.php?con=ny7) (Zou et al., 2019). Hidden Markov Model (HMM) and BLASTP program were applyed for preliminary identification of BnSnRK proteins. Local BLASTP (E-value-20) searches were performed based on Hidden Markov Model (HMM) profiles of SnRK gene domains from the Pfam database (http://pfam.janelia.org/). All candidate sequences of SnRK genes were reconfirmed through the SMART database (http://smart.embl-heidelberg.de/) (Letunic et al., 2012), the NCBI Conserved Domain database (Marchler-Bauer et al., 2011) and the Pfam database (Finn et al., 2016). Moreover, the number of amino acids, molecular weights (MW) and isoelectric point (pI) of each SnRK protein were calculated using tools from ExPASy (http://www.expasy.ch/tools/pi_tool.html).

Phylogenetic analysis and classification of SnRK gene family in B.napus

The ClustalW with default parameters was used for multiple sequence alignment of 114 BnSnRK non-redundent amino acid sequences (Saitou and Nei, 1987; Kumar et al., 2016). Based on the alignments, a phylogenetic tree was constructed using MEGA 7.0 by the Neighbor-Joining (NJ) method (Kumar et al., 2016), with the following parameters: poisson model, pairwise deletion and 1000 bootstrap replications. Unrooted NJ tree of all SnRK protein sequences from A.thaliana and B.napus was also constructed using MEGA 7.0.

Identification of motif compositions and gene structures

To identify conserved motifs of BnSnRK proteins, the Multiple Expectation Maximization for Motif Elicitation (MEME) online program (Bailey et al., 2009) (http://meme.sdsc.edu/meme/itro.html) was performed with the following parameters: number of repetition = any, maximum number of motifs = 10; and optimum motif length = 6 to 100 residues. The exon-intron structures of BnSnRK family genes were analyzed using the Gene Structure Display Server online program (GSDS: http://gsds.cbi.pku.edu.ch) (Hu et al., 2015).

Chromosomal location and gene duplication

The chromosomal positions of all BnSnRK genes were mapped to 19 chromosomes of the rapeseed genome according to the physical location information from the NY7 genome database using MapChart version 3.0 and Circos (Voorrips, 2002; Krzywinski et al., 2009). To identify gene duplication, all B.napus gene sequences were aligned using BLASTP, with an e-value of 1e-10. MCScanX with default values was used to classify the duplication patterns of the SnRK into segmental, tandem duplications (Wang et al., 2012). The definition of tandem duplication is a chromosomal region within 200 kb containing two or more genes (Holub et al., 2001). To exhibit synteny relationships of the orthologous SnRK genes between B.napus and A.thaliana, the syntenic analysis maps were constructed using python script written by ourselves. Non-synonymous (Ka) and synonymous (Ks) substitution of each duplicated BnSnRK gene were calculated using KaKs_Calculator 2.0 (Wang et al., 2010). Divergence time was estimated using the formula T = Ks/2R, where R is 1.5×10−8 synonymous substitutions per site per year (Wang et al., 2010).

Cis-elements in promoter regions of BnSnRKs

Upstream sequences (1500 bp) from the start codon of each BnSnRK gene was extracted from the genome sequence of NY7, and then used for cis-element distributions in promoter regions using PlantCARE software (http://bioinformatics.psb.ugent.be/webtools/plantcare/html/) (Lescot et al., 2002).

Expression patterns of BnSnRKs

Transcriptome data of 12 different B.napus tissues were obtained from a previous study by Sun et al. (2017) under the project ID of PRJNA394926 in the NCBI. Moreover, transcriptome data of B.napus under dehydration, salt, ABA induce and cold stress condition were described in Zhang et al. (2019). The data was available at (https://bigd.big.ac.cn/) under the project ID of CRA001775. Differential expression analysis of BnSnRKs was performed using the DSEeq2 R package and the heatmaps were created by TBtools software (Chen et al., 2018).

SNP distribution analysis of BnSnRKs in B.napusaccessions

The SNPs in the coding regions of BnSnRK genes were extracted from the total SNPs of 300 B.napus accessions, which were determined by the genome re-sequencing in our previous studies (Wu et al., 2019; Xuan et al., 2020). High-quality SNPs with MAF greater than 5% and missing rate less than 50% were used for further analysis. All high-quality SNPs were mapped to the “Darmor-bzh” genome (B.napus v4.1 genome, http://www.genoscope.cns.fr/brassicanapus/data/).

Results

Identification and phylogenetic analysis of SnRKgenes in B.napus.

A total of 114 proteins with Ser/Thr kinase domain were identified as the members of SnRK family in the B.napus genome (Table S2). The longest amino acid sequence of each protein was selected for further analyses. The information of gene names, ID, chromosomal locations, amino acid numbers, molecular weights (MW), isoelectric points (pI) and domains were listed in Table S1 and S2. The amino acid length of 114 BnSnRKs protein is ranged from 190 to 1241, correspondingly the melocular weight ranges from 22.0 to 140.6 kDa. The coding regions and sequences of each gene were listed in Table S3.

To analyze evolutionary relationships of SnRK genes in B.napus and A.thaliana, an unrooted phylogenetic tree was constructed using the full-length amino acid sequences of all SnRKs. Totally, 39 SnRKs from A.thaliana and 114 SnRKs from B.napus were identified and used in the study (Figure 1). It was reported that 39 AtSnRKs could be clustered into three groups (Hrabak et al., 2003). In this study, based on the phylogenetic analysis, 114 BnSnRKs were also classfied into three groups. In details, 10 proteins are in BnSnRK1 subfamily with Pkinase (PF00069 of Pfam), UBA (PF00627) and KA1 (PF02149) domains, excepting BnSnRK1.3, while 31 proteins contain pkinase domain in BnSnRK2 subfamily with high similarity to AtSnRK2 subfamily, and 73 proteins with Pkinase and NAF (PF03822) domains are in SnRK3 subfamily, respectively (Figure 1).

Gene structure and protein motif analysis of BnSnRK genes

The web server GSDS (Gene Structure Display Server) analysis was done to determine gene structures of BnSnRKs. As shown in Figure 2, one to fifteen exons were found in these BnSnRK genes. The BnSnRK1 subfamily genes have more than 6 exons, and 6 to 13 exons for the BnSnRK2 subfamily genes. While it is various for the BnSnRK3 subfamily with ranged from 1 to 15 exons. Furthermore, there are two subgroups in the BnSnRK3 family. The genes in the SnRK3 subgroup 1 contained more than 8 exons, while the genes in the subgroup 2 contain less than 4 exons.

The schematic structures of all BnSnRK proteins were constructed using the MEME (Multiple Em for Motif Elicitation) motif analysis results (Figure 2C). The sequence and length information of the conserved motif were shown in Supplementary Table S4. In this study, it was found that all BnSnRK genes retained the conserved Pkinase domain containing the motif 1, 2, 3, 4, 5 (Figure 2C). Moreover, BnSnRK genes within the same subfamily share a similar motif composition, while BnSnRK genes between diverse subfamilies showed distinct differences in motif composition. In detail, BnSnRK1 subfamily genes contain motif 3, 4, 5, 6, 7; BnSnRK2 genes have motif 1, 2, 5, 6, 10 or motif 1, 4, 5, 6, 10; BnSnRK3-1 genes contain motif 1, 2, 5, 6, 7, 8 or motif 1, 4, 5, 6, 7, 8; BnSnRK3-2 genes have motif 1, 5, 6, 7, 8 or motif 1, 6, 7, 8, 9 (Figure 2C). In summary, the conserved motif compositions and similar gene structures of the SnRK genes within the same subfamily, strongly support the reliability of the subfamily classifications by phylogenetic analysis.

Chromosomal distribution, genome synteny and gene duplication of BnSnRK genes

Chromosomal analysis showed that 110 BnSnRK genes were distributed over 19 chromosomes, excepting 4 genes (BnSnRK1.10, BnSnRK2.30, BnSnRK2.31, BnSnRK3.73) could not be mapped into chromosomes of the rapeseed genome (Figure 3). Among them, 56 BnSnRK genes were located in the AA subgenome, including 4 BnSnRK1 genes, 14 BnSnRK2 genes, 38 BnSnRK3 genes; while 54 genes were located in the CC subgenome, containing 5 BnSnRK1genes, 15 BnSnRK2 genes, 34 BnSnRK3 genes. In addition, some BnSnRK2 and BnSnRK3 subfamily genes were formed as clusters in diverse chromosomes, such as BnSnRK3.22 and 3.23 (Figure 3). These results indicated that BnSnRK genes were randomly distributed in the chromosomes of the AA or CC subgenome.

Using BLAST and MCScanX methods, 106 segmental duplication events were identified (Figure 4 and Table S5). Among these events, 104 events were happened between diverse chromosomes, while only 2 duplication events were detected within same chromosome. The result suggested that part of BnSnRK genes were possibly generated by gene duplication, and the segmental duplication events played vital roles in the expansion of SnRK genes in the B.napus genome. We also analyzed the occurrence of the tandem duplication events.

Here, 25 tightly linked BnSnRK genes (BnSnRK2.7/2.8, BnSnRK3.3/3.4, BnSnRK3.18/3.19, BnSnRK3.22/3.23, BnSnRK3.24/3.25, BnSnRK3.33/3.34, BnSnRK3.36/3.37, BnSnRK3.42/3.43, BnSnRK3.49/3.50, BnSnRK3.58/3.59/3.60, BnSnRK3.69/3.70, BnSnRK3.71/3.72) were found located less than 200 kb. However, the identities of these gene pairs were less than 70%, indicating they were not tandem duplication events.

Futhermore, the synteny of SnRK gene pairs between B.napus genome and A.thaliana genome was performed. The result showed that 65 SnRK genes of B.napus exhibiting syntenic relationship with AtSnRK genes (Figure 5 and Table S6), suggesting that these genes might have contributed to the evolution of BnSnRK gene family. To test the evolutionary constraints acting, Ks values, Ka values, Ka/Ks ratios and divergence time of paralogues and orthologues on SnRK family genes were calculated (Table S5, Table S6). The majority of segmental duplicated BnSnRK gene pairs had Ka/Ks ratios less than 1, and the mean value of BnSnRK3 gene pairs (Ka/Ks = 0.30) was lower than BnSnRK1 (Ka/Ks = 0.41) and BnSnRK2 (Ka/Ks = 0.75). In addition, the Ka/Ks ratios of most orthologous SnRK gene pairs were less than 1, and the mean value of SnRK2 gene pairs (Ka/Ks = 0.09) was much lower than that of SnRK1 (Ka/Ks = 0.75) and SnRK3 (Ka/Ks = 0.23). These results suggested that the BnSnRK gene family may have suffered robust purifying selective pressure during evolution.

Stress-related cis-elements in the promoters of BnSnRK genes

To understand the potential function and regulatory mechanisms of BnSnRK genes, cis-elements (1.5-kb upstream from ATG) were analyzed by using PlantCARE. Totally 104 of 114 BnSnRKs were identified with cis-elements including DRE, ABRE and LTRE, involving in dehydration responses, ABA responses and low-temperature responses (Figure 6, Table S2). In detail, 88 BnSnRK genes (77.19%) had ABRE cis-elements, 33 BnSnRK genes (28.95%) carried DRE cis-elements, and 44 BnSnRK genes (38.60%) owned LTRE cis-elementsIt was also found that most genes had more than one kind of cis-element types. FBnSnRK1 (1.80) family was higher than BnSnRK2 (1.17) and BnSnRK3 (1.21) families (Table S2). In conclusion, the cis-elements analysis suggested that most BnSnRK genes could be responsed to diverse environmental stresses, and different subfamily genes may be regulated diversely.

Expression profiles of BnSnRKs in different tissues

The expression profiles of 114 BnSnRK genes were compared among 12 tissues of the rapeseed cultivar ZS11 (Figure 7 and Table S7). The expression pattern of BnSnRKs was divided into three groups. A total of 67 genes with high expression level in all examined tissues (log2-based values > 6) were classified into the group 1. For example, BnSnRK3.40 was highly expressed in all vegetative organs with average log2-based values of 10.85. 21Twenty one BnSnRK genes belonging to group2 exhibited relatively low expression levels across the detected tissues (log2-based values > 2). The group3 included 26 BnSnRK genes, with lowest expression levels in all tissues (log2-based values < 2). Furthermore, 6 genes belonging to the group3 (BnSnRK2.26/2.27 and BnSnRK3.18/3.22/3.43/3.59) were even not expressed in all tissues. Meanwhile, the group1 contained 6 BnSnRK1, 21 BnSnRK2, 40 BnSnRK3 genes; the group 2 had 2 BnSnRK1, 4 BnSnRK2, 15 BnSnRK3 genes and the group 3 consisted of 2 BnSnRK1, 6 BnSnRK2, 18 BnSnRK3 genes. These results demonstrated that BnSnRKs displayed diverse expression patterns, and genes within the same subfamily also expressed differently.

Expression profiles of BnSnRKs under different abiotic stresses

The expression pattern of BnSnRK genes under various abiotic stress were studied using the published transcriptome data of B.napus under drought, salinity, ABA induction and cold stresses (Zhang et al, 2019). Overall, BnSnRK genes significantly changed the expression level under different abiotic stresses (Figure8, Table S8). Multiple BnSnRK genes were significantly induced by diverse treatments. For example,BnSnRK2.1 was extremely induced by all treatments. BnSnRK3.39 incerased expression levels responding to dehydration and ABA treatments. In contrast, several BnSnRKs did not respond to any abiotic stresses. for instance, BnSnRK2.12 and BnSnRK2.13 showed almost no expression changes during all treatments. Interestingly, many genes showed opposing expression profiles under diverse treatments. For instance, BnSnRK3.54 was highly up-regulated by dehydration, whereas was repressed by ABA and cold treatments.

SNPs analysis of BnSnRKs in a core germplasm of B.napus

The polymorphism of BnSnRK genes was determined using our previous resequencing data of a core germplasm resource containing 300 accessions (Wu et al., 2019; Xuan et al., 2020). The SNPs with MAF more than 5% in BnSnRK genes were selected out (Table S9). It was shown that each BnSnRK genes contained an average of 11.76 SNPs. In detail, each BnSnRK1 subfamily genes included 14.00 SNPs, each BnSnRK2 and BnSnRK3 subfamily genes contained an average of 9.67 and 12.53 SNPs. Furthermore, the SNP density of each BnSnRK genes within the same subfamily was also different. For instance, there was no SNP identified in BnSnRK1.1, whereas there were 62 SNPs in BnSnRK1.7.

Moreover, a detailed SNP distribution analysis of BnSnRK2.10 and BnSnRK3.39 was conducted because the two genes significantly changed their expression level under drought stress (Figure 8). It was found that there were 9 SNP loci in the 1500 bp promoter region , 19 SNPs in the exon/intron region and 1 SNP in the 3’UTR region of BnSnRK2.10; while we identified 6 SNP loci in promoter region, 8 SNPs in exon/intron region and 1 SNP in the 3’UTR region of BnSnRK3.39. Combined with phenotypic data of the core germplasm resource under drought stress (data not shown), it was found that SNP16, SNP26, SNP27 in BnSnRK2.10, and SNP4, SNP7, SNP15 in BnSnRK3.39 were significantly associated with their drought tolerance.

Discussion

In this study, 114 members of BnSnRK genes were identified in the B.napus genome, which were designated as BnSnRK1.1 to BnSnRK3.73 on the basis of their subfamily classification. Systematically investigation on the BnSnRK gene family were carried out, inculding phylogenetic relations, gene structures, protein motifs, chromosome distributions, gene duplication and cis-elements in the promoters. Furthermore, the deep analysis was done on the expression pattern and SNP determination of BnSnRK family genes using published data. This study provides a basic information for further functional characterization of SnRK genes to enhance plant adaptive capacity under abiotic stress.

Previous studies identified 39, 48, 44 and 34 SnRK genes in Arabidopsis thaliana (Hrabak et al., 2003), Oryza sativa (Kobayashi et al., 2004), Brachypodium distachyon (Wang et al., 2015) and Eucalyptus grandis (Wang et al., 2019), respectively. In B.napus genome, BnSnRK gene number is much higher than diplont plants. The 114 BnSnRK genes were identified and classified into three subfamilies, including 10 BnSnRK1 genes, 31 BnSnRK2 genes and 73 BnSnRK3 genes. More specific description showed similar member proportions of each subfamilies between B.napus and other species. B.napus originates from natural crossing of B.rapa (AA) and B.oleracea (CC) (Chalhoub et al., 2014). 56 and 54 BnSnRK genes were located in the AA subgenome and the CC subgenome, which indicated that SnRK genes played a similarly important functional role in both ancestral species.

Various SnRK subfamily genes contain different conserved domains, but all genes retained a N-terminal protein kinase domain. For example, SnRK3 subfamily genes were reported interact with CBLs in a calcium-dependent manner because of the NAF domain permits. In addition, the NAF domain defines a set of heterologous kinases involved in diverse signaling processes, as targets of CBL calcium sensor proteins (Albrecht et al., 2001). In this study, it was also found that different BnSnRKs subfamily genes shared the different type of conserved domains. It may suggest there is functional diversification of the BnSnRK gene family according to their domains.

In AtSnRK and BnSnRK gene families, different subfamily genes exhibited significant gene length and exon-intron structural divergences. In previous studies, genes with less introns were considered to have higher expression levels in plants (Chung et al., 2006; Jeffares et al., 2008). A compact gene structure with less introns allowed rapid gene activation and timely response to diverse environmental stresses (Jeffares et al., 2008). However, combined with transcriptome data used in this study, we did not detect that BnSnRK genes with less introns showed higher expression levels than the other BnSnRKs.

Accumulating evidence suggested that gene activities were generally correlated with disparities in promoter regions (Stephen et al., 2008). The cis-elements located in the promoter regions of genes played key roles in regulating gene expression during growth and environmental changes (Ali et al., 2006; Sarda et al., 1997). The promoter analysis revealed that BnSnRKs contained various types of cis-elements, such as DRE, ABRE and LTRE. Most of gene promoters contained at least one of these components, indicating that many of BnSnRKs were able to respond to diverse abiotic stresses. For instance, by combining gene expression profiles under 4h-cold stress, the expression level of BnSnRKs with LTRE and ABRE cis-elements increased by an average of 2.71-fold, while BnSnRKs with no cis-elements only showed 1.47-fold increase. Therefore, the cis-elements analysis provide a clue for gene function study, especially for gene expression pattern under different stress.

Gene expression profiles provided imperative clues to map out gene functionality. Firstly, we used pre-published transcriptome data to investigate BnSnRK genes expression levels in diverse tissues and organs of B.napus (Sun et al., 2017). The analysis results revealed that the expression patterns of these genes were divided into three groups (Figure7). Meanwhile, We found BnSnRKs in group2 contained less cis-elements in promoter than group1 and group3. Specifically, each gene in group1 contains an average of 2.84 ABRE, 0.40 DRE, 0.93 LTRE, and each gene in group3 contains 2.81 ABRE, 0.38 DRE, 0.30 LTRE, whereas each of group2 gene only has 2.00 ABRE, 0.19 DRE, 0.57 LTRE. These evidences suggested that BnSnRKs activities were correlated with disparities in promoter regions.

The roles function of some BnSnRKs reponding to diverse abiotic stresses were also derived. According to the drought stress data (table S8), the responsive gene BnSnRK2.10 was orthologous to AtSnRK2.3, which regulated ABA synthesis and signaling responding to drought in A.thaliana (Tan et al., 2018), indicating identical function of BnSnRK2.10 in responding to drought stress.BnSnRK2.24, exhibited extremely expression changes in drought, salt stresses and ABA induce, while its orthologs AtSnRK2.2, could also respond to osmotic stress and ABA induction in A.thaliana, indicating that BnSnRK2.24 may share similar functions in B.napus (Kazuo et al., 2009; Yoshida et al., 2015). However, BnSnRK3.39 gene, which could be significantly induced under diverse stresses, was only existed in B.napus. Moreover, Combined with resequencing data of the core germplasm resource, SNPs of each BnSnRK genes were identified. And the SNP distribution analysis of BnSnRKs such as BnSnRK2.10 and 3.39, could further provide functional molecular markers and alleles for these BnSnRK family genes.

This study provides a comprehensive knowledge of SnRK gene family in B.napus. It gives an important implication for further understanding the biological functions of individual BnSnRK genes in B.napus. However, large functional characterization work need to be done in further work to understanding the roles of BnSnRK family.

Conclusion

SnRK genes play important roles in signaling pathways including responses to biotic and abiotic stresses in plants. In this study, a comprehensive study of SnRK gene family in B.napus was performed. A total of 114 BnSnRK genes were characterized and divided into three subfamilies, which showed high similarity in gene structure and motif composition within the same subfamily. Phylogenetic comparison and synteny analysis of SnRK genes between A.thaliana and B.napus provide valuable clues for the evolutionary characteristics of the BnSnRK genes. Moreover, the cis-acting elements, gene expression and SNPs distrubution of BnSnRK family were also determined. These results provide an important information for further understanding biological functions of BnSnRK genes in B.napus.

Declarations

Acknowledgments

This work was funded by the National Natural Science Foundation of China (31961143008, 31701411), the Science and Technology Program of Zhejiang Province of China (LGN20C130007), Jiangsu Collaborative Innovation Center for Modern Crop Production, and the 111 project for introduction of foreign experts (B17039). We thank Dr. Guoping Zhang for his insightful advising and contribution in manuscript revision.

Availability of data and materials

All data analyzed during this study are included in this article and its Additional files.

Authors’ contributions

WZ, DW and LY conceived and designed the research. WZ and LY performed the experiments and data analyses. WZ, DW, LJ and LY wrote the article; all authors read and approved the final article.

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Abbreviations

ABA: Abscisic Acid; ABRE: ABA responses; AMPK: AMP-activated protein kinase; AtSnRK: Arabidopsis thaliana SnRK; BnSnRK: Brassica napus SnRK; BLASTP: Basic local alignment search tool-protein; CBLs: calcium sensor calcineurin B-like proteins; CIPK: CBL-interacting protein kinases; DRE: dehydration responses element; GSDS: Gene structure display server; HMM: Hidden markov mode; LTRE: low-temperature responses element; ROS: reactive oxygen species; SAPK: stress-activated protein kinase; SnRK: sucrose non-fermenting 1 (SNF1)-related protein kinases; SOS: salt overly sensitive.

References

[1] Bohnert HJ, Gong QQ, Li PH, Ma SS. Unraveling abiotic stress tolerance mechanisms-getting genomics going. Current Opinion in Plant Biology. 9 (2006) 180-188.

[2] Hunter H. Protein kinases and phosphatases: the yin and yang of protein phosphorylation and signaling. Cell. 80 (1995) 225-236.

[3] Hrabak EM, Chan CW, Gribskov M, Harper JF, Choi JH, Halford N, Kudla J, Luan S, Nimmo HG, Sussman MR, Thomas M, Walker-Simmons K, Zhu JK, Harmon AC. The Arabidopsis CDPK-SnRK superfamily of protein kinases. Plant Physiol. 132 (2003) 666-680.

[4] Halford NG, Hardie DG. SNF1-related protein kinases: Global regulators of carbon metabolism in plants? Plant Molecular Biology. 37 (1998) 735-748.

[5] Celenza JL, Carlson M. A yeast gene that is essential for release from glucose repression encodes a protein kinase. Science. 233 (1986) 1175-1180.

[6] Kulik A, Wawer I, Krzywi ´nska E, Bucholc M, Dobrowolska G. SnRK2 Protein Kinases—Key Regulators of Plant Response to Abiotic Stresses. OMICS: A Journal of Integrative Biology. 15 (2011) 859-872.

[7] Albrecht V, Ritz O, Linder S, Harter K, Kudla J. The NAF domain defines a novel protein-protein interaction module conserved in Ca2+-regulated kinases. EMBO Journal. 20 (2011) 1051-1063.

[8] Masaru O, Yan Y, Ursula H, Zhu JK. A novel domain in the protein kinase SOS2 mediates interaction with the protein phosphatase 2C ABI2. Proc Natl Acad Sci USA. 100 (2003) 11771-11776.

[9] Purcell PC, Smith AM, Halford NG. Antisense expression of a sucrose non-fermenting-1-related protein kinase sequence in potato results in decreased expression of sucrose synthase in tubers and loss of sucrose-inducibility of sucrose synthase transcripts in leaves. Plant Journal. 14 (1998) 195-202.

[10] Ramon M, Dang TVT, Broeckx T, Hulsmans S, Crepin N, Sheen J, Rolland F. Default Activation and Nuclear Translocation of the Plant Cellular Energy Sensor SnRK1 Regulate Metabolic Stress Responses and Development. Plant Cell. 31 (2019) 1614-1632.

[11] Halford NG, Hey SJ. Snf1-related protein kinases (SnRKs) act within an intricate network that links metabolic and stress signalling in plants. Biochemical Journal. 419 (2009) 247-259.

[12] Maszkowska J, Dębski J, Kulik A, Kistowski M, Bucholc M, Lichocka M , Klimecka M, Sztatelman O, Szymanska KP, Dadlez M, Dobrowolska G. Phosphoproteomic analysis reveals that dehydrins ERD10 and ERD14 are phosphorylated by SNF1-related protein kinase 2.10 in response to osmotic stress. Plant Cell And Environment. 42 (2019) 931-946.

[13] Diédhiou CJ, Popova OV, Dietz KJ, Golldack D. The SNF1-type serine-threonine protein kinase SAPK4 regulates stress-responsive gene expression in rice. BMC Plant Biology. 8 (2008) 49.

[14] Zhang HY, Jia HF, Liu GS, Yang SN, ZHANG ST, Yang YX, Yang PP, Cui H. Cloning and characterization of SnRK2 subfamily II genes from Nicotiana tabacum. Molecular Biology Reports. 41 (2014) 5701-5709.

[15] Anna-Chiara M, Sylvain M, Alain V, Francesca F, Giraudat J. Arabidopsis OST1 Protein Kinase Mediates the Regulation of Stomatal Aperture by Abscisic Acid and Acts Upstream of Reactive Oxygen Species Production. Plant Cell. 14 (2002) 3089-3099.

[16] Ali A, Kim JK, Jan M, Khan HA, Khan IU, Shen MZ, Park J, Lim CJ, Hussain S, Bae, D, Wang K, Chung WS, Rubio V, Lee SY, Gong ZZ, Kim WY, Bressan RA, Pardo JM, Yun DJ. Rheostatic Control of ABA Signaling through HOS15-Mediated OST1 Degradation. Molecular Plant. 12 (2019) 1447–1462.

[17] Kim KN, Lee JS, Han H, Choi SA, Go SJ, Yoon IS. Isolation and characterization of a novel rice Ca2+-regulated protein kinase gene involved in responses to diverse signals including cold, light, cytokinins, sugars and salts. Plant Molecular Biology. 52 (2003) 1191-1202.

[18] Guo Y, Xiong L, Song CP, Gong D, Halfter U, Zhu JK. A calcium sensor and its interacting protein kinase are global regulators of abscisic acid signaling in Arabidopsis. Developmental Cell. 3 (2002) 233-244.

[19] Liu J, Ishitani M, Halfter U, Kim CS, Zhu JK. The Arabidopsis thaliana SOS2 gene encodes a protein kinase that is required for salt tolerance. Proc Natl Acad Sci USA. 97 (2000) 3730-3734.

[20] Ma QJ, Sun MH, Kang H, Lu J, You CX, Hao YJ. A CIPK protein kinase targets sucrose transporter MdSUT2.2 at Ser (254) for phosphorylation to enhance salt tolerance. Plant Cell And Environment. 42 (2019) 918-930.

[21] Ma QJ, Sun MH, Lu J, Kang H, You CX, Hao YJ. An apple sucrose transporter MdSUT2.2 is a phosphorylation target for protein kinase MdCIPK22 in response to drought. Plant Biotechnology Journal. 17 (2019) 625-637.

[22] Chen L, Ren F, Zhou L, Wang QQ, Zhong H, Li XB. The Brassica napus calcineurin B-like 1/CBL-interacting protein kinase 6 (CBL1/CIPK6) component is involved in the plant response to abiotic stress and ABA signalling. Journal of Experimental Botany. 63 (2012) 6211-6222.

[23] Chalhoub B, Denoeud F, Liu S, Parkin IAP, Tang H, Wang X, Chiquet J, Belcram H, Tong C, Samans B, et al. Early allopolyploid evolution in the post-Neolithic Brassica napus oilseed genome. Science 345 (2014) 950-953.

[24] Sun FM, Fan GY, Hu Q, Zhou YM, Guan M, Tong CB, Li JN. The high-quality genome of Brassica napus cultivar ‘ZS11’ reveals the introgression history in semi-winter morphotype. The Plant Journal. 92 (2017) 452-468.

[25] Zou J, Mao LF, Qiu J, Wang M, Lei J, Wu DY, He ZS, Chen MH, Shen YF, Shen EH, Huang YJ, Li RY, Hu DD, Shi L, Wang K, Zhu QH, Ye CY, Ian B, Graham JK, Meng JL, Fan LJ. Genome-wide selection footprints and deleterious variations in young Asian allotetraploid rapeseed. Plant Biotechnology Journal. 17 (2019) 1998-2010.

[26] Letunic I, Doerks T, Bork P. SMART 7: Recent updates to the protein domain annotation resource. Nucleic Acids Reserch. 40 (2012) 302-305.

[27] Marchler-Bauer A, Lu S, Anderson JB, Chitsaz F, Derbyshire MK, DeWeese-Scott C, Fong JH, Geer LY, Geer RC, Gonzales NR, Gwadz M, Hurwitz DI, Jackson JD, Ke Z, Lanczycki CJ, Lu F, Marchler GH, Mullokandov M, Omelchenko MV, Robertson CL,Song JS, Thanki N, Yamashita RA, Zhang D, Zhang N, Zheng C, Bryant SH. CDD: a Conserved Domain Database for the functional annotation of proteins. Nucleic Acids Research. 39 (2011) 225-229.

[28] Finn RD, Coggill P, Eberhardt RY, Eddy SR, Mistry J, Mitchell AL, Potter SC, Punta M, Qureshi M, Sangradorvegas A. The Pfam protein families database: Towards a more sustainable future. Nucleic Acids Research. 44 (2016) 279-285.

[29] Saitou N, Nei M. The neighbor-joining method: A new method for reconstructing phylogenetic trees. Molecular Biology and Evolution. 4 (1987) 406-425.

[30] Kumar S, Stecher G, Tamura K. MEGA7: Molecular evolutionary genetics analysis version 7.0 for bigger datasets. Molecular Biology and Evolution. 33 (2016) 1870-1874.

[31] Bailey TL, Boden M, Buske FA, Frith M, Grant CE, Clementi L, Ren J, Li WW, Noble WS. MEME SUITE: Tools for motif discovery and searching. Nucleic Acids Res. 37 (2009) 202-208.

[32] Hu B, Jin J, Guo AY, Zhang H, Luo J, Gao G. GSDS 2.0: An upgraded gene feature visualization server. Bioinformatics. 31 (2015) 1296-1297.

[33] Voorrips RE. MapChart: Software for the graphical presentation of linkage maps and QTLs. Journal of Heredity. 93 (2002) 77-78.

[34] Krzywinski M, Schein J, Birol I, Connors J, Gascoyne R, Horsman D, Jones SJ,Marra MA. Circos: an information aesthetic for comparative genomics. Genome Research. 19 (2009) 1639–1645.

[35] Wang Y, Tang H, DeBarry JD, Tan X, Li J, Wang X, Lee T, Jin H, Marler B, Guo H. MCScanX: a toolkit for detection and evolutionary analysis of gene synteny and collinearity. Nucleic Acids Research. 40 (2012) 49.

[36] Wang D, Zhang Z, Zhang Z, Zhu J, Yu J. KaKs_Calculator 2.0: a toolkit incorporating gamma-series methods and sliding window strategies. Genomics Proteomics Bioinformatics. 8 (2010) 77-80.

[37] Lescot M, Déhais P, Thijs G, Marchal K, Moreau Y, Yves VDP, Pieree R, Stephane R. Plantcare, a database of plant cis-acting regulatory elements and a portal to tools for in silico analysis of promoter sequences. Nucleic Acids Research. 30 (2002) 325-7.

[38] Zhang YT, Ali U, Zhang GF, Yu LQ, Fang S, Iqbal S, Li HH, Lu SP, Guo L. Transcriptome analysis reveals genes commonly responding to multiple abiotic stresses in rapeseed. Molecular Breeding. 39 (2019).

[39] Chen CJ, Xia R, Chen H, He YH. TBtools, a Toolkit for Biologists integrating various HTS data handling tools with a user-friendly interface. bioRXiv. (2018).

[40] Wu DZ, Liang Z, Yan T, Xu Y, Xuan LJ, Tang J, Zhou G, Lohwasser U, Hua SJ, Wang HY, Chen XY, Wang Q, Zhu L, Maodzeka A, Hussain N, Li ZL, Li XM , Shamsi IH, Jilani G, Wu LD, Zheng HK, Zhang GP, Chalhoub B, Shen LS, Yu H, Jiang LX. Whole-Genome Resequencing of a Worldwide Collection of Rapeseed Accessions Reveals the Genetic Basis of Ecotype Divergence. Molecular Plant. 12 (2019) 30-43.

[41] Xuan LJ, Yan T, Lu LZ, Zhao XZ, Wu DZ, Hua SJ, Jiang LX. Genome-wide association study reveals new genes involved in leaf trichome formation in polyploid oilseed rape (Brassica napus L.). Plant Cell And Environment. 2019.

[42] Holub EB. The arms race is ancient history in Arabidopsis, the wildflower. Nat Rev Genet. 2 (2001) 516.

[43] Kobayashi Y, Yamamoto S, Minami H, Kagaya Y, Hattori T. Differential activation of the rice sucrose nonfermenting1-related protein kinase2 family by hyperosmotic stress and abscisic acid. Plant Cell. 16 (2004) 1163-1177.

[44] Wang LZ, Hu W, Sun JT, Liang XY, Yang XY, Wei SY, Wang XT, Zhou Y, Xiao Q, Yang GX, He GY. Genome-wide analysis of SnRK gene family in Brachypodium distachyon and functional characterization of BdSnRK2.9. Plant Science. 237 (2015) 33-45.

[45] Wang YJ, Yan HF, Qiu ZF, Hu B, Zeng BS, Zhong CL, Fan CJ. Comprehensive Analysis of SnRK Gene Family and their Responses to Salt Stress in Eucalyptus grandis. International Journal of Molecular Sciences. 20 (2019) 11.

[46] Chung BY, Simons C, Firth AE, Brown CM, Hellens RP. Effect of 5’utr intronson gene expression in Arabidopsis thaliana. BMC Genomics. 7 (2006) 120.

[47] Jeffares DC, Penkett CJ, Bähler J. Rapidly regulated genes are intronpoor. Trends Genet. 24 (2008) 375.

[48] Stephen J, Fei N, Juergen E, Zhang Stephen J, Fei N, Juergen E, Zhang S, Dong W, Xue T, Zheng C, Yuan Z. Genome-wide and expression analysis of protein phosphatase 2C in rice and Arabidopsis. BMC Genomics. 9 (2008) 1-21.

[49] Ali GM, Komatsu S. Proteomic analysis of rice leaf sheath during drought stress. Journal of Proteome Research. 5 (2006) 396-403.

[50] Sarda X, Tousch D, Ferrare K, Legrand E, Dupuis J, Casse-Delbart F, Lamaze T. Two TIP-like genes encoding aquaporins are expressed in sunflower guard cells. Plant Journal. 12 (1997) 1103-1111.

[51] Tan WR, Zhang DW, Zhou HP, Zheng T, Yin YH, Lin HH. Transcription factor HAT1 is a substrate of SnRK2.3 kinase and negatively regulates ABA synthesis and signaling in Arabidopsis responding to drought. PLOS Genetics. 14 (2018) .

[52] Nakashima K, Fujita Y, Kanamori N, Katagiri T, Umezawa T, Kidokoro S, Maruyama K, Yoshida T, Ishiyama K, Kobayashi M, Shinozaki K, Yamaguchi-Shinozaki K. Three Arabidopsis SnRK2 Protein Kinases, SRK2D/SnRK2.2, SRK2E/SnRK2.6/OST1 and SRK2I/SnRK2.3, Involved in ABA Signaling are Essential for the Control of Seed Development and Dormancy. Plant and Cell Physiology. 50 (2009) 1345-1363.

[53] Yoshida T, Fujita Y, Maruyama K, Mogami J, Todaka D, Shinozaki K, Yamaguchi-Shinozaki K. Four Arabidopsis AREB/ABF transcription factors function predominantly in gene expression downstream of SnRK2 kinases in abscisic acid signalling in response to osmotic stress. Plant Cell and Environment. 38 (2015) 35-49.