WRKY genes family study reveals tissue-specific and stress-responsive TFs in wild potato species

Wild potatoes, as dynamic resource adapted to various environmental conditions, represent a powerful and informative reservoir of genes useful for breeding efforts. WRKY transcription factors (TFs) are encoded by one of the largest families in plants and are involved in several biological processes such as growth and development, signal transduction, and plant defence against stress. In this study, 79 and 84 genes encoding putative WRKY TFs have been identified in two wild potato relatives, Solanum commersonii and S. chacoense. Phylogenetic analysis of WRKY proteins divided ScWRKYs and SchWRKYs into three Groups and seven subGroups. Structural and phylogenetic comparative analyses suggested an interspecific variability of WRKYs. Analysis of gene expression profiles in different tissues and under various stresses allowed to select ScWRKY045 as a good candidate in wounding-response, ScWRKY055 as a bacterial infection triggered WRKY and ScWRKY023 as a multiple stress-responsive WRKY gene. Those WRKYs were further studied through interactome analysis allowing the identification of potential co-expression relationships between ScWRKYs/SchWRKYs and genes of various pathways. Overall, this study enabled the discrimination of WRKY genes that could be considered as potential candidates in both breeding programs and functional studies.

cultivated and wild species are increasing. Such investigations may pave the way into exploiting these regulators for breeding purposes. A recent study carried out in the sweet potato wild ancestor Ipomoea trifida, highlighted how investigations on WRKY gene family in wild relatives can boost the molecular breeding of cultivated species 13 . However, our knowledge is still not complete and therefore WRKY gene biodiversity remains unlocked in many species.
The potato, Solanum tuberosum, is one of the most cultivated non-cereal crop in the world. Its cultivation is often hampered by the fact that it is susceptible to a wide range of stressors causing severe yield losses. Sources of resistance can be found in its tuber-bearing wild relatives, that are highly used as rootstock for cultivated Solanaceae 14 but poorly used in breeding programs. However, recent technologies can be implemented to enhance this precious source of genes/alleles. Among them, genome sequences are opening new paths for both basic research and varietal development. Nowadays, the genome sequence of two wild potato species, S. commersonii and S. chacoense, are available 15,16 . These species are excellent sources of tolerance to both biotic stressors, such as Ralstonia solanacearum 17 , Phytophthora infestans 18 and Pectobacterium carotovorum 19 , and abiotic constraints, such as cold 15 and drought 20 . Despite this, to date no studies have examined WRKY gene family components and their different characteristics in wild potato species. A few data on this gene family are available only in the cultivated potato, where Zhang et al. 21 , Liu et al. 22 and Cheng et al. 12 identified 79, 82 and 81 StWRKYs, respectively. Previously, Dellagi et al. 23 identified StWRKY1 as a good candidate for functional studies, and Shahzad et al. 24 overexpressed it in potato. They provided evidence that StWRKY1 acts as positive regulator of biotic and abiotic stress resistance through the activation of basal defence networks. Here, for the first time, we report a detailed analysis of WRKY genes in the genome of S. commersonii and S. chacoense, providing subGroup classification, gene structure and conserved motif composition. We analysed the patterns of ScWRKYs and SchWRKYs expression in flowers, leaves and tubers to determine whether some WRKYs own tissue-specificity. Furthermore, we used S. commersonii to highlight expression changes of selected ScWRKY genes after wounding and biotic (Potato Virus Y and P. carotovorum) stresses. Through the data here presented, the work aims to give a picture of the potato wild WRKY members, their nature and the complexity of their responses to unfavourable situations.

Materials and Methods
Identification of WRKY in S. commersonii and S. chacoense and phylogenetic analysis. The well-known WRKY protein sequences of S. tuberosum 22 and A. thaliana 25 were used as queries to build an HMM profile through HMMER as reported by Esposito et al. 26 and to search orthologs in S. commersonii (cmm1T clone of PI243503) and S. chacoense (M6 clone) genomes. Only sequences with an e-value lower than 10 −5 and an identity higher than 55% were regarded as putative WRKYs and further analyzed. The full-lenght WRKY candidate proteins were then manually confirmed by checking the WRKY domain using the NCBI search domain online tool 26 and used for the phylogenetic analysis. Names were assigned based on S. tuberosum orthologs using bootstrap replicates of the Maximum Likelihood (ML) phylogenetic tree (values higher than 50). Briefly, MEGAX 27 was first used to establish the best-fit model of evolution through the option "Find best DNA/Protein Models" implemented in the program and then for phylogenetic tree building using the appropriate options. In the phylogenetic analysis were integrated seven AtWRKY proteins randomly selected as representative of each WRKY Group, as already reported by Karanja et al. 28 . One-to-one orthologs were considered when candidate proteins allocated on the same clade in the phylogenetic tree with S. tuberosum. The exon-intron organization of WRKY genes was determined using the online GSDS tool (http://gsds.cbi.pku.edu.cn). Finally, the on-line tool Phenogram (http://visualization.ritchielab.org/phenograms/plot) was used to determine the location of the WRKY genes on S. chacoense chromosomes.

Plant materials and stress treatments.
In-vitro plantlets of S. commersonii clone cmm1T, accession PI243503, derived from the Inter-Regional Potato Introduction Station (Sturgeon Bay, Wisconsin), were micro-propagated as described by D' Amelia et al. 30 . Four-week-old vitroplants were transplanted into 14-mm plastic pots containing sterile soil and grown in a greenhouse under long-day conditions (16- 31 . For assessing bacterial resistance, the protocol of Melito et al. 32 was used with few modifications. The stem base of vitroplants (one injection per plant) was inoculated with 20 μl of P. carotovorum strain Ecc 009 sospension under greenhouse conditions (with temperatures ranging from 20 to 30 °C during the day and from 12 to 17 °C during the night). The bacterial culture was adjusted to 10 6 CFU·mL −1 in MgCl 2 solution. The whole plant was then covered with a transparent plastic bag. For both treatments (viral and bacterial), plants inoculated with buffer were considered as mock control. At each time point, leaves were collected from three biological replicates, both for treated and untreated samples. Each biological replicate consisted of a pool of three plants. Young leaf samples were collected from treated and mock control plants following the time course and stored at − 80 °C before RNA extraction. Wounding stress was induced according to the protocol of Vannozzi et al. 33 with few modifications. Leaf discs (15 mm diameter) were punched from healthy leaves detached from glasshouse-grown plantlets and incubated upside down on 3MM moist filter paper in large Petri dishes at 22 °C under 12 h light / 12 h dark conditions until harvest. Collected discs were immediately frozen in liquid nitrogen and stored at −80 °C for subsequent RNA extraction. Five discs were randomly chosen per each time point. No treated leaves were used as control. Each treatment consisted of three biological replications.
RNA extraction, cDNA synthesis and quantitative Real-Time PCR (RT-qPCR). Total RNA was isolated from 100 mg of grinded leaves as reported by Rinaldi et al. 34 and Villano et al. 35 . The Spectrum TM Plant Total RNA kit (Sigma-Aldrich, St. Louis, MO, USA) was used following the manufacturer's protocol with some modifications. Quantity and quality of the isolated RNA was measured using the NanoDrop ND-1000 spectrophotometer (Thermo Scientific, Wilmington, DE, USA) and Qubit 2.0 fluorometer (Life Technologies, Carlsbad, CA). For cDNA synthesis, 1 µg of each RNA sample was reverse transcribed using the SuperScript III cDNA Synthesis Kit (Life Technologies, Paisley, UK) following the manufacturer's protocol. Specific primers were designed using the website Primer3 as reported by Koressaar et al. 36 (Supplementary Table 1). Expression analysis was conducted by RT-qPCR as reported by Di Meo et al. 37 and Brulè et al. 38

Protein-protein interaction in silico analyses. An interactome analysis was carried out to investigate
the function of tissue-specific and stress responsive ScWRKYs and SchWRKYs selected in the expression study through the analysis of direct ortholog of StWRKY genes. The protein-protein interaction networks STRING database was used (http://string-db.org/). It reports protein associations based on various sources, such as experimental results, pathway understanding, text-mining and genomic information 40 . The interactome was constructed using a medium confidence score (0.400).

Phylogenetic analysis and classification of ScWRKYs and SchWRKYs.
A total of 79 and 84 candidates corresponding to the Pfam WRKY family were distinguished in S. commersonii and S. chacoense, respectively (Table 1). Based on phylogenetic analysis, 71 ScWRKYs and 80 in SchWRKYs were identified as direct orthologs of StWRKYs, while the remaining were classified as not direct orthologs and named with the suffix -like. Two paralog genes of SchWRKY080, SchWRKY086 and ScWRKY087 were also identified and named with the suffix -a and -b ( Table 1). The phylogenetic analysis of seven AtWRKY proteins randomly selected as representative of each WRKY Group and all S. commersonii and S. chacoense WRKY proteins revealed ScWRKY and SchWRKY classification in three large Groups corresponding to Group I, II and III ( Fig. 1), with the exception of nine proteins in S. commersonii (ScWRKY047,ScWRKY051, ScWRKY052, ScWRKY055, ScWRKY085, ScWRKY087a, ScWRKY087b, ScWRKY088 and ScWRKY089) and eight proteins in S. chacoense (SchWRKY047, SchWRKY051, SchWRKY052, SchWRKY056, SchWRKY057, SchWRKY085, SchWRKY088 and SchWRKY089), that were not assigned to any Group (Table 1). In S. commersonii, 12 ScWRKY proteins belonged to Group I, 47 to Group II, and 10 to Group III. Group II proteins were further categorized into subGroups. Group IIa, IIb, IIc, IId and IIe included 5, 8, 13, 7 and 14 ScWRKYs respectively ( Table 1). As far as S. chacoense is concerned, 14 proteins belonging to Group I, 45 to Group II, and 15 to Group III were identified. Those of Group II were classified in subGroup IIa (5 SchWRKYs), IIb (5), IIc (15), IId (7) and IIe (12) ( Table 1). Gene and protein features, including the length of the protein sequence, the WRKY domain motif composition and the exons/introns number were analyzed and reported in Supplementary Table 2. In S. commersonii, the "WRKYGQK" pattern was highly conserved in 69 ScWRKYs, while five variations were observed in the other proteins ("WGKYGQK", "WRWLKCG", "WSKYGQK", "WRKCGQK", "WRKYGMK"). In S. chacoense, 74 SchWRKYs contained the "WRKYGQK" domain, while the other proteins contained one of the following variations: "WIKYGEN", "WHKYGQK", "WRKYGMK", "WKKHGSN", "WHKCGQK". Concerning the Zinc Finger motif, the most common pattern in both species was "C-X 4-5-7 -C-X 22-23-24 -H-X-H/C". The only exceptions were ScWRKY068 with "C-X 8 -C-X 27 -H-X 2 -H", ScWRKY074 with "C-X 1 -C-X 26 -H-X-C", and SchWRKY074 with "C-X 8 -C-X 24 -H-X-C". Regarding the number of WDs in the studied proteins, out of 12 members belonging to Group I in S. commersonii, eight contained two WDs, two had two WDs and other two possessed three WDs. All Group I members in S. chacoense harbored two WDs, except SchWRKY014 (one WD). Seven ScWRKYs belonging to Group II and two of Group III contained two WDs, while all other members had only one WD. In S. chacoense, Group II and III proteins harbored one WD. All Group III members contained the H X C Zinc Finger domain (Supplementary Table 2). Our analysis pointed out that the number of amino acids of ScWRKYs varied from 107 (ScWRKY30) to  Table 2). The exon-intron organization of our WRKY genes was examined to gain more insight into the evolution of the WRKY family in potato. As shown in Supplementary Table 2, all ScWRKY genes possessed from one to eight exons. A similar trend was observed in S. chacoense. Concerning the genomic localization of WRKY genes, due to the unavailability of S. commersonii physical map, we plotted genes only on S. chacoense chromosomes using the Phenogram on-line tool (http://visualization.ritchielab.org/phenograms/plot) ( Figure S1). Out of 84 SchWRKY genes identified, 83 were mapped. As represented in Figure S1, most of the genes were located on chromosome 3 (11 genes; 13.1%), followed by chromosome 5 (10; 11.9%), Unknown (8; 9.5%) and 2 (7; 8.3%). A total of 25 SchWRKY genes (5 on each chromosome) were localized on chromosomes 7 to 12, whereas no one was mapped on chromosome 11.
Four ScWRKY genes (ScWRKY016, ScWRKY023, ScWRKY045 and ScWRKY055) distributed in different Groups were selected to further investigate WRKYs behaviour in response to biotic (wounding) and abiotic (PVY and P. carotovorum) stressors using qRT-PCR (Fig. 3). The expression trend of our WRKYs was variable among and during treatments. In particular, wounding stress caused ScWRKY023 and ScWRKY045 overexpression during the whole treatment and ScWRKY055 overexpression at 4-and 6-hours post treatment (hpt). As for viral infection response, ScWRKY016 and ScWRKY045 were always downregulated, while the other genes were upregulated only at one of the five hpt. The bacterial inoculation with P. carotovorum did not activate ScWRKY016 Figure 1. Phylogenetic analysis WRKY proteins in S. commersonii, S. chacoense and seven representative proteins of Arabidopsis. Multiple sequence alignments of WRKY amino acid sequences were performed using ClustalX, and the phylogenetic tree was constructed using MEGAX by the Maximum Likelihood (ML) method and 1000 bootstrap replicates. The tree was divided into seven phylogenetic subGroups and distinguished by colours: dark purple for Group I, light blue for subGroup IIa, orange for subGroup IIb, light purple for subGroup IIc, dark blue for subGroup IId, green for subGroup IIe, red for Group III. The bootstrap values were ≥85%.

Scientific RepoRtS |
(2020) 10:7196 | https://doi.org/10.1038/s41598-020-63823-w www.nature.com/scientificreports www.nature.com/scientificreports/ and ScWRKY045, while the other two genes were upregulated at 2-and 6-hpt. Given the involvement of WRKYs in several biological processes, we wondered whether they might play roles under other stresses. Since WRKY expression data on wild potato species exposed to any stress are not available, we retrieved WRKYs RPKM values www.nature.com/scientificreports www.nature.com/scientificreports/ from S. tuberosum experiments involving several treatments and stressors. As shown in Figure S2, the transcription of most WRKY genes was affected by various treatments. Only StWRKY61 to StWRKY67 did not change their transcriptional activity upon stress. The late blight infection did not perturbate the expression of StWRKYs. StWRKY023, StWRKY044, StWRKY054 and StWRKY055 increased their expression following mannitol treatment, whereas ABA, IAA and GA3 hormonal treatments affected the transcriptional activity of 3 (StWRKY027, StWRKY028 and StWRKY046), 1 (StWRKY035) and 4 (StWRKY023, StWRKY054, StWRKY068, StWRKY070) S. tuberosum WRKYs, respectively. BABA and BTH treatments induced an overexpression of 18 and 15 StWRKYs respectively, of which StWRKY042, StWRKY075, StWRKY078 and StWRKY080 were in common. Concerning heat stress, 12 StWRKYs were overexpressed. Finally, under in-vitro culture conditions, 10 StWRKYs were overexpressed in shoots and one (StWRKY004) in roots.

In silico protein interaction network of selected ScWRKYs and SchWRKYs.
A network of interaction was studied for WRKYs showing either tissue-specific or stress-induced expression ( Figure S3a and S3b). The S. commersonii flower-specific expressed WRKY002 formed a node with the anthocyanins and cell differentiation regulatory proteins. STRING analyses provided evidence that ScWRKY002 interacts, among the others, with JAF13 and TTG1, two well-characterized potato anthocyanins bHLH and WD40 TFs 27,39 . Both the leaf-specific expressed ScWRKY042 and ScWRKY080 formed a cluster of interaction with a Leucine Rich Repeat (LRR) protein (an evolutionarily conserved protein associated with innate immunity in plants). The two wounding-responsive www.nature.com/scientificreports www.nature.com/scientificreports/ ScWRKY023 and ScWRKY045 established two independent nodes of interaction. The former set a cluster with "Wound-responsive Apetala2 like factor 2 (WRAF2)" (annotation for transcript PGSC0003DMT400021314 on SpudDB database), while ScWRKY045 interacted with a cluster of proteins linked to a class of glycosyltransferase. Concerning S. chacoense, SchWRKY030 (found to be flower-specific) interacted directly with eIF2B_5, a key protein involved in mRNA translation mechanisms. On the counterpart, the leaf-specific SchWRKY017, SchWRKY043, SchWRKY059 and SchWRKY077, together with the flower-specific SchWRKY028, showed the same interaction with LRR proteins already described for ScWRKY042 and ScWRKY080.

Discussion
Due to its importance in the regulation of several processes in plants 5 , WRKY family has been studied in more than 60 plant species. In Solanaceae, data are available in some important crops, such as S. tuberosum (79 21 , 82 22 and 81 12 WRKYs), S. lycopersicum (83 WRKYs 41 ) and S. melongena (50 WRKYs 42 ). However, no information is available on the number and structural variability of WRKY TFs in Solanaceae wild species, which represent an important reservoir of genetic variation for breeding. This study was set up with the aim to profile WRKY encoding genes in S. commersonii and S. chacoense, two noteworthy tuber-bearing potato species used in potato breeding programs 20,43-45 . Structural analysis of ScWRKYs and SchWRKYs revealed interspecific diversification. The recently published genome annotation of S. commersonii 15 and S. chacoense 16 enables a comprehensive investigation of the WRKY family. We detected 79 and 84 genes encoding putative WRKY TFs in S. commersonii and S. chacoense, respectively. These results indicate that, compared to the cultivated potato 22 , S. commersonii possesses a lower number and S. chacoense a higher number of WRKY genes. Both species displayed a number of WRKYs greater than that of barley (45) 46 , castor bean (58) 47 , cucumber (55) 48 , rapeseed (43) 49 and grapevine (59) 50 , and lower than that of cotton (120) 51 , maize (136) 52 , soybean (131) 53 and rice (100) 25 . From this comparison, it appears that the number of WRKY encoding genes is not proportional to the genome size of the respective plant species, as also reported by Waqas et al. 54 . ScWRKY and SchWRKY proteins were primarily divided into three main phylogenetic Groups with Group II further classified into five subGroups (IIa-IIe). Most of WRKYs found in the two wild species belonged to Group II and this is in line with results obtained in S. tuberosum 22 . As known, WRKY proteins are characterized by one or more WRKY domain. In this study, we found that ScWRKYs and SchWRKYs had either one or two WDs. Interestingly, two ScWRKYs (ScWRKY010 and ScWRKY002) carried three WDs. This might be the result of the acquisition of a WRKY domain during evolution, supporting findings of Aversano et al. 15 and Esposito et al. 31,55 , who reported that S. commersonii prosper lineage-specific segmental duplications during evolution. Not only WDs number, but also WDs structural divergences identified in S. commersonii and S. chacoense might be the consequence of mutations during the process of evolution. Almost all ScWRKYs and SchWRKYs contained the highly conserved heptapeptide WRKYGQK motif, except for eight variants. Among them, WGKYGQK of ScWRKY014, WRWLKCG of ScWRKY061, WHKYGQK of SchWRKY014 and WKKHGSN in SchWRKY057 were not found in any other species. On the counterpart, the remaining variants were identified also in S. tuberosum 22 , S. lycopersicum 56 , H. vulgare 57 and C. annum 58 . Zhou et al. 59 hypothesized that these variations may change the DNA targets' binding specificity. The structural diversity has been investigated also at the genomic level through the identification of exons and introns. As reported by Shiu and Bleecker 60 , this can highlight events of diversification and neo-functionalization of WRKY genes. In contrast to findings by Wang et al. 61 , our results did not reveal a conservation of gene structure among the members of the same Group, even though they allowed the discrimination of eight intron-lacking WRKYs (two ScWRKYs and six SchWRKYs). This is in agreement with results reported in the cultivated potato, where StWRKY23 and StWRKY24 had no introns 21 . Lynch et al. 62 hypothesized that the intron turnover can be the result of reverse transcription of the mature mRNA followed by homologous recombination with intron-containing alleles.
Identification of tissue-specific and stress responsive WRKYs in wild potatoes. WRKY TFs have been found to play important roles under abiotic stresses, such as drought 8 , heat 63 , wounding 50 , and biotic constraints, such as bacteria 59 and viruses 64 . Tissue-specificity of WRKY genes has also been highlighted in different crops, such as pepper 58 , cotton 65 and soybean 66 , elucidating their role in developmental and functional processes. Our study investigated for the first time the stress response and tissue-specificity of WRKY genes in two wild potato species. Six and 11 WRKY genes were identified as flower-and leaf-specific, respectively. Zhang and collaborators 21 considered that the known protein-protein interaction network can provide important clues to better understand gene expression regulation. Basing on this, we investigated the interactome of tissue-specific and stress-responsive WRKYs identified here and found potential co-expression relationships between ScWRKYs/SchWRKYs and genes of various pathways. From our analyses, interesting observations and different clues for future functional studies have emerged. For example, ScWRKY002 could be in some way involved in anthocyanin activation in flowers of S. commersonii: it interacts with anthocyanin bHLHs and the flower of this wild species strongly accumulates anthocyanins 14,67 . Previous studies reported that some WRKYs can be involved in the coordination of multiple biological processes. For example, AtWRKY33 regulates disease resistance, NaCl tolerance and thermotolerance [68][69][70] , while GhWRKY40 modulates tolerance to wounding stress and resistance to R. solanacearum 43 . This suggests that some WRKY proteins provide important nodes of crosstalk between different physiological processes. However, the putative members of WRKY family and their possible roles in signalling crosstalk are still barely known. To the authors' best knowledge, no expression data are available on ScWRKYs and SchWRKYs; by contrast, StWRKYs have previously received attention. Among them, only Shahzad et al. 71 and Yogendra et al. 72 found StWRKY010 (PGSC0003DMP400029302) and StWRKY020 (PGSC0003DMP400028763) to be active in P. infestans-potato interaction. Consistently with these data, our results indicated that the same genes increased their expression after BABA treatment, known to (2020) 10:7196 | https://doi.org/10.1038/s41598-020-63823-w www.nature.com/scientificreports www.nature.com/scientificreports/ confer protection against several biotic threats. Furthermore, we focused our attention on a group of proteins (ScWRKY016, ScWRKY023, ScWRKY045 and ScWRKY055) which were reported to be stress-responsive 21,73 .
For these genes, we tested the transcriptional activity of wild S. commersonii alleles after wounding and bacterial infection and investigated on their direct orthologs expression following various treatments. Among them, StWRKY016, StWRKY045 and StWRKY055 appeared to be required by plants to face damages by heat stress, while StWRKY023 was reported to be active under mannitol and GA3 treatments as well as drought stress 73 . Our results suggested that the wild alleles of ScWRKY023 and ScWRKY045 might represent promising candidates for multiple stress responses as they are leaf-specific and constantly expressed after wounding in S. commersonii but not in the cultivated potato. In addition, WRKY023 is also induced by bacterial infection and it is suggested to interact with both a WRAF2-like protein and with the LRR mediated immunity system 74,75 .

Conclusions
The present study identified 79 and 84 genes encoding putative WRKY TFs in S. commersonii and S. chacoense, respectively. Their protein structure and data from the comparative analyses suggested an interspecific variability of WRKY genes. Most of them were up-regulated under stress conditions and across different tissues, hinting a possible role in the cross-talk between plant and environmental cues in potato species. Taken as hole, these analyses will help to hasten the determination of the function of WRKY TFs especially in response to biotic and abiotic stresses. Candidate ScWRKY and SchWRKY genes identified here can be employed in potato breeding programs.