Abstract Background Leonurus japonicus (L. japonicus) is a herbaceous flowering plant, widely distributed in Asia. Drought is one of the primary environmental stress factors affecting L. japonicus growth. Previous studies have demonstrated that WRKY transcription factors (TFs) play a crucial role in plant responses to drought stress. So far, there has been no research on the function of WRKY genes in L. japonicus. Results The physiological experiment results showed that drought stress significantly increased the malondialdehyde (MDA), proline, and hydrogen peroxide (H₂O₂) content of L. japonicus. Transcriptome analysis revealed significant changes in the expression levels of the WRKY gene family. Based on bioinformatics analysis, 67 WRKY genes (LjWRKYs) were identified in the genome of L. japonicus, with amino acid lengths ranging from 85 to 574. The LjWRKYs can be divided into three subfamilies. Among them, the expression of LjWRKY (1/4/23/44) were significantly up-regulated under drought stress, whereas the expression of LjWRKY (21/25/65) were significantly down-regulated. Additionally, Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis indicated that after drought stress, differentially expressed genes (DEGs) were enriched in plant hormone signal transduction pathway, the MAPK signaling pathway and biosynthesis of secondary metabolites pathway. In the MAPK pathway, there were 19 DEGs, 9 of which contained W-box regions, suggesting that they may be potential regulatory targets of LjWRKY TFs under drought stress. Conclusion These findings suggested that WRKY gene family may participate in the response to drought stress in L. japonicus. This study provides a scientific basis for the further development and functional validation of the WRKY gene family in L. japonicus. Clinical trial number Not applicable. Supplementary Information The online version contains supplementary material available at 10.1186/s12870-025-06606-7. Keywords: Leonurus japonicus, WRKY family, Genome-wide analysis, Drought stress Background Leonurus japonicus (L. japonicus), a member of the Lamiaceae family, is an annual or biennial herbaceous plant widely used in traditional Chinese medicine. L. japonicus is renowned for its effects on women’s physiological health. L. japonicus contains various bioactive components, such as alkaloids, diterpenes, and flavonoids, making it effective for treating multiple diseases [[38]1]. Contemporary research has demonstrated its pharmacological effects, including improvements in hemodynamics, uterine stimulation, anti-inflammatory action, analgesia, and diuretic properties [[39]2]. L. japonicus, is native to China, Korea, India, and other parts of Asia [[40]3]. Nowadays, L. japonicus is widely cultivated in various regions of China, including Xinjiang, Zhejiang, Henan, Guangdong, Shandong, Jiangsu, and Anhui [[41]4]. Drought stress has an adverse impact on plants, triggering changes in gene expression and adjustments in various metabolic pathways. Key genes affected include those involved in stress responses, such as DREB (Dehydration-Responsive Element-Binding proteins) and RD22 (Responsive to Dehydration 22) genes, which play crucial roles in protecting cellular structures and maintaining osmotic balance [[42]5, [43]6]. Metabolic pathways such as the biosynthesis of osmoprotectants and the abscisic acid (ABA) signaling pathway [[44]7, [45]8]. Additionally, genes encoding antioxidant enzymes like SOD (Superoxide Dismutase), CAT (Catalase), and APX (Ascorbate Peroxidase) are up-regulated to mitigate oxidative damage caused by reactive oxygen species (ROS) [[46]9]. In medicinal plants, drought stress has an influence on the synthesis, accumulation, and transformation of active compounds, ultimately affecting the plants’ medicinal value and quality [[47]10]. Reducing soil moisture from 80 to 15% decreases the production of flavonoid compounds in the tubers of Polygonatum kingianum [[48]11]. Similarly, longer irrigation intervals decline antioxidant activity and total phenolic content in Thymus vulgaris [[49]12]. Increasing the interval between irrigation days reduces the total phenolic content, antioxidant activity, and essential oil concentration in Mentha piperita [[50]13]. However, the impact of drought on L. japonicus remains unclear. Several transcription factors (TFs) play a crucial role in plant respond to drought stress, such as AREB/ABF, NAC, WRKY, and MYB TFs [[51]14–[52]18]. The WRKY TFs family is an important group of proteins involved in drought stress response and widely distributed across various plant species [[53]19, [54]20]. WRKY TFs name is derived from the highly conserved WRKY domain, which includes the conserved WRKYGQK sequence and a zinc-finger domain [[55]21–[56]23]. The core sequence of a WRKY motif is WRKYGQK, with variants such as WRKYGKK, WRKYGEK, WRKYGRK, WKKYGQK, WKRYGQK, and WSKYEQK [[57]22, [58]24, [59]25]. WRKY gene families have been identified at the genome level in 234 species, including Arabidopsis thaliana (A. thaliana), Hordeum vulgare, Pinus monticola, Carica papaya, Gossypium aridum, Phoebe bournei and Citrullus lanatus [[60]26–[61]33]. WRKY TFs typically regulate downstream gene expression by specifically binding to W-box elements (TTGACC/T) on promoters [[62]34]. Based on their structural characteristics, WRKY TFs are classified into three major subgroups: Group I contains two WRKY domains and a C2H2-type zinc finger; Group II contains one WRKY domain and a C2H2-type zinc finger; and Group III contains one WRKY domain and a C2HC-type zinc finger [[63]35, [64]36]. Furthermore, Group II is divided into five subgroups (II-a, II-b, II-c, II-d, and II-e) based on their phylogenetic relationships [[65]26, [66]36, [67]37]. The functional importance of WRKY TFs in stress responses has been demonstrated in various plant species. For example, in maize, silencing ZmWRKY79 leads to reduced drought resistance [[68]38]. In Triticum aestivum, knocking out TaWRKY24 results in slow growth of seedlings under salt and drought stress [[69]39]. Similarly, in Gossypium hirsutum, GhWRKY33 acts as a negative regulator mediating the plant’s response to drought stress [[70]40]. These studies highlight the diverse roles of WRKY TFs in enhancing plant tolerance to abiotic stresses, underscoring their potential as targets for improving crop resilience. There are currently no studies on the drought response of L. japonicus, nor on the involvement and of WRKY genes in this process or their identification. Here, bioinformatics analyses and expression pattern analysis of the WRKY gene family under drought stress were performed in L. japonicus. These findings speculated that LjWRKYs may participate in the response to drought stress in L. japonicus. This study provides a scientific basis for the further development and functional validation of the WRKY gene family in L. japonicus. Results Simulation of drought stress in L. japonicus using PEG6000 Drought stress in L. japonicus was simulated using different concentrations of polyethylene glycol (PEG6000). L. japonicus seedlings were treated with 20%, 40%, and 60% PEG6000 solutions for 0, 4, 8, and 12 h, while a control (CK) group received an equivalent amount of deionized water. Following treatment with different PEG6000 concentrations, the L. japonicus leaves displayed varying degrees of wilting (Fig. [71]1A). Under the 20% PEG6000 treatment, the leaves showed almost no noticeable changes throughout the experiment and remained relatively healthy. In the 40% PEG6000 treatment, slight changes began at 4 h, with noticeable curling at 8 h, and wilting starting at 12 h. In contrast, under the 60% PEG6000 treatment, noticeable changes were observed at 4 h, with more severe curling and wilting at 8 and 12 h, indicating a stronger drought stress damage. Fig. 1. [72]Fig. 1 [73]Open in a new tab Physiological responses of L. japonicus to simulated drought stress induced by PEG6000 treatment. (A) Dynamic effects of different concentrations of PEG 6000 treatment on the growth of L. japonicus at 0, 4, 8, and 12 h. (B) MDA content under simulated drought stress. (C) Proline content under simulated drought stress. (D) H[2]O[2] content under simulated drought stress. Asterisks indicate significant differences with respect to values of 0 h (Student’s t-test): **, P < 0.01 Under drought stress induced by varying concentrations of PEG6000, the levels of MDA, proline, and H₂O₂ exhibited significant changes (Fig. [74]1B-D). MDA, a key indicator of membrane lipid peroxidation, exhibited a significant upward trend with increasing PEG6000 concentrations and extended treatment durations. MDA levels increased slightly under 20% PEG6000, while they rose significantly under 40% and 60% PEG6000, particularly at 8 h and 12 h. This indicated that drought stress significantly damages the cell membrane system. Proline, an important osmotic regulator in plants, significantly increased under drought stress conditions. Proline levels increased slightly under 20% PEG6000, while significant increases were observed at 40% and 60% PEG6000, especially at 12 h. This suggested that proline may be involved in osmotic regulation and stress resistance under drought conditions. H₂O₂, a representative molecule of ROS, showed a significant increase with higher PEG6000 concentrations and extended drought stress durations. H₂O₂ accumulation was relatively low under 20% PEG6000, while it increased significantly under 40% and 60% PEG6000, particularly at 8 and 12 h, indicating a strong correlation between ROS accumulation and the intensity of drought stress. Based on phenotypic observations and physiological analyses, 12 h of 40% PEG6000 was shown to produce significant effects. Consequently, 12 h of 40% PEG6000 treatment was chosen as the drought stress condition for this study. Transcriptome analysis of L. japonicus reveals gene expression changes under drought stress To examine gene expression in L. japonicus under drought stress, transcriptome sequencing was performed before and after drought stress treatment. The sequencing results revealed significant alterations in the L. japonicus transcriptome following this treatment, with consistent changes observed across replicate samples (Fig. [75]2A). Principal component analysis (PCA) demonstrated clear separation between the CK and PEG, the plot showed two primary components, PC1 (99.5%) and PC2 (0.4%), the separation between the two groups PC1 reflected significant differences in gene expression between the drought stress and CK conditions. Fig. 2. [76]Fig. 2 [77]Open in a new tab Transcriptome analysis of L. japonicus under simulated drought stress. (A) PCA analysis. (B) Differentially Expressed Genes (DEGs) under simulated drought stress. UP, up-regulated DEGs, DOWN, down-regulated DEGs. (C) Top 20 bubble plot of KEGG enrichment analysis for DEGs. (D) Number of the top 20 transcription factor families identified in the transcriptome analysis of L. japonicus under simulated drought stress A total of 4,140 genes exhibited significant differential expression under drought stress (Fig. [78]2B), with 1,711 genes up-regulated and 2,429 genes down-regulated. This indicated that drought stress substantially impacts gene expression, leading to widespread changes in expression levels. KEGG pathway enrichment analysis showed the DEGs mainly enriched in pathways including the Mitogen-Activated Protein Kinase (MAPK) signaling, plant hormone signal transduction, biosynthesis of secondary metabolites, and diterpenoid biosynthesis, a pathway closely associated with the biosynthesis of isoleosibirin, a terpenoid compound with significant medicinal properties in L. japonicus (Fig. [79]2C). Transcription factor families with the highest expression levels, such as bHLH, MYB, ERF, C2H2, WRKY and NAC, were widely involved in regulating responses to drought and other abiotic stresses (Fig. [80]2D). Additionally, 67 WRKY TFs were identified, and transcriptome analysis confirmed the expression of these genes, highlighting the significant representativeness of WRKY in drought stress responses. Therefore, studying the WRKY family can contribute to a deeper understanding of the molecular mechanisms underlying L. japonicus drought resistance and provide a theoretical basis for breeding drought-tolerant plants. Identification and physicochemical analysis of LjWRKYs To systematically identify the WRKY gene family in L. japonicus, a bioinformatics approach was employed utilizing the Hidden Markov Model (HMM) profile corresponding to the WRKY domain (PF03106). This search against the whole genome of L. japonicus initially yielded 123 candidate genes. Subsequent domain validation using the SMART web tool culminated in the identification of 67 LjWRKYs. Based on their physical locations on the chromosomes (from top to bottom), the LjWRKYs were designated as LjWRKY1 to LjWRKY67 (Fig. [81]3). The 67 LjWRKYs were distributed across all 10 chromosomes of L. japonicus, but their distribution was uneven. Chromosome 7 (Chr 07) contained the most LjWRKYs (n = 11), followed by Chr 01 (n = 10), while Chr 08 had the fewest LjWRKY genes (n = 1). Fig. 3. [82]Fig. 3 [83]Open in a new tab Chromosomal localization of LjWRKYs in L. japonicus The physicochemical characteristics of these proteins, including the number of amino acids, molecular weight (MW), isoelectric point (pI), aliphatic index, instability index, and grand average of hydropathicity, are presented in Table [84]1. The molecular weights of the LjWRKY proteins ranged from 10.2 to 77.4 kDa, with the highest the theoretical isoelectric points (pI) value being 9.72 for LjWRKY41 and the lowest being 4.77 for LjWRKY2, suggesting that each LjWRKY protein possesses a unique functional microenvironment within the plant. The sizes of the encoded proteins ranged from 85 to 574 amino acids. Subcellular localization analysis indicated that most LjWRKY proteins were primarily located in the nucleus, while LjWRKY15 and LjWRKY28 were prediction in the peroxisome, only LjWRKY24 was located in the cytoplasm. Table 1. Genomic information and protein characteristics of 67 LjWRKYs in L. japonicus Gene ID Gene name Chromosome Amino acid number Molecular Weight (Da) pI Instability index Aliphatic index Grand average of hydropathicity Subcellular localization Lj_00070 LjWRKY1 1 252 28516.70 5.66 57.21 55.36 -0.954 Nucleus Lj_00084 LjWRKY2 1 259 29691.56 4.77 62.52 47.37 -1.077 Nucleus Lj_00527 LjWRKY3 1 310 34710.49 9.34 55.60 69.19 -0.697 Nucleus Lj_01078 LjWRKY4 1 525 58123.58 8.39 57.34 57.60 -0.798 Nucleus Lj_01371 LjWRKY5 1 187 21424.14 9.56 36.84 54.76 -0.852 Nucleus Lj_01861 LjWRKY6 1 547 61708.15 6.08 57.62 59.69 -0.773 Nucleus Lj_02230 LjWRKY7 1 364 39873.70 5.84 53.16 62.75 -0.604 Nucleus Lj_02231 LjWRKY8 1 301 34384.16 6.19 52.01 67.71 -0.814 Nucleus Lj_02427 LjWRKY9 1 316 34182.49 5.70 50.11 49.53 -0.778 Nucleus Lj_02523 LjWRKY10 1 157 18690.22 9.58 49.21 50.83 -1.078 Nucleus Lj_03033 LjWRKY16 3 266 30216.71 6.46 39.39 57.56 -0.897 Nucleus Lj_03381 LjWRKY17 3 574 62931.51 6.53 58.24 56.83 -0.840 Nucleus Lj_03759 LjWRKY18 3 286 31704.46 8.76 38.86 68.22 -0.792 Nucleus Lj_04754 LjWRKY19 3 175 20115.39 5.97 42.33 43.31 -0.992 Nucleus Lj_04953 LjWRKY20 3 424 46455.95 6.51 48.61 69.27 -0.646 Nucleus Lj_05080 LjWRKY21 3 354 39457.46 5.25 50.28 54.94 -0.758 Nucleus Lj_05274 LjWRKY22 3 337 36432.16 9.73 35.65 66.59 -0.530 Nucleus Lj_05348 LjWRKY23 3 218 24572.94 8.44 50.63 75.96 -0.590 Nucleus Lj_05577 LjWRKY24 3 85 10281.63 9.17 43.33 46.82 -1.347 Cytoplasm Lj_05852 LjWRKY34 5 271 30751.26 5.10 50.94 68.41 -0.702 Nucleus Lj_06091 LjWRKY35 5 377 42354.14 7.29 54.47 48.12 -1.030 Nucleus Lj_07379 LjWRKY36 5 327 35711.64 9.60 49.45 61.53 -0.676 Nucleus Lj_08060 LjWRKY37 5 282 30958.15 5.66 56.98 57.34 -0.650 Nucleus Lj_08068 LjWRKY38 5 524 58070.67 8.86 58.74 55.25 -0.940 Nucleus Lj_08605 LjWRKY25 4 408 45467.91 8.07 58.27 64.31 -0.774 Nucleus Lj_08856 LjWRKY26 4 314 34672.54 6.08 56.62 54.36 -0.824 Nucleus Lj_09523 LjWRKY27 4 204 23336.78 8.79 46.45 48.24 -1.078 Nucleus Lj_09904 LjWRKY28 4 191 22085.67 8.92 36.69 51.52 -1.024 Peroxisome Lj_09906 LjWRKY29 4 314 35741.61 5.67 63.18 54.01 -0.854 Nucleus Lj_10612 LjWRKY30 4 441 49139.32 5.98 49.40 58.19 -0.821 Nucleus Lj_10675 LjWRKY31 4 505 55113.41 7.73 50.14 56.61 -0.757 Nucleus Lj_11458 LjWRKY32 4 168 19939.36 8.95 46.93 48.63 -1.056 Nucleus Lj_11639 LjWRKY33 4 309 33341.64 9.61 53.95 59.13 -0.557 Nucleus Lj_11964 LjWRKY39 6 307 34338.45 6.97 40.84 68.31 -0.804 Nucleus Lj_12409 LjWRKY40 6 209 24480.82 9.08 62.64 58.23 -0.735 Nucleus Lj_13516 LjWRKY41 6 334 36550.47 9.72 58.16 60.78 -0.645 Nucleus Lj_14182 LjWRKY11 2 372 40958.40 6.01 52.68 66.10 -0.491 Nucleus Lj_14376 LjWRKY12 2 270 30733.71 6.00 64.38 33.96 -1.116 Nucleus Lj_15342 LjWRKY13 2 497 54587.03 9.00 55.80 56.28 -0.787 Nucleus Lj_15627 LjWRKY14 2 289 31890.60 5.45 50.31 65.78 -0.546 Nucleus Lj_16479 LjWRKY15 2 416 46598.15 9.19 49.80 57.16 -0.902 Peroxisome Lj_17528 LjWRKY42 7 347 38331.41 5.78 46.82 54.58 -0.746 Nucleus Lj_17605 LjWRKY43 7 185 20566.38 6.22 58.15 41.84 -1.145 Nucleus Lj_17864 LjWRKY44 7 368 40503.35 6.26 45.74 68.97 -0.561 Nucleus Lj_18115 LjWRKY45 7 282 31751.08 6.89 61.38 50.18 -0.879 Nucleus Lj_18194 LjWRKY46 7 411 45057.10 5.54 58.57 58.49 -0.727 Nucleus Lj_18280 LjWRKY47 7 581 63291.83 6.34 45.92 63.10 -0.694 Nucleus Lj_18384 LjWRKY48 7 714 77422.88 6.05 53.78 49.08 -0.864 Nucleus Lj_18624 LjWRKY49 7 538 58156.66 6.73 47.08 59.39 -0.648 Nucleus Lj_18829 LjWRKY50 7 92 10252.34 7.97 67.85 37.07 -1.280 Nucleus Lj_18934 LjWRKY51 7 255 28546.39 6.66 51.60 49.02 -0.981 Nucleus Lj_19002 LjWRKY52 7 291 32988.43 4.87 85.53 53.99 -0.927 Nucleus Lj_19526 LjWRKY53 8 323 36272.02 5.80 70.10 48.02 -0.991 Nucleus Lj_21218 LjWRKY54 9 535 59624.24 6.88 47.45 61.53 -0.990 Nucleus Lj_21857 LjWRKY55 9 482 52974.56 7.56 51.48 55.19 -1.011 Nucleus Lj_22027 LjWRKY56 9 474 51872.89 6.54 51.87 63.65 -0.689 Nucleus Lj_22312 LjWRKY57 9 404 45985.51 8.78 57.85 41.51 -1.194 Nucleus Lj_22909 LjWRKY58 9 496 55027.86 8.45 54.23 63.12 -0.622 Nucleus Lj_22996 LjWRKY59 9 453 50277.66 7.00 61.64 61.99 -0.814 Nucleus Lj_23020 LjWRKY60 9 261 29093.27 5.94 61.21 59.77 -0.815 Nucleus Lj_23411 LjWRKY61 10 269 30992.91 8.50 76.65 47.88 -1.004 Nucleus Lj_24308 LjWRKY62 10 202 22467.07 6.96 36.53 68.12 -0.621 Nucleus Lj_24569 LjWRKY63 10 344 38401.44 9.52 56.93 62.67 -0.795 Nucleus Lj_25000 LjWRKY64 10 208 24137.86 8.97 45.91 46.78 -1.100 Nucleus Lj_25289 LjWRKY65 10 284 31702.03 6.24 58.27 47.11 -0.844 Nucleus Lj_25293 LjWRKY66 10 166 19566.97 9.19 65.86 53.43 -0.922 Nucleus Lj_25340 LjWRKY67 10 305 34393.58 5.26 61.77 67.11 -0.635 Nucleus [85]Open in a new tab Analysis of intron-exon structure and motifs in LjWRKYs of L. japonicus To elucidate the key role of exon/intron patterns in the evolution of the L. japonicus gene family, the structure of LjWRKYs was determined through exon/intron boundary analysis. The 67 LjWRKYs contained 2 to 7 exons (7 with 2 exons, 34 with 3 exons, 10 with 4 exons, 11 with 5 exons, 3 with 6 exons, and 2 with 7 exons) (Fig. [86]4). This suggested that exon loss and gain events occurred during the evolution of the LjWRKY family, potentially contributing to the functional diversity among LjWRKYs. Fig. 4. [87]Fig. 4 [88]Open in a new tab Intron-Exon structure and motifs in LjWRKYs of L. japonicus. Schematic representation of the WRKY genes, including 3’ and 5’ UTR (untranslated regions) and CDS (coding sequences) Conserved motifs in the LjWRKY proteins were predicted using the MEME program. A total of 10 conserved motifs were identified among the 67 LjWRKY protein sequences, with each motif comprising 15 to 50 amino acids. Each WRKY protein contained 2 to 6 motifs. Notably, three motifs, designated as motif_1, motif_3, and motif_8, formed the core structure of the WRKY domain. Motif_1 was present in 66 LjWRKY proteins, while LjWRKY50 contained only motif_8. According to the phylogenetic tree classification, Group I had 3 to 5 motifs, with most members containing 5 motifs. In Group II, the number of motifs ranged from 1 to 6, with significant differences observed between the subfamilies. The number of motifs in subfamilies Group II-a and Group II-b ranged from 1 to 6, while in subfamily Group II-c, only LjWRKY43 contained two motifs, and the others had three motifs. Both subfamilies Group II-d and Group II-e contained three motifs. Overall, almost all LjWRKY proteins within the same subfamily exhibited very similar motif compositions, suggesting that they have similar functions. Phylogenetic analyze and multiple sequence alignment of LjWRKYs To analyze the evolutionary relationships of LjWRKY proteins, a phylogenetic tree (Fig. [89]5) was constructed using the maximum likelihood method based on 67 LjWRKY proteins from L. japonicus and 21 AtWRKY proteins from A. thaliana. Based on the characteristics of AtWRKYs in A. thaliana, the 67 LjWRKY proteins were classified into three groups: I, II, and III. 13 LjWRKYs were assigned to Group I, 48 to Group II, and the remaining 6 to Group III. Group II, which contains 48 LjWRKYs, was further divided into five subgroups: II-a (n = 5), II-b (n = 9), II-c (n = 19), II-d (n = 7), and II-e (n = 8). Fig. 5. [90]Fig. 5 [91]Open in a new tab Phylogenetic tree of WRKY proteins from L. japonicus and A. thaliana Multiple sequence alignment of the LjWRKY proteins revealed that members of Group I contained two complete WRKY domains, including the conserved WRKYGQK sequence, and possess C2H2-type zinc fingers. All 48 Group II proteins contain one WRKY domain, with 45 exhibiting the conserved WRKYGQK sequence. However, some LjWRKYs, such as LjWRKY32, LjWRKY19, and LjWRKY24, display variations in this sequence, like WRKYGKK. Group II LjWRKYs typically exhibit the C2H2-type zinc finger, although some, such as LjWRKY43, LjWRKY60, and LjWRKY50, showed zinc finger structure loss. All Group III WRKYs contain the WRKYGQK sequence and a C2HC-type zinc finger (Fig. [92]6). Fig. 6. [93]Fig. 6 [94]Open in a new tab Multiple sequence alignment of the WRKY domains of 67 LjWRKY proteins. The red box represents the WRKY domain, and the blue box represents the zinc finger structure LjWRKYs duplication and synteny analyses To gain insight into the expansion of the LjWRKYs, an investigation of gene duplication events was conducted. Five pairs of tandem repeats were identified: LjWRKY7 and LjWRKY8, LjWRKY28 and LjWRKY29, LjWRKY37 and LjWRKY38, and LjWRKY65 and LjWRKY66. In addition, gene duplication occurred across different chromosomes, such as between LjWRKY41 and LjWRKY36 or between LjWRKY53 and LjWRKY61 (Fig. [95]7A). These results suggested that fragment duplication events have played more important role in the expansion of LjWRKYs compared to tandem duplication events. Fig. 7. [96]Fig. 7 [97]Open in a new tab Syntenic analysis of LjWRKYs in L. japonicus. (A) Circos plot showing the gene duplication events of LjWRKYs within the L. japonicus genome. The red lines indicate duplicated gene pairs. (B) Syntenic relationships of WRKY genes among L. japonicus, A. thaliana, and S. bowleyana. Blue lines represent syntenic WRKY gene pairs between the species To further understand the evolutionary relationships among WRKY members in different plant species, we analyzed the syntenic relationships between LjWRKYs and their homologs in other species. We performed synteny analysis of L. japonicus with the model plant A. thaliana and the Lamiaceae species Salvia miltiorrhiza (S. miltiorrhiza) (Fig. [98]7B). Among the 67 identified LjWRKYs, 46 and 41 were collinear with WRKYs in S. miltiorrhiza and A. thaliana, respectively. A total of 46 homologous gene pairs were identified between L. japonicus and S. miltiorrhiza, while 41 homologous gene pairs were identified between L. japonicus and A. thaliana. This high gene homology between L. japonicus and S. miltiorrhiza may be attributed to the conserved nature of WRKY genes during evolution. Analysis of cis-regulatory elements in LjWRKYs promoters To further investigate the functions and expression patterns of the LjWRKY genes in L. japonicus, this study conducted a comprehensive analysis of cis-acting elements in the promoter regions of 67 LjWRKY genes using the PlantCARE database(Table [99]S1). The study illustrated the distribution and abundance of 14 types of cis-acting elements within the promoter regions of 67 LjWRKY genes (Fig. [100]8, Fig [101]S1). Each type of cis-acting element is linked to specific functions, such as stress responses, plant hormone responses, and growth and developmental regulation. Fig. 8. [102]Fig. 8 [103]Open in a new tab Cis-acting elements in the promoter regions of LjWRKYs in L. japonicus Among the stress response-related cis-acting elements, the W-box was the most prominent, which was identified in 40 LjWRKYs and known to be strongly associated with WRKY TFs binding. A total of 48 LjWRKYs were found to contain the MBS cis-acting element, which is associated with drought responsiveness. The LTR (low-temperature response element) was identified in 27 LjWRKYs, with relatively higher frequencies in LjWRKY5 and LjWRKY23. TC-rich repeats were found in 36 LjWRKYs, with relatively higher frequencies in LjWRKY12 and LjWRKY44. Among the plant hormone response-related cis-acting elements, ABRE (abscisic acid-responsive element) was identified in 52 LjWRKYs, with notable enrichment in LjWRKY23, LjWRKY35, and LjWRKY53. The TGACG-motif (jasmonic acid-responsive element) was identified in 32 LjWRKYs, with relatively higher frequencies in LjWRKY3 and LjWRKY19. The CGTCA-motif (methyl jasmonic acid-responsive element) was present in 39 LjWRKYs and showed the highest enrichment in LjWRKY44. The TGA-element (auxin-responsive element) was found in 18 LjWRKYs, including LjWRKY33 and LjWRKY56. The AuxRR-core (cytokinin-responsive element) was detected in 12 LjWRKYs, with a sparse distribution in genes such as LjWRKY16 and LjWRKY25.Among the growth and development-related cis-acting elements, the CAT-box was present in 19 LjWRKYs, with relatively higher occurrences in LjWRKY5 and LjWRKY19. The O2-site (a regulatory element of endosperm expression) was found in 5 LjWRKYs, mainly in LjWRKY23. The P-box was identified in 10 LjWRKYs, with a scattered distribution, particularly in LjWRKY44 and LjWRKY61. Overall, the figure clearly reveals the potential functional differences of each LjWRKYs in stress responses, hormone signaling, and growth and development, providing a solid data foundation for further research on the functions of the LjWRKYs. Expression patterns of LjWRKYs in L. japonicus across tissues and under drought stress To better understand the functional differentiation of the LjWRKYs, we conducted a comprehensive analysis of the tissue expression patterns and stress-responsive expression levels. The heatmap illustrated the expression patterns of LjWRKYs across different organs of L. japonicus (leaf, stem, flower, and root) (Fig. [104]9A, Table [105]S2). The heatmap clearly showed that most LjWRKYs exhibit higher expression levels in leaf and root, particularly genes such as LjWRKY4, LjWRKY5, LjWRKY22, and LjWRKY23. These genes may play important roles in plant growth and environmental responses. In contrast, certain genes (LjWRKY44, LjWRKY64) showed negligible or very low expression in the stem and flower, which may suggest that these genes have weak or unimportant functions in these organs. Notably, some genes, such as LjWRKY21 and LjWRKY67, have distinct tissue-specific expression patterns, indicating that they may be involved in the development and function of specific organs. Fig. 9. [106]Fig. 9 [107]Open in a new tab Expression patterns of the LjWRKYs in different tissues and under drought stress. (A) The expression patterns of LjWRKYs in different tissues, including leaf, stem, flower, and root. The color gradient ranges from red (high expression) to white (low expression), indicating variations in expression levels. (B) The expression changes of LjWRKYs under drought stress treatments, represented by three biological replicates (PEG-1, PEG-2, PEG-3), are shown as log₂ fold change (log₂FC) values relative to the control group (CK). The color gradient ranges from purple (high expression) to cyan (low expression), illustrating the magnitude of expression changes. Gene grouping is based on the phylogenetic tree classification Groups I to III In contrast, many genes showed a significant increase in expression levels under drought treatment (Fig. [108]9B), particularly LjWRKY1, LjWRKY23, and LjWRKY57, which were up-regulated under drought stress conditions, suggesting that these genes may play key roles in the plant’s response to drought. Additionally, some genes, such as LjWRKY34 and LjWRKY35, exhibit relatively low expression under drought treatment, suggesting that these genes may not be involved in drought responses or that their roles in drought stress are not significant. To validate the accuracy of the transcriptome results, we selected a representative gene from each group and performed quantitative reverse transcription PCR (qRT-PCR). LjWRKY1, LjWRKY4, LjWRKY23, LjWRKY44 expression was up-regulated after drought treatment and LjWRKY21, LjWRKY25, LjWRKY65 expression was down-regulated after drought treatment. The qRT-PCR validation results were closely consistent with the transcriptome data, further confirming the differential expression patterns of these genes under drought stress (Fig [109]S2). Overall, this figure reveals the expression characteristics of LjWRKYs in different organs and their changes under drought stress. The comprehensive analysis of the organ expression patterns and stress-responsive expression levels enhances our understanding of the functional differentiation of the LjWRKYs and its roles in plant growth, development, and environmental adaptation. Regulatory role of LjWRKYs in MAPK pathway under drought stress During the initial analysis of DEGs, the MAPK pathway was found to be significantly enriched (Fig. [110]2C). The response mechanisms of the MAPK signaling pathway under drought stress were subsequently investigated. Furthermore, a gene expression heatmap was employed to visually demonstrate variations in gene expression under drought stress (Fig. [111]10). From the Fig. [112]10, it can be observed that the expression levels of PYL/PYR (Lj_22318, Lj_08822, Lj_05087, Lj_04481, Lj_22066) genes varied under drought stress. For example, Lj_22318 exhibited fold change of 1.04, 0.68, and 1.36, indicating a significant up-regulation under drought stress. PP2C (Lj_20896, Lj_05704, Lj_20999, Lj_02051) genes showed a substantial increase in expression levels after treatment, especially Lj_02051, which reached as high as 3.32 under drought stress, indicating its sensitivity to drought signals. SnRK2 (Lj_08942, Lj_08623) showed a slight up-regulation. In the MAPKKK17/18 pathway, Lj_22021 expression significantly increased, reaching 3.18 and 3.40 after drought stress treatments, respectively. Lj_04014, associated with MKK3, showed minimal changes after drought stress. Among genes related to CAT1 metabolism, Lj_22257 expression significantly increased after drought treatment, with fold change of 0.65, 0.75, and 0.75 under drought stress, while Lj_13701 showed lower variations. Fig. 10. [113]Fig. 10 [114]Open in a new tab Expression heatmap of MAPK pathway genes in L. japonicus under drought stress. The expression changes of LjWRKYs under drought stress treatments, represented by three biological replicates (PEG-1, PEG-2, PEG-3), are shown as log₂ fold change (log₂FC) values relative to the CK. Purple indicates up-regulation, and blue indicates down-regulation. The genes marked with a five-pointed star contain W-box elements In summary, these gene expression changes reveal the activity levels of different signaling pathways under drought stress. Notably, the promoters of genes containing W-box structures (Fig. [115]S3), such as Lj_22257 (CAT1), Lj_20896 (PP2C), Lj_02051 (PP2C), and Lj_22021 (MAPKKK), exhibited significant up-regulation, suggesting they may be the target genes regulated by LjWRKYs. Through the regulation of these genes, the MAPK signaling pathway effectively responds to drought stress, enabling plants to adapt to adverse conditions. Discussion Drought stress has a significant impact on plant growth and development, particularly on medicinal plants [[116]41, [117]42]. It affects various physiological processes and medicinal properties of plants [[118]43]. Drought stress triggers a series of physiological responses, including alterations in metabolic pathways and the accumulation of stress-related compounds. These compounds, such as ROS, secondary metabolites, and osmoprotectants, play a critical role in the plant’s response to drought stress [[119]44, [120]45]. However, they may also influence the medicinal characteristics of the plants. Transcriptome sequencing was performed to investigate the molecular response of L. japonicus to drought stress, identifying 4,140 DEGs (Fig. [121]2). KEGG pathway analysis emphasized the importance of the MAPK and plant hormone signaling pathways, consistent with their roles in other drought-tolerant plants, such as A. thaliana and Oryza sativa [[122]46, [123]47]. The MAPK signaling pathway plays a critical role in regulating stress responses, as observed in other species [[124]48]. Among the differentially expressed transcription factors, WRKY family members stand out due to their high expression and involvement in multiple stress pathways. WRKY TFs integrate various stress signals and regulate critical networks involved in drought tolerance, as demonstrated in A. thaliana and O. sativa [[125]49, [126]50]. The genome size of L. japonicus is 472 Mb, and 67 WRKY family members were identified within its genome (Fig. [127]3, Table [128]S1). Our findings suggest that plant species closely related to the Lamiaceae have similar WRKY gene numbers, while more distantly related species exhibit greater variation. For example, 61 WRKY genes were identified in Salvia miltiorrhiza [[129]51]. Similarly, 58 WRKY genes were found in Andrographis paniculata, a species from the Acanthaceae family, which is phylogenetically close to the Lamiaceae [[130]52]. In contrast, Miscanthus sinensis contains 179 WRKY genes [[131]60], the Brassicaceae plant Brassica napus harbors 287 WRKY genes, and the legume Glycine max has 182 WRKY genes [[132]53, [133]54]. Interestingly, the number of WRKY genes does not appear to directly correlate with genome size. For instance, Vitis vinifera has a genome size of 487 Mb but contains only 55 WRKY genes [[134]55], while the herbaceous plant Brachypodium distachyon (with a genome size of 272 Mb) contains 86 WRKY genes [[135]56]. These comparisons indicate that, in addition to genome size, other factors may influence the expansion or contraction of the WRKY gene family across different plant species. Evolutionary analysis of 67 LjWRKY genes classified them into three main groups: Group I, Group II (further subdivided into five subgroups: II-a, II-b, II-c, II-d, and II-e), and Group III (Fig. [136]5). Of these, 13 genes contain dual WRKY domains with a C2H2-type zinc finger motif, while the remaining genes, which have a single WRKY domain, are classified into Group II (48 members) or Group III (6 members). This classification aligns with previous studies on Salvia miltiorrhiza and Andrographis paniculata, suggesting structural conservation across species [[137]51, [138]52]. Group II, which contains the most members, may play a complex role in plant adaptation and growth, potentially due to frequent gene duplications. WRKY proteins typically contain a conserved WRKYGQK heptapeptide motif, however, in some Group II members, such as LjWRKY32, LjWRKY19, and LjWRKY24, this motif was replaced by WRKYGKK. This variant has been observed in species such as O. sativa and Ricinus communis [[139]47, [140]55–[141]57]. Additionally, some LjWRKY genes exhibit partial or complete loss of zinc finger structures, a rare phenomenon that could affect their DNA-binding and regulatory functions. These structural variations may represent adaptive mechanisms in L. japonicus, however, further research is needed to explore their biological significance. WRKY TFs with the WRKYGQK motif, such as AtWRKY, typically bind to W-box elements in the promoters of target genes. In contrast, WRKY TFs with the WRKYGKK variant, such as GmWRKY6 and GmWRKY21, cannot bind to the W-box [[142]58]. This suggests that LjWRKY TFs lacking the typical WRKYGQK motif may not recognize the W-box, potentially affecting their DNA-binding activity and regulatory functions. The MAPK signaling pathway plays a pivotal role in plant stress responses, particularly in osmotic regulation and oxidative defense mechanisms [[143]59]. This study highlights the complex regulatory processes underlying the plant’s response to drought via the MAPK pathway (Fig. [144]10). In L. japonicus, key components of the ABA signaling pathway, such as PYL/PYR proteins, are significantly up-regulated under drought stress, highlighting their involvement in the regulation of drought responses [[145]45]. Similarly, in Sorghum bicolor, MAPK cascade genes, particularly SbMPK14, was identified as important regulators of drought sensitivity. The up-regulation of these genes suggests a similar role for MAPK components in regulating drought tolerance in L. japonicus and Sorghum bicolor [[146]60]. In Vitis vinifera, a transcriptomic analysis revealed a complex immune signaling network in response to drought stress, where MAPK signaling interacts with immune pathways to enhance drought tolerance. This suggests that MAPK signaling not only acts as a stress response mediator but also plays a critical role in modulating the immune response, providing additional protection against drought-induced damage [[147]61]. This interaction between stress and immune signaling is also evident in L. japonicus, where MAPK signaling appears to regulate both stress responses and immune function, contributing to enhanced drought resilience. Moreover, a comprehensive transcriptomic analysis of Gossypium anomalum seedlings under drought stress identified key regulators and metabolic pathways involved in the plant’s response to drought. Among these, the MAPK signaling pathway was found to be central in regulating the expression of genes involved in drought adaptation. This aligns with findings in L. japonicus, where MAPK signaling also plays a critical role in regulating metabolic pathways that help mitigate drought stress [[148]62]. Additionally, genes involved in CAT1 metabolism are crucial for alleviating oxidative stress during drought, with sustained up-regulation of Lj_22257 and stable expression of Lj_13701 indicating differential regulation of antioxidant defenses. Several up-regulated genes, including Lj_22257, Lj_20896, Lj_02051, and Lj_22021, contain W-box motifs, suggesting that LjWRKY TFs may regulate these drought-responsive genes [[149]63]. Consistent with findings in other species, L. japonicus WRKY TFs, likely through the ABA-mediated MAPK signaling pathway, may enhance drought tolerance. This study speculates that LjWRKY TFs directly regulate PP2C gene expression and indirectly participate in the MAPK cascade, thereby contributing to enhanced drought adaptation. The activation of MAPKKK17 further propagates the drought signal, influencing downstream genes and enhancing plant tolerance to osmotic stress and oxidative damage. Conclusion This study investigated the drought stress response of L. japonicus at both physiological and molecular levels. In this study, we identified 67 LjWRKY family members from the L. japonicus genome. By analyzing their phylogenetic relationship, genes duplication and cis-acting elements in the promoter regions using bioinformatics methods and using transcriptome to analyze the expression patterns of the LjWRKYs in different tissues and under drought stress. Our findings revealed that WRKY family may work as key TFs in the drought response of L. japonicus through the ABA-MAPK signaling cascade. This study provides valuable insights into the molecular mechanisms underlying drought tolerance in L. japonicus, offering potential targets for enhancing drought stress in medicinal plants. Materials and methods Plant material and growth conditions The full name of the material used in this study is Leonurus japonicus Houtt., collected from Yiyuan, Zibo, Shandong Province, China (118.07°E, 36.08°N), identified by Jianing Xu. The sampling was conducted on public land, a sampling permit was obtained from the consent of the users of that land prior to the collection of experimental material. Seeds of L. japonicus were sown in plastic pots containing a mixture of perlite and silica sand (2:1 v/v). The pots were placed under long-day conditions (25 °C, 20,000 lx, 14 h light/10 h dark) for germination. For subsequent experiments, plants that were 11 weeks old and of uniform height were used. The plant samples were kept in the laboratory of the School of Life Sciences and Medicine, Shandong University of Technology. The voucher specimen of this material has been deposited in a publicly available herbarium in School of Life Sciences and Medicine, Shandong University of Technology. Drought stress simulation To investigate the expression changes of LjWRKYs under drought stress, different concentrations of polyethylene glycol (PEG6000) were used. Healthy L. japonicus seedlings were treated with PEG6000 solutions at concentrations of 20%, 40%, and 60% to simulate varying levels of drought stress. The control group was treated with an equal volume of distilled water. Observations of the control and treated plants were conducted at 0 h, 4 h, 8 h, and 12 h. Plant samples were collected and immediately frozen in liquid nitrogen. Total RNA was extracted using a commercial RNA extraction kit. The quality and purity of RNA were evaluated based on the OD[260/280] and OD[260/230] ratios. Three biological replicates per group. MDA content, proline and H[2]O[2] activity measurement The contents of malondialdehyde (MDA), proline, and hydrogen peroxide (H₂O₂) in fresh L. japonicus leaves were measured using assay kits (MDA: G0109W, Proline: G0111W, H₂O₂: G0112W) provided by Suzhou Grace Biotechnology Co., Ltd., following the manufacturer’s instructions. Fresh leaf samples (0.5 g) were washed, dried, and ground into a fine powder with liquid nitrogen, then homogenized with the respective extraction buffers provided in the kits. The extracted solutions were mixed with detection reagents, and absorbance values were recorded using a Thermo Fisher Scientific microplate reader (model Variskan). at specific wavelengths: MDA: 532 nm, Proline: 520 nm, H₂O₂: 415 nm. The final contents of MDA, Proline, and H₂O₂ were calculated based on the absorbance readings. Three biological replicates per group. Transcriptome sequencing and qRT-PCR validation RNA was then subjected to library preparation for transcriptome sequencing. mRNA was enriched using oligo(dT) beads, fragmented, and converted into cDNA using reverse transcriptase. The cDNA was then subjected to library construction, including adapter ligation, PCR amplification, and purification. Libraries were sequenced on a high-throughput sequencer, such as Illumina Nova-Seq or Hi-Seq, to generate raw sequence data. The obtained reads were quality checked and trimmed. The primers were done using Primer-BLAST (Primer designing tool (nih.gov); (table [150]S3) their synthesis was entrusted to GENEWIZ (Suzhou, China), and rRNA served as a reference gene. Total RNA was extracted from the samples using a Plant RNA Rapid Extract Kit (coolaber, Beijing, China), and cDNA was obtained using a TransScript^® One-Step gDNA Removal and cDNA Synthesis SuperMix (Trans Gene). The qRT-PCR was performed on a StepOne Real-Time PCR System (Applied Biosystems) using SYBR Green realtime PCR master mix (SYBR Green, Toyobo). The PCR reaction program consisted of 95 °C for 30 s, followed by 40 cycles of 95 °C for 10 s, 60 °C for 30 s, 95 °C for 15 s, and 60 °C for 1 min, finishing with 95 °C for 15 s. The experiment was repeated three times, the relative expression level was calculated by the 2^−∆∆CT method. Three biological replicates per group. KEGG annotation analysis of LjWRKYs KEGG enrichment analyses were conducted for LjWRKYs to further explore their biological functions, using a pipeline on the Omicshare Bioinformatics Cloud Platform (GENEDENOVO, Guangzhou, China). Plant transcription factor statistics Utilize the transcription factor databases ([151]https://planttfdb.gao-lab.org/) to annotate genes for transcription factors (TFs), determine whether a gene is a member of a specific transcription factor family, and summarize how many genes belong to each TF family. Sequence retrieval and LjWRKY identification First, download the seed file about WRKY gene (PF03106) in Pfam official website ([152]http://pfam.xfam.org/). The genome of L. japonicus was compared with TBtools, and the sequences with E-value > 0.001 were removed. The matching sequence id and the sequence containing WRKY domain were further verified by making Wayne diagram on the online website ([153]https://bioinformatics.psb.ugent.be/webtools/Venn/). Then 123 sequences were extracted by [154]https://www.tbtools.com/. Then, whether the 123 sequences contain WRKY domain was further verified by SMART online website ([155]https://smart.embl.de/). Upload the sorted sequence files on the online website MEME ([156]https://meme-suite.org/meme/tools/meme) to find and identify motif. The parameter of the maximum number of motifs is set to 10, the minimum motif width is 6, and the maximum motif width is 50. By selecting the longest sequence and removing redundant sequences at the same time, 67 WRKY protein sequences remain. Additionally, the Compute pI/Mw tool on the ExPASy server ([157]http://web.expasy.org/compute_pi/) was used to predict the theoretical isoelectric point (pI) and molecular weight (Mw) of the WRKY proteins from L. japonicus. Chromosomal location and gene duplication Based on the physical location information of the L. japonicus genome, use the Gene Location Visualize from GTF/GFF function in TBtools ([158]https://github.com/CJ-Chen/TBtools) to map all LjWRKYs onto the L. japonicus chromosomes. Perform gene duplication event analysis using MCScanX with default parameters. We used the one step MCScanX in TBtools to illustrate the collinearity of WRKY genes between L. japonicus and two selected species (S. bowleyana, A. thaliana). Multiple sequence alignment phylogenetic analysis and classification of LjWRKYs Using MEGA 11.0, perform multiple sequence alignment of WRKY transcription factors from A. thaliana and L. japonicus. Subsequently, construct a phylogenetic tree based on the amino acid sequences of the two species using the maximum likelihood (ML) method with 1000 bootstrap replicates, the Jones-Taylor-Thornton (JTT) model, and Gamma distribution (G). Then, use Chiplot ([159]https://www.chiplot.online/tvbot.html) to transform the chart. Based on the phylogenetic tree and previously reported AtWRKYs classifications [[160]64], categorize LjWRKYs into different groups. Identification of gene structure and conserved motifs Analyze and visualize the exon/intron composition of LjWRKYs genes using TBtools ([161]https://github.com/CJ-Chen/TBtools). Discover motifs in LjWRKYs using the MEME 5.5.5 online tool ([162]https://meme-suite.org/meme/tools/meme) with the parameters set to: site distribution as ZOOPS, and the number of motifs set to 10. Then, use Chiplot ([163]https://www.chiplot.online/tvbot.html) to transform the chart. Promoter cis-acting element and annotation analysis of LjWRKYs The cis-acting elements of the 2000 bp promoter sequence of LjWRKYs were analyzed using the PlantCARE database ([164]http://bioinformatics.psb.ugent.be/webtools/plantcare/html/). Creating a Heatmap Using Excel. Gene expression analysis of LjWRKYs Total RNA was isolated from four distinct tissues using the Trizol reagent kit from Invitrogen (Carlsbad, CA, USA), following the manufacturer’s instructions. Sequencing libraries for mRNA were prepared and sequenced on an Illumina HiSeqTM 4000 platform by Gene Denovo Biotechnology Co. (Guangzhou, China). A reference genome index was created, and clean paired-end reads were aligned to this reference genome using HISAT2 version 2.1.0 with the “rna-strandness RF” option and default settings for other parameters. The aligned reads for each sample were then assembled using StringTie v1.3.1 in a reference-based manner. Expression levels within each transcriptional region were estimated using FPKM (fragments per kilobase of transcript per million mapped reads) values, which were computed with RSEM software to assess expression abundance and variation. Each sample was analyzed in triplicate to ensure data reliability. To explore the specific expression patterns of LjWRKYs in different L. japonicus organs, RNA-seq data from L. japonicus mature roots, leaves, stems and flowers were were from our lab and used for gene expression analysis. The heatmap was drawn in TBtools using heatmap. Statistical analysis All data are presented as means with standard deviations (SDs). Statistical analysis and bar graphs were created by GraphPad Prism 8.3.0 software. * p < 0.05, ** p < 0.01. Electronic supplementary material Below is the link to the electronic supplementary material. [165]12870_2025_6606_MOESM1_ESM.jpg^ (2.3MB, jpg) Supplementary Material 1: Fig. S1: Distribution of cis-acting elements in the LjWRKYs [166]12870_2025_6606_MOESM2_ESM.png^ (480.1KB, png) Supplementary Material 2: Fig. S2: Expression levels of LjWRKYs in response to drought stress as determined by qRT-PCR and RNA-seq. [167]12870_2025_6606_MOESM3_ESM.png^ (1.1MB, png) Supplementary Material 3: Fig. S3: Distribution of W-box motifs in MAPK Pathway-Associated genes. [168]12870_2025_6606_MOESM4_ESM.xlsx^ (12.4KB, xlsx) Supplementary Material 4: Table S1: Cis-Acting elements in the promoter regions of LjWRKYs in L. japonicus [169]12870_2025_6606_MOESM5_ESM.xlsx^ (9.9KB, xlsx) Supplementary Material 5: Table. S2: Expression patterns of the LjWRKYs in different tissues and under drought stress. [170]12870_2025_6606_MOESM6_ESM.xlsx^ (43.5KB, xlsx) Supplementary Material 6: Table. S3: Primers Used for qRT-PCR Analysis. [171]Supplementary Material 7^ (12.8KB, docx) Acknowledgements