Computational Survey Identifies 35 WRKY Transcription Factor Genes in Chia

Chia, the ancient Mesoamerican grain that has become one of the world's most celebrated functional foods, has received a detailed computational inventory of one of its most important gene families. In a study published in Discover Plants, researchers report the first genome-wide identification and characterization of the WRKY transcription factor family in Salvia hispanica, identifying 35 genes, designated CsWRKY1 through CsWRKY35, and classifying them into the seven canonical phylogenetic groups that define this family across the plant kingdom. Even so, the study delivers something chia researchers have lacked until now: a systematic framework for locating, naming and prioritizing the genes most likely to govern how this crop responds to its environment and orchestrates its own development.

The significance of the WRKY family in plants is difficult to overstate. WRKY transcription factors are among the largest families of gene regulators in the plant kingdom, acting as molecular switches that reprogram gene expression during defense responses, abiotic stresses such as drought and salinity, senescence and normal developmental processes. They are defined by one or two conserved WRKY domains, each carrying the invariant seven-amino-acid sequence WRKYGQK together with a zinc-finger motif of either the C2H2 type, found in groups I and II, or the C2HC type, characteristic of group III.

These domains bind to DNA sequences called W-boxes, with the core motif TTGACC/T, in the promoters of target genes. Genome-wide WRKY inventories have previously been compiled for Arabidopsis thaliana with 74 genes, rice with 109, maize with 86, tomato with 67, grapevine with 55, cotton with 97 and sugar beet with 72, but chia had never been surveyed in this way. Of the 35 genes recovered, 12 belong to group I, which carries two WRKY domains and C2H2 zinc fingers; 15 belong to group II, with a single domain and C2H2 motif, subdivided into classes IIa through IIe; and 8 belong to group III, distinguished by the C2HC zinc finger.

Protein lengths ranged from 143 amino acids for CsWRKY21 to 606 amino acids for CsWRKY20, with a mean of 332.8 amino acids. The mean molecular weight was 37,499 daltons and the mean theoretical isoelectric point was 7.62, spanning a range from 5.16 to 9.67. All proteins were predicted to be hydrophilic, with uniformly negative GRAVY values and a mean aliphatic index of 76.25.

Subcellular localization predictions, generated with three separate algorithms, WoLF PSORT, DeepLoc and CELLO, placed 18 of the 35 proteins, or 51.4 percent, primarily in the nucleus, which is consistent with their role as transcriptional regulators that must physically contact DNA. Another 14 proteins, 40.0 percent, were predicted to reside in the cytoplasm, and 3 were assigned to the chloroplast. Gene-structure analysis revealed between 2 and 7 exons per gene, with the two-domain group I members generally carrying the more complex architectures of 4 to 7 exons, while group IIa and IIb members predominantly showed a streamlined two-exon, one-intron organization.

Motif analysis with the MEME suite identified 15 conserved motifs across the family. Motif 1 corresponds to the core WRKY domain containing the WRKYGQK heptapeptide, which was present in all 35 proteins without exception, while Motif 2 matched the C2H2 zinc finger of groups I and II, Motif 3 the C2HC motif of group III, and Motifs 4 and 5 the second WRKY domain characteristic of group I. The 35 genes were distributed unevenly across the six chromosomes of chia's compact genome: chromosome 3 carried the largest share with 12 genes, chromosome 1 held 7, chromosomes 2 and 5 each carried 5, and chromosomes 4 and 6 each carried 3.

Duplication analysis identified 16 duplicated pairs, of which 13 arose from segmental duplication, the large-scale duplication of chromosomal blocks, and 3 from tandem duplication, where genes sit side by side. For every one of the 16 pairs, the ratio of nonsynonymous to synonymous substitution rates, Ka/Ks, was below 1, ranging from 0.207 to 0.642 with a mean of 0.4176. This uniform pattern indicates that the entire duplicated set has evolved under strong purifying selection, which removes harmful mutations, and no pair showed any evidence of positive selection.

The authors calculated a mean divergence time for the duplicated pairs of 23.05 million years, with a standard deviation of 13.54 million years. Secondary-structure predictions showed that the alpha helix is the dominant structural element in all 35 proteins, averaging 41.65 percent of residues, followed by random coil at 38.56 percent, extended strand at 13.38 percent and beta turn at 6.41 percent. The extended-strand component is functionally meaningful, because it reflects the four-stranded beta sheet that forms the DNA-binding surface through which WRKY domains recognize W-box elements.

A maximum-likelihood phylogenetic tree, built with 1000 bootstrap replicates, resolved the 35 CsWRKY genes together with orthologs from the five reference species into the seven established groups. In most groups, the chia genes clustered with their fellow eudicots, Arabidopsis, soybean and Medicago, rather than with the monocot rice, a pattern consistent with the botanical position of Salvia hispanica within the eudicot family Lamiaceae. Analysis of 2000-base-pair promoter regions upstream of each gene uncovered 20 classes of cis-acting regulatory elements, revealing that the CsWRKY genes are embedded in dense hormone- and stress-signaling networks.

The promoters carried W-box elements, which enable WRKY autoregulation, along with MYB-binding sites, stress-responsive STRE elements, TC-rich repeats, and hormone-responsive motifs including ABRE for abscisic acid, ERE for ethylene, GARE for gibberellin, TGA for auxin and the TCA element for salicylate, as well as light-responsive and circadian elements. These promoter features are predictive rather than demonstrative, but they suggest candidate roles in integrating hormonal and environmental cues. The expression results showed that most CsWRKY genes are constitutively expressed with tissue-preferential patterns, and one-way analysis of variance confirmed statistically significant differences in expression among tissues.

Principal component analysis and hierarchical clustering of the tissue-expression matrix separated these tissue-preferential patterns along the first two components. Notably, the authors state that no salinity, drought, cold, heat or abscisic-acid stress RNA-Seq dataset is currently available for chia in the public archives, so no stress-responsive expression analysis was performed, and earlier draft claims of stress-responsive expression were removed because they were not supported by real data. Statistical checks reinforced the internal consistency of the dataset: protein length and molecular weight were almost perfectly correlated, with a Pearson coefficient of 0.9995, and a one-sample t-test confirmed that the mean Ka/Ks of 0.4176 deviates significantly from the neutral expectation of 1, supporting purifying selection.

Homology models were built with I-TASSER for one representative CsWRKY protein from each of the seven phylogenetic groups, and model quality was validated with TM-scores ranging from 0.72 to 0.90, where values above 0.5 indicate a correct fold, and more than 89 percent of residues occupying favoured regions of the Ramachandran plot. Docking predicted binding in the WRKY-domain cleft for all representatives, with CsWRKY3 showing the strongest predicted binding at a free energy of minus 7.00 kilocalories per mole and CsWRKY17 the weakest. The authors present this as an exploratory structural hypothesis, not evidence of a physiological interaction, and recommend confirmation by techniques such as isothermal titration calorimetry or fluorescence binding assays.

The study's limitations are stated with unusual candor. No wet-laboratory experiments were performed to validate gene identity, expression patterns, protein-DNA interactions or protein-ligand binding. All analyses rest on a single reference genome assembly, so cultivar-level genetic variation in chia was not captured.

What the study delivers, therefore, is a preliminary but coherent genomic framework: 35 named and classified WRKY genes, their chromosomal addresses, duplication history, promoter architecture, tissue-expression profiles and a shortlist of candidates, led by CsWRKY3, for future functional work.

Subject of Research: Genome-wide in silico identification and expression profiling of the WRKY transcription factor gene family in Salvia hispanica (chia)

Article Title: Genome wide identification and in silico mRNA expression profiling of the WRKY transcription factor gene family in Salvia hispanica

Article References: Tiwari, S. (2026). Genome wide identification and in silico mRNA expression profiling of the WRKY transcription factor gene family in Salvia hispanica. Discover Plants, 3(1), Article 482. https://doi.org/10.1007/s44372-026-00972-8

Image Credits: AI Generated

DOI: 10.1007/s44372-026-00972-8

Keywords: chia, Salvia hispanica, WRKY transcription factors, genome-wide identification, in silico analysis, gene family, purifying selection, phylogenetics, RNA-Seq expression profiling, molecular docking, apigenin, Lamiaceae