Genome-Wide Identification of Bilberry WRKY Transcription Factors: Go Wild and Duplicate

WRKY transcription factor genes compose an important family of transcriptional regulators that are present in several plant species. According to previous studies, these genes can also perform important roles in bilberry (Vaccinium myrtillus L.) metabolism, making it essential to deepen our understanding of fruit ripening regulation and anthocyanin biosynthesis. In this context, the detailed characterization of these proteins will provide a comprehensive view of the functional features of VmWRKY genes in different plant organs and in response to different intensities of light. In this study, the investigation of the complete genome of the bilberry identified 76 VmWRKY genes that were evaluated and distributed in all twelve chromosomes. The proteins encoded by these genes were classified into four groups (I, II, III, and IV) based on their conserved domains and zinc finger domain types. Fifteen pairs of VmWRKY genes in segmental duplication and four pairs in tandem duplication were detected. A cis element analysis showed that all promoters of the VmWRKY genes contain at least one potential cis stress-response element. Differential expression analysis of RNA-seq data revealed that VmWRKY genes from bilberry show preferential or specific expression in samples. These findings provide an overview of the functional characterization of these proteins in bilberry.


Introduction
The bilberry (Vaccinium myrtillus L.), also known as the European wild blueberry, is a fruit that grows in the wild from the forests of Northern Europe to the Caucasus and towards the northern Asia-Pacific [1][2][3].V. myrtillus is one of the 400 species belonging to the genus Vaccinium of the Ericaceae family, which also includes cultivated blueberry species and wild fruits such as V. corymbosum, V. angustifolium, V. virgatum, and V. macrocarpon [4].Blueberries are usually diploid (2n = 2x = 24), but polyploid subspecies can be found in North America [2].Due to their high levels of beneficial nutrients and bioactive phytochemicals present in the skin, pulp, and leaves, especially anthocyanins, the demand for blueberries has been increasing in the market [3,5].Anthocyanins are responsible for their colors of red, blue, purple, and black fruits [5,6].In the last two decades, several studies have linked the polyphenols and high natural anthocyanin levels of blueberries to the prevention of chronic diseases such as cancer, coronary heart disease, and neurological disorders [5].
In addition, there has been an increased interest in the Vaccinium genus, and as a result, genomic and transcriptomic resources for this genus have been developed.However, Plants 2023, 12, 3176 2 of 20 despite this progress, there is still incomplete genome assembly and annotation data, as is the case of blueberries (V.corymbosum) [7].On the other hand, efforts to understand the genomic regions that control the anthocyanin composition of the berry have led to the generation of the complete genome, as is the case with bilberry (V.myrtillus) [4].Anthocyanins and flavonoid compounds present in blueberry fruits have been associated with structural genes and transcription factors (TFs) in many species [2].Advances in genomic resources and the understanding of Vaccinium spp.have made it possible to describe the MYBA TF loci and identify the main regulatory genes of this family that determine anthocyanin production [4].Recent studies have shown that anthocyanin profiles in Vaccinium spp.are regulated by prevailing light and temperature conditions [8,9].Blue light induces anthocyanin accumulation in pear fruit (Pyrus communis x pyrifolia cv Red Zaosu) and activates anthocyanin-production-responsive genes [10].In blueberry, the light intensity regulation of anthocyanin accumulation represents a valuable data set to guide future functional and crop improvement studies [11].In bilberry, it was demonstrated that blue and red/blue lights have the ability to promote anthocyanin biosynthesis by inducing the expression of key structural genes and accumulation of metabolites involved in the anthocyanin synthesis pathway, as well as the relationship between photosynthesis under different light qualities in blueberry leaves [12].Several regulatory genes, including WRKY family TFs, have been proposed to control fruit development and ripening processes, as well as their relationship with responses to various stresses [13,14].
The WRKY TFs are one of the largest families of transcriptional regulators in plants and are present in several plant species [15,16].These TFs are characterized by the WRKY domain, which contains about 60 amino acids that may be located near the N-or C-terminal regions and is composed of a highly conserved WRKYGQK heptapeptide DNA-binding sequence, followed by a zinc finger motif that binds to specific cis-regulatory elements called W-boxes (TTGACT/C) [17][18][19][20][21]. Evidence suggests that W-boxes are present in the promoter region of genes related to plant innate immunity [20].Based on the number of WRKY domains and the pattern of zinc finger motifs, WRKY genes can be divided into four groups (I, II, III, and IV) [22][23][24].Group I WRKYs have two domains containing a C2H2 zinc finger motif.Group II representatives have only one WRKY domain and one C2H2 zinc finger motif and can be divided into five subgroups (IIa, IIb, IIc, IId, and IIe).Group III WRKYs also have a single domain, but their zinc finger motif is C2HC [22,23] Recently, Group IV was identified, composed of incomplete WRKYs that do not fit into the previous groups because they lack the zinc finger motif [23,24].
We can find several examples of WRKY genes that act in the biosynthesis of anthocyanins, plant pigments that are present in leaves, flowers, fruits, and roots.Ultraviolet-B (UV-B) radiation promotes the synthesis of anthocyanins in many plants, and several transcription factors respond to UV-B radiation [25].In M. domestica, MdWRKY72 TF increases anthocyanin synthesis by direct and indirect mechanisms when induced by UV-B radiation [26].The overexpression of MdWRKY75 led to an increase in anthocyanin levels by binding to the MYB transcription factor promoter, MdMYB1 [27].Other examples are the MscWRKY12 and MscWRKY19 TFs from M. sylvestris, which are involved in leaf pigmentation during the development of this organ [24], and the involvement of WRKY TFs in the synthesis of anthocyanins in Raphanus sativus L., the carmine radish [28].Especially in fruits, such as Pyrus spp., their red skin accumulates anthocyanin, and the genes involved can contribute to the improvement of their appearance, as is the case with the combination of PyWRKY26 and PybHLH3 capable of co-directing the PyMYB114 promoter to generate anthocyanin accumulation in red pears [29].Eight PcWRKYs also participate in color development in red fruits of P. communis, and color fading in some fruits is due to reduced biosynthesis, increased degradation, and the suppression of anthocyanin transport [30].In addition, by having a better understanding of the genes involved in this biosynthesis, it is possible to overexpress them.In M. domestica, the overexpression of MdWRKY11 [31] and MdWRKY71-L [32] plays a role in the synthesis of anthocyanins in red fruits and demonstrates that TFs participate in the UV-B signaling pathway to regulate the accumulation of anthocyanin in apple.Just like in the apple, the red coloration of mango (Mangifera indica L.) peels results from anthocyanin accumulation, where it was found that MiWRKY1 and MiWRKY81 were upregulated during the light induction accumulation of anthocyanin in mango, indicating that these genes may regulate anthocyanin biosynthesis [33].
The growing interest in the nutraceutical properties of blueberries makes it essential to deepen our understanding of fruit ripening regulation and anthocyanin biosynthesis.In this context, WRKY transcription factors have been identified as a promising tool for manipulating stress tolerance and the synthesis of nutritionally desirable compounds.Therefore, it is crucial to identify the WRKY proteins, from their structural classification to their phylogenetic relationships, as well as their gene structure, conservation of domains and motifs, cis elements, chromosomal mapping, tandem and segmental duplication, and genetic divergence in the complete bilberry genome (Vaccinium myrtillus L.).The detailed characterization of these proteins will provide a comprehensive view of the evolution and modification of the WRKY protein family in the crop and will help to determine the functional features of VmWRKY genes in different plant organs and in response to different intensities of light.

Analysis of Cis Elements in VmWRKY Gene Promoters
The analyses of promoter sequences of 69 VmWRKY genes generated 108 types of putative cis-regulatory elements, which were categorized into four known and one unknown cis-regulatory action elements (Table S4).The largest category among the known elements comprises light-responsive elements (25%), predominantly represented by cis-action G-Box, Box 4, and GT1-motif.Following that, hormone-responsive elements (9%) with cis-action ABRE, CGTCA-motif, and TGACG-motif were identified.Elements related to development (9%) were found with cis-action CAT-box, O2-site, and GCN4_motif.Furthermore, environment-responsive elements (8%) exhibited cis-action AREM, LTR, and MBS.The elements associated with promoters and binding sites (7%) were predominantly represented by cis-action CAAT-box, TATA-box, and W-box.Notably, the WRKY transcription factors are bound to the W-box to initiate transcription.Additionally, the analysis revealed several other unknown functional elements (40%) with cis-action MYB, MYC, and STRE.The unrooted phylogenetic tree of 69 VmWRKY protein sequences (Figure 1A) displays their relationships, forming eight subgroups (I, IIa, IIb, IIc, IId, IIe, III, and IV).The largest subgroups are group I (14 members) and group IIc ( 14), followed by group III (12), group IIe (10), group IIb (9), group IId (4), group IIa (4), and group IV (2).In group IV, the protein sequences contain the conserved WRKY domain but lack the zinc finger motif.The primary structure of the 69 VmWRKY gene sequences (Figure 1B) shows variation in the number of introns, ranging from 0 to 9. By conducting the analysis of conserved functional motifs in the protein sequences, a total of 20 functional motifs were identified, distributed across each VmWRKY subgroup (Figure 1C).Each subgroup has specific motif patterns (Figures S1 and S2).

Chromosomal Localization, VmWRKY Gene Duplication, and Divergence
A total of 64 VmWRKY genes were mapped onto twelve chromosomes of the bilberry genome (Figure 2), and the genes VmWRKY65, VmWRKY66, VmWRKY67, VmWRKY68, and VmWRKY69 are sequenced at the scaffold level, and their location on the chromosome is unknown.The highest number of VmWRKY genes was found on chromosome 12 (11 genes), followed by chromosome 8 (9 genes), chromosome 4 (8 genes), chromosomes 1 and 6 (7 genes), chromosome 3 (5 genes), chromosomes 11 and 9 (4 genes), chromosome 2 (3 genes), and chromosomes 10, 5, and 7 (2 genes).Tandem and segmental duplications were identified (Table 2).Four pairs of VmWRKY genes were identified as tandem duplications (<100 Kb and >70% similarity) (Figure S3), and 12 pairs of VmWRKY genes were identified as segmental duplications (>70% similarity) (Figure S4).This information provides valuable insights into the genetic evolution of the gene family and its domains regarding other model species such as Arabidopsis thaliana (Figure S5).The estimation of divergence time for paralogous VmWRKY gene pairs, based on their synonymous substitution rates, ranged from 0.13 to 83.40 million years for the VmWRKY40/VmWRKY41 and VmWRKY61/VmWRKY63 gene pairs.Other duplicated gene pairs did not yield divergence time estimates with the methods applied using DNA sequences, coding sequences (CDS), and the transcriptome annotation of the genome.Another analysis was applied, and the Nei and Gojobori (NG) method was utilized to determine the divergence between proteins (Figure S6).

Identification of the WRKY Protein Family in Bilberry
Transcription factors exist in the form of a superfamily of genes and play a critical regulatory role in plants' growth, development, and response to the environment [34].In recent years, there has been a significant increase in interest in investigating the WRKY protein family in different plant species, including the bilberry crop.With the publication of the complete genome of several species, the identification and analysis of TF families at the whole genome level has become one of the focuses of genomic research.The WRKY TF family plays an important role in plant growth and development and defense mechanisms [35].

Analysis of Cis Elements in VmWRKY Gene Promoters
A particular feature of the WRKY TF family is the ability to specifically bind to W-box ((C/T)TGAC(T/C)).Nevertheless, they bind also to other cis-acting elements located in the promoter region of their target genes [19,36,39].In this study, it was possible to observe that the promoter region of VmWRKY genes contains a total of 108 conserved cis-regulatory elements with diverse functions.Strikingly, 40 VmWRKYs presented one or more W-boxes, indicating the putative autoregulation of these TFs, a number that is higher than the 26 FaWRKYs observed in strawberry [23,[40][41][42].Common promoter and enhancer regions, such as A-Box and CAAT-box, were also identified, as well as TATA-box, which are regions located around −30 of transcription initiation.In addition, the presence of HD-Zip 1 and HD-Zip 3, the elements involved in mesophyll cell differentiation and the site where protein binding occurs, were detected.
The presence of a large number of phytohormone-responsive elements, such as TGA (responsive to auxin), TATC-box (response to gibberellin), SARE (salicylic acid-responsive element), TCA (response to salicylic acid), ABRE (response to abscisic acid), AuxRR-core (responsive to auxin), CGTCA motif (response to MeJA), TGACG motif (response to MeJA), GARE (responsive to gibberellin), and P-box (gibberellin-responsive element), was detected.Additionally, 27 promoter regions associated with the response to light were identified, with G-Box, Box 4, GT1-motif, and TCT-motif being the most abundant within the species.This coincides with what was reported in raspberry (Rubus occidentalis), where the most abundant promoter region was Gbox, ATC, and TCT-motif [37].
The VmWRKY gene promoter regions also showed conserved cis-regulatory elements, which are involved in a variety of functions, such as responses to biotic and abiotic stresses.Among these elements are TC-rich repeats, which act in defense and response to stress; LTRs, which respond to low temperatures; AREs, which act in anaerobic induction; GC motifs, which act in specific anoxic induction; MBSs, which act in response to water stress; the WUN motif, which responds to tissue damage; and AT-rich sequences and AT-rich elements.The presence of diverse cis-acting elements, which mediate responses to environmental stresses and plant hormones, suggests that these WRKY TFs are involved in various biological processes.

Phylogeny, Gene Structure, and Motif Analysis of WRKY Proteins in Bilberry
Plants have adaptation mechanisms to adverse environments that involve signal transduction and molecular regulation in response to stress.TFs play a key role in this process by activating or inhibiting gene transcription by their specific binding to gene promoter regions.This results in the induction of functional gene expression and contributes to signal transduction in response to stress [43,44].
In the phylogenetic analysis of TFs, the 69 VmWRKY proteins were distributed into seven clusters.Cluster 1 presents nine WRKY proteins, cluster 2 presents four, cluster 3 presents ten, cluster 4 presents four, cluster 5 presents fourteen, cluster 6 presents sixteen, and cluster 7 presents twelve.The proteins belonging to each group (I, II, III, and IV), with group II subdivided into five distinct subgroups (IIa-e), formed clusters with members of only one group, and others were made up of members of two or more groups.During the intragroup evolutionary analysis of VmWRKY genes, it was observed that genes from group IIc were more divergent from the other members within group II (Iia, Iib, Iid, and IIe).Within this branch, IIa and IIb were grouped separately in a branch within group II, similarly to what was found in plum [45] and black raspberry [37].Also, IId and IIe formed another subgroup.Group I presented almost all members of this group clustered together, with two additional members from group III and group IV.The last group formed was composed of almost all members from group III but included one genotype from group IV (Figure 1).All WRKY TFs from group I in bilberry contain two WRKY domains.Among them, 15 VmWRKYs (21.7%) were assigned to group I, and 47 VmWRKY (68.12%) genes were distributed in group II, which was further classified into five subgroups, IIa, IIb, IIc, IId, and IIe, which contained 7, 9, 14, 5, and 12 VmWRKYs, respectively.The remaining 12 VmWRKYs (17.4%) belonged to group III, and 2 VmWRKYs (2.9%) belonged to group IV.These results are similar to those found in grape (Vitis vinifera) and blackberry (Rubus occidentalis), in which the number of WRKYs genes found in group II were the most abundant with 39 VvWRKYs (66.1%) [36] and 25 RoWRKYs (41.6%), respectively [37].
Group III of WRKY gene members presents a zinc finger motif different from groups I and II and can be considered the most dynamic in terms of evolution [34,46].In this study, 12 VmWRKYs were identified as group III, which is similar to the 13 AtWRKYs found in Arabidopsis [47] and the 10 found in wild strawberry [48].However, this count is higher than the 6 VvWRKYs found in grapes [36] and lower than the 28 OsWRKYs in rice [34].Some WRKYs genes from group III are part of the signaling pathway of the plant's defense system, having a significant impact on resistance to diseases and drought [34].
Finally, we have VmWRKY68 and VmWRKY69, which were grouped in group IV for not fitting into any of the other groups because the WRKY domain has a partial or lacks a complete zinc finger motif structure.WRKY genes belonging to group IV were also reported in species such as Arabidopsis thaliana with two AtWRKY genes [22] and Pennisetum glaucum with twelve PgWRKY genes [41].Likewise, M. domestica (10% = 13 MdWRKYs) and Vitis vinifera (2% = 2 VvWRKYs) [49] presented WRKY proteins without a complete domain.This occurrence may mean the loss of the WRKY domain [21,50,51] and modifying the functional properties of the proteins [42].In studies performed in Hylocereus undatus, the absence of the WRKY domain, the zinc finger or coiled-coil sequence, did not allow the binding between the WRKY protein and the promoter region of the gene [52], suggesting that genes belonging to group IV can originate non-functional proteins.

Chromosomal Localization, VmWRKY Gene Duplication, and Divergence
Overall, 69 VmWRKY genes were identified in the bilberry genome; however, only 64 VmWRKY genes (VmWRKY1-VmWRKY 64) had a known location on the chromosomes, and they were distributed into twelve chromosomes.The WRKY genes (VmWRKY65−VmWRKY69) are sequenced at the scaffold level and have no known chromosome locations.And most of the VmWRKYs were abundant on Chr 12.The number of VmWRKY genes is similar to that of the pitaya (Hylocereus undatus), with 70 HuWRKY genes distributed on eleven chromosomes [52].The raspberry also presents a similar number of WRKY genes (60 RoWRKY); however, it has a lower number of chromosomes (7), and Chr 6 had the largest number of RoWRKY genes, representing 23.33% of the WRKY genes [37].In the case of grapes (Vitis vinifera), the number of WRKY genes present in this species is also similar (59 VvWRKYs); however, they are mapped to nineteen chromosomes [38], a greater number of chromosomes than that present in bilberry.The differential distribution of the WRKY genes present in a species may imply chromosomal rearrangements and duplication events that took place during its evolution [41].
It has been observed that gene duplication, including tandem duplication, fragment duplication, and genome duplication, is a key factor in gene family amplification in plant genomes [53].In this particular study, 15 segmental duplication events and four tandem duplication events (VmWRKY37/VmWRKY38, VmWRKY40/VmWRKY41, VmWRKY46/VmWRKY47, and VmWRKY48/VmWRKY49) were identified (Table 2).A high frequency of segmental duplication was also observed in Rubus occidentalis, with five genes containing homologous segments [37].In Fragaria vesca, segmental duplications were higher than tandem duplications, representing 84.2% and 15.8% of the total duplications, respectively.Compared to what has been observed in V. vinifera [54] and Oryza rufipogon [55], tandem duplication events can significantly contribute to the expansion of VmWRKY genes [52].

Expression Patterns of WRKY Genes by Induction of Light and Plant Organs
The processing of transcriptomic samples from bilberry plant organs induced by light showed variation in the total gene expression and specific expression of VmWRKY genes in each sample (Table S5).In the case of samples treated with red light (24,491 = 67%), gene expression was higher than in the control and blue light samples.This pattern was also observed in the specific expression of VmWRKY genes.However, 17 VmWRKY genes showed differential expression under red light compared to the control and blue light (Figure 3(A3)).On the other hand, in pear (P.communis x pyrifolia cv 'Red Zaosu'), exposure to blue light led to an increase in anthocyanin accumulation and the activation of genes responsible for anthocyanin production [10].In red mango fruit (M.indica L.), UV-B light induction positively regulates anthocyanin accumulation, and the genes MiWRKY1 and MiWRKY81 are involved in this regulation [33].In general, anthocyanin accumulation is influenced by light availability, and the specific impact of different light qualities varies among plant species [56][57][58].These results suggest that VmWRKYs are concentrated in the skin and berry of the fruit, which may indicate an association with higher levels of anthocyanins, the main source of organoleptic and antioxidant properties in blueberries, as reported in a gene expression study during bilberry development [59].However, the expression level in bilberry leaves was high compared to that observed in Tetrastigma (Tetrastigma hemsleyanum), a grape family plant with eight WRKYs, which showed high expression levels in leaves [60].
Gene expression in plant organ samples showed similar levels in the entire berry (24,268 = 67%) and the berry skin (24,449 = 67%), which were comparable to those observed in the berry pulp (24,190 = 66%) and higher than those numbers found in leaves (22,791 = 63%) and roots (23,370 = 64%) (Table S6).However, a greater number of VmWRKY genes were expressed in the root sample (62 = 0.17%), and the differential expression of 41 VmWRKY genes was high compared to the rest of the samples (Figure 3(B3)).The higher number of genes and their differential expression in the root sample could be explained by the role played by WRKY transcription factors, which are mainly involved in development and stress responses, such as salt and water stress tolerance [61].Since bilberries are sensitive to these conditions, it is important to pay close attention to the presence of VmWRKYs in directly affected organs such as roots, where VmWRKY54 and VmWRKY11 showed higher expression levels (Figure 3(B4)).This is supported by examples of the effectiveness of these transcription factors in Arabidopsis, where the overexpression of AtWRKY46, GmWRKY13, or VvWRKY11 can positively regulate salt and water stress tolerance [62][63][64][65], and in Nicotiana benthamiana, where the overexpression of GhWRKY41 conferred tolerance to water and salt stress [66].
The complete genome of bilberry (Vaccinium myrtillus L.) allowed the characterization of 69 members of the VmWRKY gene family.Segmental and tandem duplications were detected and could enhance biotic/abiotic resistance in the bilberry genome.The average ages of duplications were identified as 8.27 mya (range 013-16.41)for the tandem and 26.43 mya (range 0.48-83.40)for the segmental duplications, suggesting more recent events for tandem duplications than segmental duplications.

Multiple Alignment and Phylogenetic Analyses
A phylogenetic tree was constructed to compare bilberry WRKY proteins.Multiple alignment of WRKY protein sequences was performed with ClustalW software using standard parameters [73].The phylogenetic tree was constructed using BEAST v.2.5 software [74], with the UPGMA clustering method [75].A bootstrap analysis was conducted using 10,000,000 replicates, and evolutionary distances were calculated using the JTT matrix-based method [76].

Analysis of Gene Structure and Identification of Conservation of Motifs
To investigate the diversity and structure of VmWRKY family members, genomic sequences for their exon/intron were used and plotted on TBtools [77], based on bilberry genome annotation information (Vaccinium myrtillus isolate NK2018 v1.0 genome sequence).VmWRKY protein sequences were used to identify conserved motifs using the Expectation Maximization Tool for Motive Elicitation MEME (https://meme-suite.org/meme/tools/ meme) (accessed on 20 February 2023) [78].The parameters were as follows: number of repetitions: any; maximum number of motifs: 20; optimal motif widths: 8 to 50 amino acid residues.

Chromosomal Localization, Gene Duplication, Ka/Ks Calculation, and Divergence Time Estimation
The chromosomal location image of the VmWRKY genes was generated by the TBtools software [77], according to chromosomal position information provided in the genomic database for Vaccinium-GDV (https://www.vaccinium.org/)(accessed on 30 February 2023).To identify specific tandem duplications of VmWRKYs, the following criteria were used: genes within a 100 kbp region on an individual chromosome with a sequence similarity of ≥70% [79].The pairwise local alignment calculation of two protein sequences was performed by the Smith-Waterman algorithm of EMBOSS Water (http://www.ebi.ac.uk/Tools/psa/) (accessed on 17 February 2023) [80].
For the calculation of non-synonymous substitutions per site (Ka) and the number of synonymous substitutions per synonymous site (Ks), in addition to comparing the selection pressure, a Ka/Ks ratio greater than 1, less than 1, and equal to 1 represent positive, negative, and neutral selection, respectively.For each pair of genes, the value of Ks was used to estimate the time of divergence in millions of years based on a rate of 6.1 × 10 −9 replacements per site per year, and the time of divergence (T) was calculated as T = Ks/(2 × 6.1 × 10 −9 ) × 10 −6 million years ago (Mya) [81].The bioinformatics tool used for genetic divergence was the simple Ka/Ks Calculator (NG) from TBtools-II [77].
Data processing was performed in the following three steps.(a) Quality control and adjustment of samples: SRA toolkit [83] was used to download the data samples, FastQC [84] was employed to analyze and visualize the quality of readings, and Trimmomatic ver.0.39 [85] was applied to remove the low quality and library adaptors.(b) The reads were mapped against the bilberry (V.myrtillus) reference genome [4] using the software HISAT2 [86].(c) For the counting of total reads aligned by gene in the different libraries, FeatureCounts [87] was used.The quantification analysis (d) was performed using the packages limma, edgeR, and DESeq2 in R software [88].In this protocol, a FPKM normalization method (fragments per kilo base per million mapped reads) was used, which are counts scaled by the total number of reads and the expression of VmWRKY genes per library.The analysis of gene expression proportion was based on the number of genes annotated in the bilberry genome (36,405 genes = 100%), and the expression of each transcriptome sample and VmWRKY genes was proportional to the total annotated genome genes.A heatmap was produced using TBTools with normalized data (log2 counts per FPKM) [77].Then, a multi-dimensional scaling (MDS) plot was generated to check the repeatability of the sample and the overall difference between samples.Additionally, a MeanVar plot and Biological Coefficient of Variation (BCV) were calculated.

Conclusions
The complete genome of bilberry (Vaccinium myrtillus L.) contains 76 identified and 69 characterized members of the VmWRKY gene family, which are located on twelve chromosomes.Additionally, 12 segmental and four tandem duplications were identified in the bilberry genome.The cis-regulatory elements found in the promoters of VmWRKY genes are putatively involved in various functions related to biotic and abiotic stress responses.The differential expression of VmWRKY genes was induced by light, suggesting its involvement in anthocyanin biosynthesis.The characterization of VmWRKY genes in the bilberry genome, their phylogenetic relationships, and differential expression are important for understanding their role in regulating anthocyanin biosynthesis and their adaptation to different environmental conditions.

Figure 1 .
Figure 1.Phylogenetic relationships and structure of genes encoding the bilberry VmWRKY proteins: (A) The unrooted tree was generated with the BEAST program using the full-length amino acid sequences of the 69 bilberry VmWRKY proteins by the UPGMA method, with 1,000,000 bootstrap replications.VmWRKY protein subfamilies (I, IIa, IIb, IIc, IId, IIe, III, and IV).(B) Exon/intron organization of bilberry VmWRKY genes.Yellow boxes represent exons, and black lines represent introns.Untranslated regions (UTRs) are indicated by blue boxes.The sizes of exons and introns can be estimated using the scale at the bottom.(C) Schematic representation of conserved motifs in bilberry VmWRKY proteins, elucidated from publicly available data (NCBI CDD Domain-Pfam-18,271 PSSMs).Each colored rectangular box represents a motif with the given name and motif consensus.

Figure 1 .
Figure 1.Phylogenetic relationships and structure of genes encoding the bilberry VmWRKY proteins: (A) The unrooted tree was generated with the BEAST program using the full-length amino acid sequences of the 69 bilberry VmWRKY proteins by the UPGMA method, with 1,000,000 bootstrap replications.VmWRKY protein subfamilies (I, IIa, IIb, IIc, IId, IIe, III, and IV).(B) Exon/intron organization of bilberry VmWRKY genes.Yellow boxes represent exons, and black lines represent introns.Untranslated regions (UTRs) are indicated by blue boxes.The sizes of exons and introns can be estimated using the scale at the bottom.(C) Schematic representation of conserved motifs in bilberry VmWRKY proteins, elucidated from publicly available data (NCBI CDD Domain-Pfam-18,271 PSSMs).Each colored rectangular box represents a motif with the given name and motif consensus.

Plants 2023 ,
12, 3176 9 of 20 divergence time for paralogous VmWRKY gene pairs, based on their synonymous substitution rates, ranged from 0.13 to 83.40 million years for the VmWRKY40/VmWRKY41

Figure 2 .
Figure 2. Chromosomal map and coordinates of VmWRKY gene duplication events: The identity of each linkage group is indicated in the central part of each bar.The putative segmental duplicated genes are connected by red color lines, and the duplicated gene pair in tandem is represented by the blue color on the chromosome.

Figure 2 .
Figure 2. Chromosomal map and coordinates of VmWRKY gene duplication events: The identity of each linkage group is indicated in the central part of each bar.The putative segmental duplicated genes are connected by red color lines, and the duplicated gene pair in tandem is represented by the blue color on the chromosome.

Figure 3 .
Figure 3. (A) Gene expression of bilberry transcriptomic light samples (under control, red light, and blue light conditions, measured in FPKM): (1) multi-dimensional scaling (MDS) plot of the light samples; (2) VmWRKY gene expression common to each sample shown in a Venn diagram; (3) differential expression of VmWRKY genes in each sample (values presented in log2, ranging from 0.40 to 2.30); (4) highly expressed VmWRKY genes in all samples (values presented in log2, ranging from 0.35 to 5.66).(B) Gene expression of bilberry transcriptomic organ samples (in leaves, roots, whole berries, berry flesh, and berry skin, measured in FPKM): (1) multi-dimensional scaling (MDS) plot of the organ samples (roots, leaves, berry skin, whole berry, and berry flesh); (2) VmWRKY gene expression common to each sample shown in a Venn diagram; (3) differential expression of VmWRKY genes in each sample (values presented in log2, ranging from −2.00 to 2.00); (4) highly expressed VmWRKY genes in all samples (values presented in log2, ranging from −1.50 to 3.00).

Figure 3 .
Figure 3. (A) Gene expression of bilberry transcriptomic light samples (under control, red light, and blue light conditions, measured in FPKM): (1) multi-dimensional scaling (MDS) plot of the light samples; (2) VmWRKY gene expression common to each sample shown in a Venn diagram; (3) differential expression of VmWRKY genes in each sample (values presented in log2, ranging from 0.40 to 2.30); (4) highly expressed VmWRKY genes in all samples (values presented in log2, ranging from 0.35 to 5.66).(B) Gene expression of bilberry transcriptomic organ samples (in leaves, roots, whole berries, berry flesh, and berry skin, measured in FPKM): (1) multi-dimensional scaling (MDS) plot of the organ samples (roots, leaves, berry skin, whole berry, and berry flesh); (2) VmWRKY gene expression common to each sample shown in a Venn diagram; (3) differential expression of VmWRKY genes in each sample (values presented in log2, ranging from −2.00 to 2.00); (4) highly expressed VmWRKY genes in all samples (values presented in log2, ranging from −1.50 to 3.00).

Table 1 .
Information on VmWRKY genes in the bilberry genome (Vaccinium myrtillus L.).

Table 2 .
Ks, Ka, and Ka/Ks calculation and divergence time of duplicated bilberry VmWRKY gene pairs.

Table 2 .
Ks , Ka, and Ka/Ks calculation and divergence time of duplicated bilberry VmWRKY gene pairs.