WRKY gene family of Licorice (Glycyrrhiza uralensis): identification and expression analysis in response to abiotic stress

Abstract

WRKY transcription factors (TF) control a variety of biological processes in higher plants, including abiotic and biotic stress responses. WRKY genes have been identified and functionally characterized in many plant species during the past ten years, but despite this progress, few investigations have focused on this TF family in licorice (Glycyrrhiza uralensis Fisch.). This study identified 72 WRKY genes from the licorice genome and classified them into three groups. The gene structure analysis indicated that the GuWRKY genes contain from one to eight introns, with the majority having two introns and three exons. Conserved motif analysis identified the conserved WRKY domain in most GuWRKYs. Further analysis of the promoter region revealed that almost all GuWRKYs had plant hormone- or stress-responsive cis-regulatory elements, with differences in the number of elements among the groups. The role of 10 selected GuWRKYs in regulating abiotic stress responses was further explored via quantitative real-time PCR (qRT-PCR) after low-temperature (4 °C), salt (150 mmol·L⁻¹ NaCl), drought (15% polyethylene glycol 6000) treatments at 3, 6, 12, and 24 h. Thus, our study of the WRKY gene family in licorice provides a basis for further research on the structure and function of WRKY genes.

Keywords:

• Glycyrrhiza uralensis
• WRKY transcription factors
• gene family
• bioinformatic analysis
• abiotic stress

Introduction

In plants, the WRKY family is one of the largest families of transcription factors (TF). Their name is derived from the most conserved and characteristic WRKY domain in the protein sequence [1]. In recent years, WRKY TFs have attracted research interest because of the numerous members of each species, diverse functions in development, and their application in phylogenetics. The WRKY domain, which is 60 amino acids long and comprises a conserved WRKYGQK sequence at the N-terminus and a C₂H₂ or C₂HC zinc motif at the C-terminus, is unique to WRKY proteins [1, 2]. The WRKYGQK sequence and the zinc-finger motif are both essential for WRKY TFs to bind to the consensus W-box (TTGACC/T) cis-elements in their target genes [3, 4]. WRKY proteins are classified into three major groups based on the number of WRKY domains and the structure of the zinc-finger motifs. Two WRKY conserved domains, including a C₂H₂ motif, are found in Group I WRKY proteins. WRKY proteins containing a single WRKY domain and a C₂H₂ motif are members of Group II. A few WRKY proteins feature a WRKY domain with a finger motif that is distinct from the members of groups I and II. They are categorized as group III and have a Cx7Cx23HxC zinc-finger domain rather than a Cx4-5Cx22-23HxH pattern [5–7].

The first WRKY gene (SPF1) was cloned from sweet potato (Ipomoea batatas) [8], following which many were identified in Arabidopsis (Arabidopsis thaliana) [9], rice (Oryza sativa L.) [10], Tamarix (Tamarix hispida) [11], mulberry (Morus notabilis) [12], oil palm (Elaeis guineensis) [13], cucumber (Cucumis sativus L.) [14], sorghum (Sorghum bicolor (L.) Moench) [15] and various other plant species. WRKY TFs have been shown in studies to have an important role in biotic and abiotic stress response and development. WRKY8 controls salt tolerance in Arabidopsis [16], while other WRKY genes (WRKY1, WRKY46, WRKY54, and WRKY63) contribute in drought stress response via distinct signaling pathways [17–19]. WRKY4 enhances tolerance to salt and abscisic acid (ABA) treatment in T. hispida by increasing superoxide dismutase activity, decreasing chlorophyll loss and protecting cells [11]. In Polygonum cuspidatum, WRKY33 negatively regulates salt tolerance by downregulating stress-related genes and increasing reactive oxygen species (ROS) [20]. In Oryza sativa, OsWRKY114 increases drought sensitivity by inhibiting stomatal closure [21]. Overexpression studies in Arabidopsis demonstrated that the maize WRKY ZmWRKY17 controls the transcription of stress- and ABA-related genes, reducing plant tolerance to salt stress and ABA sensitivity [22], while IgWRKY50 and IgWRKY32 of Iris germanica enhance drought resistance by mediating the ABA signaling pathway [23].

Licorice (Glycyrrhiza uralensis Fisch.) is a perennial legume widely cultivated as an edible and medicinal crop. Its area of distribution is mainly in the northwest arid, semi-arid and desert regions of China, where the average low temperature and high temperature in winter are very low [24]. Licorice has high-stress tolerance and plays significant roles in desertification control, animal husbandry and human health [25, 26]. Unfortunately, most cultivated licorice roots do not have quality similar to the wild variety; therefore, these materials are often diverted to non-medicinal, food, or confectionary uses. Therefore, we urgently need to improve the yield and quality traits of licorice via breeding, which demands an in-depth understanding of the growth and response of this species under adverse environmental factors. Mochida et al. [27] published the draft genome of licorice and a chromosome-scale genome assembly was reported for licorice [28]. The genome sequencing and annotation provide an essential resource for the gene functional annotation and evolutionary analysis of WRKY family in licorice.

It is crucial to thoroughly examine the WRKY family in order to extend the theoretical investigation and economic utilization of licorice. The present study analyzed the phylogenetic relationships, gene and protein structures, and promoter cis-elements of 72 GuWRKYs of licorice. Further, we treated licorice seedlings with sodium chloride (NaCl), polyethylene glycol (PEG6000), and low-temperature stress and analyzed the transcription levels of 10 GuWRKY genes by quantitative real-time PCR (qRT-PCR) to evaluate the response to abiotic stress. Thus, this study functionally characterizes novel WRKY family members of licorice, facilitating the development of stress-resistant cultivars and increase in yield and quality.

Materials and methods

Identification of the GuWRKY family

Based on the licorice genome sequence published by Mochida et al. [27], the raw data of the genome sequence was accessed from the DNA Database of Japan (DDBJ; ID: PRJDB3943). The sequences of the AtWRKY proteins were downloaded from the Arabidopsis Information Resource (TAIR; https://www.arabidopsis.org) (Supplemental Table S1) [29]. Then, the AtWRKY sequence was used as the query in BLAST, and the GuWRKY members were identified based on homology. The Pfam database (http://pfam.xfam.org/) [30] was used to predict the conserved domains of the obtained sequences, and the genes of the GuWRKY family were finally selected by removing those without the characteristic conserved domain. The molecular weight (MW), isoelectric point (pI), instability index (II), aliphatic index (AI) and grand average of hydropathicity (GRAVY) of the GuWRKY proteins were determined using the ProtParam tool in the online ExPASy software (https://web.expasy.org/protparam/) [31]. The subcellular localization of these proteins was predicted using the online tool Cell-PLoc 2.0 (http://www.csbio.sjtu.edu.cn/bioinf/plant-multi/) [32].

Phylogenetic analysis

An interspecific phylogenetic tree was constructed following the maximum likelihood (ML) method in MEGA7.0 with 1000 ultrafast bootstrap replications [33].

Conserved motifs, gene structure analysis and cis-regulatory element identification

The exon/intron structures of GuWRKY genes were visualized using the online tool Gene Structure Display Server 2.0 (GSDS 2.0; http://gsds.cbi.pku.edu.cn/) [34] to clarify the structural diversity. The conserved motifs in the GuWRKY protein sequences were analyzed using the MEME online tool with the following parameters: maximum number of motifs, 15 (http://meme-suite.org/index.html), and the final result was displayed using TBtools (v1.098726) [35]. The PlantCARE (http://bioinformatics.psb.ugent.be/webtools/plantcare/html/) database was used to predict the cis-regulatory elements 2.0 kb upstream of the start codon of each gene.

Plant materials and stress treatment

The seeds of licorice (G. uralensis Fisch.) purchased from Urumqi County, Xinjiang, China (43°47′N, 87°41′E). The seeds used in this study had been formally identified by Associate Professor Wenjing Huang (Shaanxi University of Chinese Medicine). Dry seeds were stored in labeled envelopes in a refrigerator (4 °C) in Co-construction Collaborative Innovation Center for Chinese Medicine Resources Industrialization by Shaanxi & Education Ministry. Seeds were treated with 98% sulfuric acid for 30 min (to break the seed dormancy), rinsed five times with sterilized ultrapure water, placed in a petri dish with two layers of moist filter paper, and incubated in an artificial climate box (25 °C, 16 h light/8 h dark cycle, and 60%−70% relative humidity) for seven days. Then, seedlings with consistent growth were selected, transplanted into plastic hydroponic pot containers filled with 100 g of perlite. There were 15 seedlings in each pot, a total of 12 pots were transplanted. Each pot was irrigated with 250 mL of Hoagland’s nutrient solution, and supplemented once every two days. All plants were grown in a greenhouse under a 12 h light period, 25/20 °C day/night temperatures, 60%–70% relative humidity, and 3000 lx light intensity, and maintained for one month for subsequent experiments. Low temperature treatment: the plants were placed in an artificial climate box chamber at 4 °C (16 h light/8 h dark cycle, and 60%−70% relative humidity). Salt stress treatment: the plants were treated with 250 mL Hoagland’s nutrient solution containing 150 mmol·L⁻¹ NaCl per pot. Drought stress treatment: the plants were treated with 250 mL Hoagland’s nutrient solution containing 15% PEG-6000 (w/v) per pot. Control (CK): the plants received 250 mL Hoagland’s nutrient solution per pot. Three pots were used for each treatment. Plant samples, including stem, roots and leaves, were collected at 3, 6, 12, and 24 h after the beginning of the treatment. Three plants were collected per time point per treatment. Then, all roots, stems and leaves were frozen in liquid nitrogen immediately and stored at −80 °C until used for the RNA extraction. All experiments consisted of three biological replicates.

RNA isolation, cDNA synthesis and quantitative real‑time PCR analysis

The roots, stems and leaves of licorice were quickly transferred to a mortar pre-cooled with liquid nitrogen. The tissue samples were ground into powder using a pestle and mortar while being continuously sprayed with liquid nitrogen. About 0.1 g of each ground sample was taken and total RNA was extracted from the licorice roots, stems and leaves using the RNAprep Pure Plant Kit (TIANGEN, Beijing, China) and treated with RNase-Free DNase I (TIANGEN, Beijing, China) to remove DNA contamination. The quality of the extracted RNA was assessed by agarose gel (0.8%) electrophoresis, and the concentration and purity (OD₂₆₀/OD₂₈₀) were measured using NanoDrop One (Thermo Fisher Scientific, Inc., Waltham, MA, USA). Using a FastKing RT Kit (TIANGEN, Beijing, China) and a T100 Thermal Cycler (Bio-Rad, Inc., California, USA), approximately 1 μg of total RNA was reversely transcribed into first-strand cDNA by FastKing RT Enzyme. The reaction was carried out as follows: incubate at 42 °C for 15 min and then at 95 °C for 3 min.

The CDS sequences of the GuWRKY gene family members were obtained from the whole genome database of licorice (DDBJ; ID: PRJDB3943). The primers for amplification of the GuWRKYs were designed using the Primer Premier 5.0 software and synthesized by Shanghai Shenggong Bioengineering Co., Ltd (https://www.sangon.com/) (Table 1). Then, qRT-PCR was performed to determine the relative expression levels of the genes with a SuperReal PreMix Plus (SYBR Green) Kit (TIANGEN, Beijing, China) on a qTOWER2.0 real-time PCR system (Germany). The PCR mixture contained 10 μL of the SYBR Green PCR Master Mix (Without ROX), 0.6 μL of each primer, 1 μL of cDNA template (0.3 µmol·L⁻¹) and 8.2 μL of ddH₂O. The reaction was carried out as follows: pre-denaturation at 95 °C for 15 min, followed by 40 cycles of denaturation at 95 °C for 10 s, annealing at 60 °C for 30 s, extension at 72 °C for 30 s. Actin was used as the internal reference gene for specific expression analysis [36], and the assay was repeated three times. Finally, the fluorescence melting curve was analyzed. The relative expression values of the genes were calculated based on three biological replicates using the 2^−ΔΔCT method [37]. The single factor variance analysis of the real-time PCR data was performed using the SPSS 26 software.

Table 1. Primer sequences for qRT-PCR. The primers for GuWRKY amplification were designed using Primer Premier 5.0 software. The actin primer pair was used according to maroufi [36].

Gene name	Forward primer (5'→3')	Reverse primer (5'→3')
GuWRKY5	GAAGCGGCAGCATCAGAGGAAG	CTCCAGCGAAAGCCATCTCCAAG
GuWRKY6	GTCAGTGAAGAACAGTCCCAACCC	GCTTGAGTCTTCCCTTCCCTTTCC
GuWRKY19	TGTCTCCAAGTTCAAGCAGGTCATC	ACTGTGTTTGAGGTTGAGGTTGAGG
GuWRKY23	GACAAGGAATCCGAACCCGAACC	CGCCCAACCATCCCAACTCATC
GuWRKY36	GAGCAAGACTAGAGTTTCCCACAGC	GAGCAAGACTAGAGTTTCCCACAGC
GuWRKY41	ATTAGTAGCAGTAGCAGCAGCAACC	TCACAAGCAGGGTCGTTTAAGAAGG
GuWRKY42	CGAGGACGACAACAACAACAACAAC	CCTAAGTGGCATGGTGACAAGCAG
GuWRKY53	GGGCAATAAGAGGAAGCGAGTCATC	CCAACAGCAGCAGCAGTAGGATTC
GuWRKY70	TCACCACCAGAGTCCACAGTATCC	AGCCCTTGGTTTTGCGGAAGATG
GuWRKY72	GTGCTTGCCCATCATCTCCTCTG	CATGTTCTGCCCCACCCCAAAG
Actin	CCTCTCTCTTTATGCCAGTGGT	AGGGAGTGCATAACCCTCATAG

Statistical analysis

Data are represented as mean values with standard error of the mean (±SEM). All experimental data were analyzed using SPSS 26.0 software (SPSS Inc., Chicago, IL, United States). Significant differences between the treatment means were determined using the least significant difference (LSD; p < 0.05) test. All graphs were plotted using GraphPad Prism 8.0.2.263 and Origin 2021.

Results

Identification and classification of WRKY genes

The study identified 72 GuWRKY protein-coding genes from the licorice genome based on the presence of a conserved domain and designated them from GuWRKY01 to GuWRKY72. Further, the phylogenetic tree reflected the evolutionary status and grouped the 72 members of the WRKY family (Figure 1). The 72 identified members of the licorice WRKY family were sorted into three major groups, with the third group having the largest number of members (40); the first and second groups had 15 and 17 members, respectively.

Figure 1. Phylogenetic relationships among the WRKY genes of licorice. A maximum likelihood phylogenetic tree was constructed in MEGA7.0 with bootstrap testing (1000 replicates). the tree divided the WRKY genes into groups I, II and III.

Subsequently, we predicted the physicochemical properties of the GuWRKY proteins (Table 2). The analysis identified GuWRKY66 as the shortest protein [73 amino acids (aa)] and GuWRKY44 as the longest (722 aa). The predicted molecular weight (MW) of the GuWRKY proteins ranged from 8.405 kDa (GuWRKY66) to 83.685 kDa (GuWRKY44), and the pI from 4.44 (GuWRKY6) to 9.97 (GuWRKY19). The negative GRAVY values indicated that all the GuWRKY proteins are hydrophilic. Most GuWRKYs were found to be unstable, with the maximum instability index for GuWRKY2 (71.99) and the minimum for GuWRKY61 (35.74). Meanwhile, the Aliphatic index had a maximum of 86.71 (GuWRKY66) and a minimum of 38.95 (GuWRKY63). Additionally, Cell-PLoc predicted nucleus localization for the 72 GuWRKY proteins.

Table 2. Information of GuWRKY genes. The number of amino acids, molecular weight, theoretical isoelectric point, instability index, aliphatic index and grand average of hydropathicity of the GuWRKY proteins were determined using the ProtParam tool in the online ExPASy software [31]. The subcellular localization of these proteins was predicted using the online tool Cell-PLoc 2.0 [32].

Protein name	ID	Number of amino acids (aa)	Molecular weight (kDa)	Theoretical isoelectric point (pI)	Instability index	Aliphatic index	Grand average of hydropathicity	Predicted location(s)
GuWRKY1	Glyur002923s00038431.1	409	44.357	6.59	54.86	51.59	−0.733	Nucleus
GuWRKY2	Glyur000049s00003147.1	293	32.342	7.08	71.99	52.87	−0.74	Nucleus
GuWRKY3	Glyur000190s00015737.1	473	52.326	5.72	60.61	53.4	−0.777	Nucleus
GuWRKY4	Glyur000831s00034077.1	227	25.675	9.8	56.07	60.93	−0.782	Nucleus
GuWRKY5	Glyur000081s00006738.1	424	46.543	8.78	49.28	54.53	−0.881	Nucleus
GuWRKY6	Glyur000999s00024092.1	193	22.322	4.44	47.67	59.53	−0.879	Nucleus
GuWRKY7	Glyur000515s00023231.1	642	68.579	7.75	47.55	58.58	−0.667	Nucleus
GuWRKY8	Glyur000020s00001774.1	613	65.927	6.84	54.28	61.04	−0.673	Nucleus
GuWRKY9	Glyur000218s00011611.1	447	49.48	6.56	50.87	61.21	−0.952	Nucleus
GuWRKY10	Glyur000111s00011064.1	512	55.247	5.94	43.88	64.47	−0.542	Nucleus
GuWRKY11	Glyur000315s00012756.1	500	55.269	7.37	49.72	70.26	−0.59	Nucleus
GuWRKY12	Glyur000068s00006850.1	628	67.904	5.35	46.07	58.3	−0.719	Nucleus
GuWRKY13	Glyur000101s00010298.1	293	32.073	9.67	43.24	68.29	−0.56	Nucleus
GuWRKY14	Glyur002655s00040012.1	468	51.525	8.55	41.58	64.81	−0.671	Nucleus
GuWRKY15	Glyur000097s00008318.1	342	36.929	9.58	47.31	69.09	−0.439	Nucleus
GuWRKY16	Glyur005236s00044304.1	319	35.465	5.41	48.57	61.44	−0.786	Nucleus
GuWRKY17	Glyur000498s00016149.1	448	49.786	6.75	57.86	65.6	−0.754	Nucleus
GuWRKY18	Glyur000016s00003035.1	348	39.226	9.73	50.84	62.21	−0.801	Nucleus
GuWRKY19	Glyur000102s00010972.1	322	35.119	9.97	51.37	63.94	−0.619	Nucleus
GuWRKY20	Glyur000210s00014438.1	364	40.696	9.7	55.26	66.15	−0.766	Nucleus
GuWRKY21	Glyur002240s00040736.1	301	33.063	9.64	45.45	65.85	−0.625	Nucleus
GuWRKY22	Glyur000431s00018649.1	633	68.262	6.14	54.14	61.63	−0.702	Nucleus
GuWRKY23	Glyur003077s00042521.1	215	23.667	5.65	67.66	62.6	−0.726	Nucleus
GuWRKY24	Glyur000821s00020788.1	540	59.489	5.56	47.81	58.76	−0.876	Nucleus
GuWRKY25	Glyur000012s00000551.1	426	46.754	8.78	54.98	75.87	−0.47	Nucleus
GuWRKY26	Glyur000049s00003157.1	267	30.454	4.87	67.56	59.14	−0.98	Nucleus
GuWRKY27	Glyur000802s00023733.1	263	28.693	5.33	62.42	55.67	−0.711	Nucleus
GuWRKY28	Glyur000464s00022410.1	352	38.479	9.68	57.83	61.59	−0.753	Nucleus
GuWRKY29	Glyur000420s00021423.1	257	28.57	5.62	62.27	55.8	−0.755	Nucleus
GuWRKY30	Glyur000048s00002242.1	174	20.087	5.66	54.2	44.77	−1.018	Nucleus
GuWRKY31	Glyur001574s00034448.1	376	41.338	9.27	50.94	63.54	−0.614	Nucleus
GuWRKY32	Glyur000392s00023197.1	399	43.909	7.2	56.45	48.47	−1.019	Nucleus
GuWRKY33	Glyur000243s00013151.1	354	38.476	5.58	53.69	52.97	−0.714	Nucleus
GuWRKY34	Glyur000036s00004638.1	355	40.57	5.37	48.53	57.13	−0.824	Nucleus
GuWRKY35	Glyur000172s00014880.1	394	42.991	5.47	51.18	55.71	−0.613	Nucleus
GuWRKY36	Glyur000672s00022927.1	522	58.525	6.23	50.29	81.03	−0.568	Nucleus
GuWRKY37	Glyur000121s00009710.1	355	40.056	5.39	47.58	59.89	−0.719	Nucleus
GuWRKY38	Glyur000162s00009909.1	287	31.609	8.45	46.58	66.62	−0.597	Nucleus
GuWRKY39	Glyur000381s00016114.1	381	41.342	5.32	60.44	58.06	−0.729	Nucleus
GuWRKY40	Glyur001658s00045059.1	447	49.388	7.03	48.81	70.11	−0.576	Nucleus
GuWRKY41	Glyur000132s00015923.1	321	36.193	5.44	68.16	71.06	−0.774	Nucleus
GuWRKY42	Glyur000194s00016333.1	352	39.688	5.33	54.24	64.8	−0.669	Nucleus
GuWRKY43	Glyur001230s00027566.1	339	36.58	5.19	60.78	50.09	−0.811	Nucleus
GuWRKY44	Glyur000018s00003729.1	722	83.685	5.93	52.82	61.47	−0.709	Nucleus
GuWRKY45	Glyur000603s00029236.1	249	26.993	5.14	63.57	48.31	−0.717	Nucleus
GuWRKY46	Glyur001238s00031501.1	360	39.589	5.71	48.06	61.47	−0.59	Nucleus
GuWRKY47	Glyur007003s00046206.1	341	38.78	9.33	48.35	54.34	−0.996	Nucleus
GuWRKY48	Glyur000870s00029783.1	433	45.959	5.38	48.15	55.7	−0.667	Nucleus
GuWRKY49	Glyur000243s00013147.1	193	21.727	5.56	42.7	60.05	−0.905	Nucleus
GuWRKY50	Glyur001338s00034122.1	161	17.866	4.7	46.07	49.63	−1.106	Nucleus
GuWRKY51	Glyur002123s00028291.1	380	42.56	6.8	44.98	61.42	−0.669	Nucleus
GuWRKY52	Glyur000521s00021483.1	312	35.263	5.89	52.57	43.75	−1.065	Nucleus
GuWRKY53	Glyur002177s00034179.1	341	39.194	6.67	53.74	56.63	−0.948	Nucleus
GuWRKY54	Glyur001238s00031502.2	199	22.392	8.84	48.23	55.93	−0.846	Nucleus
GuWRKY55	Glyur001238s00031502.1	492	55.118	6.75	52.68	56.91	−0.828	Nucleus
GuWRKY56	Glyur000040s00007334.1	311	35.834	6.56	54.9	73.02	−0.643	Nucleus
GuWRKY57	Glyur000151s00010735.1	155	17.267	5.24	50.69	57.94	−0.577	Nucleus
GuWRKY58	Glyur000092s00007537.2	377	40.262	6.3	58.86	50.19	−0.749	Nucleus
GuWRKY59	Glyur000092s00007537.1	401	42.985	6.34	60.75	48.4	−0.764	Nucleus
GuWRKY60	Glyur002177s00034178.1	261	29.493	6.24	70.44	73.22	−0.538	Nucleus
GuWRKY61	Glyur000105s00014183.1	127	14.549	4.52	35.74	45.2	−1.131	Nucleus
GuWRKY62	Glyur000085s00008894.1	579	62.925	4.88	53.49	66.55	−0.583	Nucleus
GuWRKY63	Glyur005843s00044722.1	95	11.103	9.86	50.31	38.95	−1.106	Nucleus
GuWRKY64	Glyur001814s00036625.1	372	40.684	5.91	56.93	53.95	−1.048	Nucleus
GuWRKY65	Glyur000144s00012202.1	265	29.145	6.17	64.71	59.51	−0.804	Nucleus
GuWRKY66	Glyur008961s00047170.1	73	8.405	6.89	39.95	86.71	−0.178	Nucleus
GuWRKY67	Glyur000917s00030711.1	359	39.525	6.45	53.45	55.96	−0.849	Nucleus
GuWRKY68	Glyur000183s00013026.1	448	49.098	8.69	51.44	50.78	−0.905	Nucleus
GuWRKY69	Glyur000005s00001082.1	443	49.891	9.89	59.97	41.6	−1.374	Nucleus
GuWRKY70	Glyur000577s00027226.1	309	35.041	5.66	51.07	62.52	−0.732	Nucleus
GuWRKY71	Glyur000392s00023192.1	297	34.341	5.45	63.15	74.78	−0.648	Nucleus
GuWRKY72	Glyur000858s00035416.1	236	25.838	4.87	59.61	68.18	−0.465	Nucleus

Gene structure and motif composition of GuWRKYs

We further mapped the genetic structure of GuWRKYs to gain in-depth knowledge on the evolution of the WRKY gene family in licorice. Gene structure analysis was performed to construct a detailed exon-intron map (Figure 2). Among the 72 GuWRKYs, GuWRKY49, GuWRKY61, GuWRKY66, GuWRKY70 and GuWRKY72 had the least introns (one), while GuWRKY44 had the highest number (eight). Among the remaining GuWRKY genes, 32 had 2 introns, and another 5, 12, 10, 5 and 2 genes had 3, 4, 5, 6 and 7 introns, respectively. Meanwhile, nearly half of the genes (32) had three exons. These observations indicated that exon loss and gain occurred in GuWRKY genes during evolution, suggesting functional diversity among the genes.

Figure 2. Gene structure of GuWRKY genes. The exon/intron structures of GuWRKY genes were visualized using gene structure Display Server 2.0 (GSDS 2.0) [34]. the yellow and green boxes represent the coding sequences (CDS) and the untranslated regions (UTR), respectively. The black lines indicate introns.

Then, we analyzed the motifs of the 72 GuWRKY proteins using the MEME program to assess the conservation and diversification among the genes. A total of 10 conserved motifs were identified, namely, motifs 1–10. The size of the conserved motifs ranged from 11 (motif 3) to 45 (motif 5) amino acids (Figure 3). Among these, motif 1 was identified as the conserved WRKY domain and contained a WRKYGQK sequence, and motif 2 as the conserved zinc-finger structure. Additionally, other motifs (motifs 3–10) were also predicted within the GuWRKYs. The number of conserved motifs in each GuWRKY protein varied from 2 to 10. As shown in Figure 3, motifs 1, 2 and 3 were present in most genes, motif 8 was unique to group I, motif 9 existed only in group II, and motifs 5, 6, 7 and 10 were present in group III. Thus, our observations suggested that the GuWRKY members of each group had similar motif composition, consistent with the clustering based on phylogenetic analysis.

Figure 3. Conserved motifs of GuWRKY proteins of groups I, II and III. MEME predicted a total of ten motifs shown in different colors. Gray lines represent non-conserved sequences.

Analysis of cis-acting elements in GuWRKY promoters

Using the PlantCARE website, the cis-acting elements were found in the 2.0 kb upstream sequence from the translation start site (ATG) of each GuWRKY gene in order to infer the functions of the GuWRKYs. The analysis identified 25 cis-acting elements, including phytohormone-responsive, defense- and stress-responsive, development-related and light-responsive elements, in the promoter regions of all GuWRKY genes (Figure 4). The promoter regions contained numerous phytohormone-responsive cis-acting elements, such as those that responded to abscisic acid (ABRE), auxin (AuxRR-core and TGA-element), gibberellin (GARE-motif, TATC-box and P-box), MeJA (CGTCA-motif and TGACG-motif) and salicylic acid (TCA-element). The ABA and MeJA (ABRE, CGTCA-motif and TGACG-motif) comprised the majority of these components. In addition, five defense- and stress-responsive cis-acting elements were detected in the WRKY promoter regions, containing the TTGACC sequence (W-box), the low-temperature-responsive (LTR), defense- and stress-responsive (TC-rich repeats), drought-inducible (MBS) and wound-responsive elements (WUN-motif). The analysis detected LTR in 19 (26%), TC-rich repeat elements in 23 (32%), MBS in 27 (38%), WUN-motif in 24 (33%) and W-box in 47 (65%) WRKY gene promoters. Also, we discovered that 48 (67%) of the gene promoters had two or more stress-related elements, which suggests a role for GuWRKYs in a variety of stress responses. Endosperm expression-specific elements (GCN4_motif), anaerobic induction elements (ARE), meristem expression-specific elements (CAT-boxes), seed-specific regulatory element (RY-element), circadian control-related elements (circadian) and flavonoid biosynthetic gene-regulating MYB binding site (MBSI) are cis-acting elements related to plant development. In licorice, MBSI alone existed in GuWRKY53, 54 and 55, suggesting the role of these three genes in regulating flavonoid metabolism. In addition, our analysis showed the presence of light-responsive elements, such as GT1-motif, G-box, Box 4, AE-box and GATA-motif, essential for light induction, in most GuWRKY genes.

Figure 4. Cis-acting elements of GuWRKY genes. The PlantCARE database was used to predict the cis-regulatory elements 2.0 kb upstream of the start codon of each gene. Here, the x-axis shows the type of cis-acting elements, and the y-axis indicates the number of WRKY family members.

Expression analysis of GuWRKYs under abiotic stress

Finally, to explore the role of GuWRKY genes in drought, low temperature and salt stress responses, we analyzed the expression of ten genes specifically in leaf, stem and root tissues of licorice by qRT-PCR. We found that most GuWRKY genes responded to at least one treatment. In the roots (Figure 5), the expression levels of eight genes (GuWRKY5/6/19/23/42/53/70/72) were significantly upregulated after NaCl treatment compared to control (CK), but were downregulated considerably or remained unchanged at 12 and 24 h after cold treatment. GuWRKY6/41/53 were upregulated at 12 h after low-temperature treatment. A few GuWRKY genes were upregulated only at a specific time; GuWRKY36/41 were significantly upregulated at 12 h after NaCl treatment, GuWRKY19/72 at 12 h after PEG treatment, and GuWRKY42 at 24 h after PEG treatment. In the stem (Figure 6), GuWRKY19/23 and GuWRKY19/70 were significantly upregulated at 3 and 6 h after NaCl treatment and PEG treatment, respectively, compared to the levels in the non-treated control plants (CK). In addition, GuWRKY5/42/53/70/72 were significantly upregulated at 12 and 24 h after NaCl treatment. PEG treatment was associated with significantly upregulated GuWRKY5/41/53 at 3 h; while cold treatment, with significantly upregulated GuWRKY6/70 only at 24 h. In leaves (Figure 7), NaCl treatment was associated with upregulated GuWRKY19/70/72 at 3 and12 h and GuWRKY5/6/23/42 at 12 and 24 h compared to CK. In addition, the expression of GuWRKY36/41 was only significantly upregulated at 12 h. After PEG treatment, the expression levels of six genes, GuWRKY5/6/41/42/53/72, were significantly upregulated at 12 h, while only GuWRKY70 was significantly upregulated at 6 h. Cold treatment upregulated GuWRKY23/70 at 24 h. Thus, the qRT-PCR-based expression profiles showed that the WRKY family members in licorice responded to various abiotic stresses, providing a theoretical basis for promoting licorice adaptation under adverse conditions.

Figure 5. Relative expression of 10 GuWRKYs in root under low temperature, salt and drought treatments. The y-axis indicates the relative expression levels of GuWRKY genes using the 2^−ΔΔCT method. The x-axis represents the RNA samples from the roots in different treatments at four time points (3 h, 6 h, 12 h, 24 h), from left to right: CK (normal growth), cold (4 °C low temperature treatment), NaCl (150 mmol·L⁻¹ NaCl salt treatment), PEG (15% PEG-6000 drought treatment). Actin was used as the internal reference gene for specific expression analysis. Data shown are mean values ± SEM (n = 3). different capital letters indicate significant differences among the different time points under the same treatment at 0.05 level, and different lowercase letters indicate significant differences among the treatments at the same time point at 0.05 level (LSD test).

Figure 6. Relative expression of 10 GuWRKYs in stem under low temperature, salt and drought treatments. The y-axis indicates the relative expression levels of GuWRKY genes using the 2^−ΔΔCT method. The x-axis represents the RNA samples from the stems in different treatment at four time points (3 h, 6 h, 12 h, 24 h), from left to right: CK (normal growth), cold (4 °C low temperature treatment), NaCl (150 mmol·L⁻¹ NaCl salt treatment), PEG (15% PEG-6000 drought treatment). Actin was used as the internal reference gene for specific expression analysis. Data are shown as mean vales ± SEM (n = 3). different capital letters indicate significant differences among the different time points under the same treatment at 0.05 level, and different lowercase letters indicate significant differences among the treatments at the same time point at 0.05 level (LSD test).

Figure 7. Relative expression of 10 GuWRKYs in leaves under low temperature, salt and drought treatments. The y-axis indicates the relative expression levels of GuWRKY genes using the 2^−ΔΔCT method. The x-axis represents the RNA samples from the leaves in different treatment at four time points (3 h, 6 h, 12 h, 24 h), from left to right: CK (normal growth), cold (4 °C low temperature treatment), NaCl (150 mmol·L⁻¹ NaCl salt treatment), PEG (15% PEG-6000 drought treatment). Actin was used as the internal reference gene for specific expression analysis. Data are shown as mean values ± SEM (n = 3). different capital letters indicate significant differences among the different time points under the same treatment at 0.05 level, and different lowercase letters indicate significant differences among the treatments at the same time point at 0.05 level (LSD test).

Discussion

The WRKY gene family is ubiquitous and is essential for regulating growth, development and stress responses in plants [1, 38]. This family of TFs has received attention in many species, such as rice (O. sativa) [10], cucumber (Cucumis sativus) [39], Scutellaria baicalensis [40], Panax ginseng [41] and Taxus chinensis [42] but not in licorice. Therefore, we conducted a systematic bioinformatic analysis and identified 72 members of the WRKY gene family in the licorice genome. Licorice has fewer WRKY genes than Populus trichocarpa (104) [43] and P. ginseng (137) [44] but more than Akebia trifoliata (63) [45], T. chinensis (61) [42], Isatis indigotica (64) [46], and Salvia miltiorrhiza (61) [47]. Its comparatively high number of members points to the WRKY gene family’s excellent conservation in licorice, which is probably related to gene duplication throughout the formation and evolution of species. The WRKY gene family was increased by two gene replication events (tandem and phased replication), which resulted in gene recombination and amplification, according to Zhang et al.'s findings in S. baicalensis [40]. Furthermore, pertinent research indicated that the primary factor promoting plant WRKY evolution was gene duplication. Also, the results of this study’s subcellular localization revealed that all GuWRKYs were found in the nucleus, indicating that WRKY gene activities are connected to the control of target gene expression in licorice, similarly like they are in Xanthoceras sorbifolium [48] and Juglans regia [49].

In evolutionary trees, TF genes with similar functions tend to cluster together. Phylogenetic analysis showed that group III, accounting for 56% of all GuWRKYs, had significantly more gene members than groups I and II, consistent with Panicum miliaceum [50] and Triticum aestivum [51]. While in Arabidopsis, group II was the largest group, accounting for 62% [9]. These observations indicate both diversification and conservation in the WRKY gene family among the land plants. Specifically, the variation in group III WRKY genes probably accounts for the variation in the WRKY TF family [52].

The diversity in the exon-intron organization is an important aspect of the evolution of gene families [53]. Exon-intron structures can provide valuable information about gene evolution. The number of introns in GuWRKYs ranged from one to eight, which was consistent with the study of WRKY genes in other plants such as peanut (Arachis hypogaea L.) [54] and Akebia (Akebia trifoliata) [45]. These findings showed that the GuWRKY gene family underwent exon gain and loss over its evolution, which most likely increased the functional variety of the GuWRKY gene family members.

Among the ten functional motifs discovered in GuWRKY proteins, motifs 1–3 corresponded to WRKY domains containing the zinc-finger motif in the majority of GuWRKY members. The nuclear localization signal (NLS) was represented by Motif 8, which was unique to Group I. The number and type of conserved motifs were different among the three groups but were roughly the same among the members of a specific group, indicating highly conserved motifs among the GuWRKY proteins. These observations indicated that the WRKY family members had similar structures and biological functions. We also found that the grouping of GuWRKY genes based on phylogeny was consistent with that based on the conserved motifs, indicating that the WRKY genes, to a certain extent, were functionally conserved among plants.

In this study, the analysis of the promoter regions of 72 GuWRKY genes revealed a number of conserved cis-acting regulatory elements involved in a variety of functions, including as abiotic and biotic stress responses (MBS, TC-rich, LTR, ARE, WUN-motif, MBSI, CAT-box, RY-element, circadian and GCN4_motif) and phytohormone regulation (ABRE, AuxRR-core, TGA-element, TCA-element, TGACG-motif, CGTCA-element and GARE-motif). The existence of several cis-acting elements mediating responses to environmental stress and phytohormones indicated the role of licorice WRKY members in different biological processes, consistent with the reports in other plant species. Furthermore, the W-box containing the TGAC sequence, which is the critical recognition site for the WRKY gene, was found in the majority of GuWRKY members. Typically, WRKY TFs regulate target gene expression by binding to W-box on gene promoters [55] or by regulating other WRKY TFs by cross-regulation [56]. For example, overexpression of wheat TaWRKY2 and TaWRKY19 in Arabidopsis demonstrated their role in regulating the expression of target genes by binding to their promoter regions [57]. The present study detected one or more W-boxes in the promoters of 47 GuWRKY genes, confirming this hypothesis.

WRKY genes have been reported to play pivotal roles in regulating transcriptional programming and signal transduction associated with various stresses [1, 4, 58]. Furthermore, the functions of genes in stress and response are tightly connected to the cis-acting elements in gene promoters [59]. Group III AtWRKY53 and AtWRKY70 are known to play essential roles in leaf senescence, whereas AtWRKY70 is vital in osmotic stress signaling and defensive responses [60, 61]. Our results showed that the expression levels of GuWRKY5/6/36, belonging to group III, were upregulated in leaves treated with NaCl and PEG. Various WRKY TFs are known to regulate abiotic stress tolerance through the ABA signaling pathway [62]. In this study, the expression levels of GuWRKY41/70/72 genes were significantly upregulated in the roots, stems and leaves under salt treatment. Consistent with these observations, we detected ABRE, CGTCA-motif and TGA-element in the promoter regions of GuWRKY41/70/72, which indicates their roles in regulating the response to salt stress through ABA, MeJA or auxin signaling pathways. In addition, GuWRKY5/19/53, which were upregulated in response to drought, have the MBS element in their promoter region, which suggests that the encoded proteins control gene expression in response to drought by interacting with particular regulatory elements.

Typically, cold-adapted plants exhibit cold resistance. Yet, it is unclear exactly how WRKYs function in the cold stress response. Our study found that under low temperature stress, especially after 24 h, GuWRKY6/70 were significantly upregulated in the stem. These observations suggested that within 24 h, licorice plants adapt to low temperature stress by expressing genes related to cold stress resistance. In contrast, GuWRKY36/42 had varying degrees of downregulation following cold stress, demonstrating that both genes were negatively regulated. These observations on the variable expression patterns of the GuWRKY genes suggested that the GuWRKY TFs carry out complex regulation to accommodate the physiological changes in licorice under abiotic stress. However, further more in-depth studies about the GuWRKY genes’ regulation mechanism are required.

Conclusions

In this study, we performed a genome-wide analysis and identified 72 GuWRKY genes of licorice. We comprehensively analyzed the WRKY gene family of licorice and classified the 72 GuWRKY genes into three groups. Detailed analysis revealed that GuWRKYs, including GuWRKY5, 6, 19, 23, 36, 41, 42, 53, 70 and 72, regulate the response to abiotic stresses, such as low-temperature, salt and drought stress. Our discoveries on the GuWRKY gene family are useful for future functional investigations to explain their regulation mechanism and generate high-quality licorice cultivars through molecular breeding.

Authors’ contributions

Jing Gao designed the study; Rui Li, Yonggang Yan and Jiakun Yan performed the bioinformatics analysis and qRT-PCR analysis; Nan Wang and Gang Zhang assisted with the analysis and interpretation of the data. Rui Li wrote the manuscript with the contributions of Jing Gao and Nan Wang. All authors read and approved the final manuscript.

Acknowledgements

We thank Shanghai Shenggong Bioengineering Co., Ltd. for their technical support in primer design.

Disclosure statement

The authors declare no conflict of interest.

Data availability statement

The data that support the findings reported in this study are available from the corresponding author (GJ) on reasonable request.

Additional information

Funding

This research was supported by the National Natural Science Foundation of China (NSFC) (31600320), Natural Science Basic Research Program of Shaanxi (2022JM-455), the Key Research and Development Program of Shaanxi Province (2019TSLNY03-1), and Shaanxi Qinchuangyuan "Scientist + Engineer" Team Construction Project (2022KXJ-015).