Genome-wide in silico identification of glutathione S-transferase (GST) gene family members in fig (Ficus carica L.) and expression characteristics during fruit color development

Glutathione S-transferase (GSTs), a large and diverse group of multi-functional enzymes (EC 2.5.1.18), are associated with cellular detoxification, various biotic and abiotic stress responses, as well as secondary metabolites transportation. Here, 53 members of the FcGST gene family were screened from the genome database of fig (Ficus carica), which were further classified into five subfamilies, and the tau and phi were the major subfamilies. These genes were unevenly distributed over all the 13 chromosomes, and 12 tandem and one segmental duplication may contribute to this family expansion. Syntenic analysis revealed that FcGST shared closer genetic evolutionary origin relationship with species from the Ficus genus of the Moraceae family, such as F. microcarpa and F. hispida. The FcGST members of the same subfamily shared similar gene structure and motif distribution. The α helices were the chief structure element in predicted secondary and tertiary structure of FcGSTs proteins. GO and KEGG indicated that FcGSTs play multiple roles in glutathione metabolism and stress reactions as well as flavonoid metabolism. Predictive promoter analysis indicated that FcGSTs gene may be responsive to light, hormone, stress stimulation, development signaling, and regulated by MYB or WRKY. RNA-seq analysis showed that several FcGSTs that mainly expressed in the female flower tissue and peel during ‘Purple-Peel’ fig fruit development. Compared with ‘Green Peel’, FcGSTF1, and FcGSTU5/6/7 exhibited high expression abundance in the mature fruit purple peel. Additionally, results of phylogenetic sequences analysis, multiple sequences alignment, and anthocyanin content together showed that the expression changes of FcGSTF1, and FcGSTU5/6/7 may play crucial roles in fruit peel color alteration during fruit ripening. Our study provides a comprehensive overview of the GST gene family in fig, thus facilitating the further clarification of the molecular function and breeding utilization.


INTRODUCTION
The glutathione S-transferase (GST ) gene family is an ancient and multigene group that is found in almost all living organisms. Two binding domains, the glutathione binding site (G-site) in the N-terminal (GST-N) and the substrate binding site (H-site) in the Cterminal, are well conserved in typical GST enzymes (Edwards & Dixon, 2005). The main function of GSTs is to catalyze the conjugation of an array of electrophilic compounds of exogenous and endogenous origins to reduced glutathione (GSH) (Cummins et al., 2011).
According to the gene structure, protein sequence similarity, gene function, and immunological characteristics within the family, the GST gene family can be further categorized into 14 subfamilies. Of these, four subfamilies (tau (U), phi (F), lambda (L), and TCHQD) are unique to the plant (Lallement et al., 2014). As an ancient protein superfamily with multipurpose roles, GST genes have diverse functions in plants (Vaish et al., 2020). Previous studies have emphasized the detoxification roles of GSTs via transformation, conjugation, and compartmentation phases, which thereby decrease toxic damage (Edwards & Dixon, 2005). They also showed the widely indispensable regulatory functions when plants encounter various stressors, such as cold (Kayum et al., 2018;Song et al., 2021a), salinity (Xu et al., 2016), heavy metals (Gao et al., 2020;Jing et al., 2021), drought (Chen et al., 2012;Xu et al., 2016), and chemical toxicity (Xu et al., 2016). Another outstanding role of several GST genes is governing the sequestration of secondary metabolic compounds, such as flavonoids, anthocyanin, proanthocyanidins (PAs), and porphyrins (Dixon, Skipsey & Edwards, 2010;Liu et al., 2019a). Other studies have indicated that GST genes play important roles in coloration. Several GSTF subfamily genes are required for anthocyanin transport and accumulation, such as TT19 (GSTF12) in Arabidopsis (Kitamura, Shikazono & Tanaka, 2003;Sun, Li & Huang, 2012), An9 in petunia (Mueller et al., 2000), CkmGST3 in cyclamen (Kitamura et al., 2012), DcGSTF2 in carnations (Sasaki et al., 2012), LcGST4 in litchi (Hu et al., 2016), FvRAP in strawberry (Luo et al., 2018), CsGSTa in tea (Liu et al., 2019a), CsGSTF1 in purple tea (Wei et al., 2019a), AcGST1 in kiwifruit (Liu et al., 2019b), MdGSTF6 in apple (Jiang et al., 2019), IbGSTF4 in sweet potato (Kou et al., 2019), PpGST1 in peach (Zhao et al., 2020), GST1 in Japanese gentian (Tasaki et al., 2020), GhGSTF12 in cotton (Shao et al., 2021), LhGST in lilies (Cao et al., 2021), and PcGST57 in pear (Li et al., 2022). Several members of the GSTU subfamily are also involved in the regulation of anthocyanin accumulation, including Bz2 in maize (Marrs et al., 1995), CsGSTb and CsGSTc in tea (Liu et al., 2019a), and MdGSTU12 in apple (Zhao et al., 2021). Although it is known that GST is involved in anthocyanin accumulation in many species, there is a lack of information about GST genes' role in fig fruit color formation. Fig (Ficus carica L.), belonging to the Moraceae family, originated in western Asia, and is one of the earliest domesticated fruit crops. Ripe fruits contain a variety of beneficial bioactive ingredients, such as dietary fiber, sugars, minerals, carotenoids, and anthocyanins (Wang et al., 2021). The common fig can be classified into green, yellow, and red cultivar categories according to the different coloration levels caused by anthocyanin content in the fruit peel. Anthocyanin content in fig contributes to its nutritional quality and attractiveness, with red fruit having a high market value. To evaluate the potential roles of GST in anthocyanin accumulation, we carried out comprehensive identification and expression characteristic analysis of the GST gene family in fig. In this study, we used the public high-quality reference genome sequence of fig (Usai et al., 2020) to identify a total of 53 GST gene members. Gene features of FcGSTs were determined based on bioinformatics analysis, and we found that FcGSTF1 and FcGSTU5/6/7 may contribute to fruit peel color. Our results may facilitate the further functional investigation of FcGSTs' role in anthocyanin accumulation in fig.

FcGST candidate gene identification
Genome sequences and annotation information of fig (Ficus carica 'Dottato') were collected from the National Center for Biotechnology Information (NCBI) (https://www.ncbi. nlm.nih.gov/genome/?term=Ficus+carica) (Usai et al., 2020). The protein sequences of the fig genome were predicted using the Batch translator CDS to protein function of TBtools (Chen et al., 2020) and by protein sequences of fig obtained from EnsemblPlants (http://ftp.ebi.ac.uk/ensemblgenomes/pub/release-54/plants/fasta/ficus_carica/pep/). The protein sequences of 53 AtGSTs in Arabidopsis were downloaded from TAIR (https: //www.arabidopsis.org/browse/genefamily/gst.jsp) (Sappl et al., 2004). The protein sequences of 79 OsGSTs in rice were obtained from the Rice Genome Annotation Project Database (RGAP) (http://rice.uga.edu/downloads_gad.shtml) according gene id (Jain, Ghanashyam & Bhattacharjee, 2010), but the protein sequences of OsGSTU3 (LOC_Os10g38501) and OsGSTU4 (LOC_Os10g38495) were not found in RGAP. Therefore, a total of 130 sequences were taken as queries in the local BLASTP search against the protein sequences of fig to retrieve all GST family candidates with E-value (≤ e −10 ) and identify ≥ 50% as thresholds. Moreover, the hidden Markov model (HMM) profiles of the GST-N (PF02798) and GST-C (PF00043) domain, obtained from the Pfam database (http://pfam.xfam.org/), were also employed as queries to search the FcGST candidates using HMMER 3.1 with a cut off E-value (≤ e −5 ). After redundant sequences were removed, the rest FcGST candidates were validated using SMART (http://smart.embl-heidelberg.de/) and NCBI conserved domain database (CDD) (https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi) to examine the presence of the conserved domains with the default parameters. The protein sequences and gene id of 53 AcGSTs, 77 OsGSTs, and 53 FcGSTs are slisted in Data S1.

Phylogenetic analysis and subfamily categorization
Multiple-sequence alignments of GST protein sequences from Arabidopsis, rice, and fig were performed using the default settings of Clustal Omega (https://www.ebi.ac.uk/ Tools/msa/clustalo/). A phylogenetic tree was generated using the default parameters of neighbor-joining method in MEGA X software with 1,000 bootstrap replicates (http://www.megasoftware.net/). According to the classification records of subfamily members in Arabidopsis and rice, these FcGST genes were further categorized into five different subfamilies. To precisely predict their molecular function, GST genes related to anthocyanin accumulation were used to construct a phylogenetic tree with 53 FcGST genes via MEGA X software according to the procedures mentioned above.
The length of cDNA and CDS were retrieved from fig genome annotation data. Physicochemical parameters (such as polypeptide length, molecular weight, and isoelectric point) of FcGST proteins were evaluated using the ExPasy website (http://web.expasy.org/ protparam/). The online soft tools of PSORT (http://psort.hgc.jp/form.html) and BUSCA (http://busca.biocomp.unibo.it/) were used to predict these FcGSTs subcellular localization.

Chromosomal locations and collinearity analysis
A genome annotation file was used to collect 53 FcGST genes' chromosomal location information. The positions of the FcGSTs were plotted to chromosomes using TBtools, and only gene matched chromosomes were shown (Chen et al., 2020).
The genome data of Malus domestica and Prunus persica were downloaded from Genome Database for Rosaceae (GDR) (https://www.rosaceae.org). The genome data of Vitis vinifera were obtained from EnsemblPlants (http://ftp.ebi.ac.uk/ensemblgenomes/pub/release-52/plants/fasta/vitis_vinifera/). The genome data of Ficus hispida and F. microcarpa were collected from the database of the National Genomics Data Center (NGDC) with BioProject number PRJCA002187. The collinear relationships among the FcGST genes in fig and A. thaliana, rice, apple, peach, and grape were performed using the Multiple Collinearity Scan toolkit (MCScanX) with default parameters (Wang et al., 2012). The collinear relationships among F. hispida, F. microcarpa, and fig were also analyzed using the default parameters of MCScanX. Two or more genes located within the present 100 kb region on the same chromosome were defined as a tandem repeat event, while those located beyond the 100 kb region were considered a segmental duplication event.

Functional analysis
GO annotation analysis of all fig protein sequences was performed using BLASTP with Swiss-UniProt database with E-value (≤ e −5 ). KEGG annotation analysis was performed using the KofamKOALA online program with default parameters (https://www.genome.jp/tools/kofamkoala/). TBtools was used to obtain the GO and KEGG enrichment terms with corrected p-value (≤ 0.05) and to plot the enrichment results (Chen et al., 2020). A protein-protein interaction (PPI) network of 53 FcGST proteins was generated using STRING database V11.5 (https://cn.string-db.org/) (Szklarczyk et al., 2021) and visualized in Cytoscape 3.6.1 software program (Shannon et al., 2003).

Cis-acting element analysis of FcGST genes
For cis-element analysis, the putative promoter sequences (2,000 bp sequence upstream of the translation site) were acquired from all FcGST genomic DNA sequences using TBtools (Chen et al., 2020). The putative cis-acting elements were predicted using PlantCARE program with default parameters (http://bioinformatics.psb.ugent.be/webtools/plantcare/ html/).

RNA-seq data analysis
Three published raw RNA-seq data downloaded from the NCBI Sequence Read Archive (SRA) database (https://www.ncbi.nlm.nih.gov/sra/) were mined to investigate the expression characteristics of FcGST genes during fig fruit development. In common purple-peel fig cultivar 'Zibao', the development process is divided into six stages. Stages 1 and 2 belong to a rapid growth stage. In stages 3 and 4, fruit size and hardness remain almost unchanged. Stages 5 and 6 belong to the fruit mature stage, where stage 5 corresponds to commercial ripeness (Zhai et al., 2021). The expression patterns of FcGSTs were analyzed at six fruit development stages of the internal female flower (F1-F6) and syconia peel (P1-P6) (accession number SRP315833) (Zhai et al., 2021). RNA-seq data of 'Green Peel' and its bud mutation cv. 'Purple Peel' at young and mature fruit development stages were used to estimate the expression level of FcGSTs in peel (accession number SRP114533) (Wang et al., 2017). Unlike typical climacteric fruit, such as banana, figs cannot ripen by themselves or with exogenous treatment once harvested before commercial ripeness stage. Ethephon treatment has been reported to stimulate fig fruit growth and shorten the ripening period with no negative effect on fruit quality. Ethephon solution was injected through the fruit ostiole in order to determine the role of ethylene in fig-fruit ripening using RNA sequencing (Cui et al., 2020). Transcriptome data of fig fruit after ethephon application at two, four, and six days were used to investigate the FcGST gene expression profile in ethylene-regulated pigmentation changes of fig female flowers (accession number SRP249501) (Cui et al., 2020).
The quality of raw reads were checked using FastQC (version 0.11.8, https://www. bioinformatics.babraham.ac.uk/projects/fastqc/). After adapter removal by Trimmomatic (Bolger, Lohse & Usadel, 2014), the clean reads were mapped to fig genome data by HISAT2 with default parameters (Kim et al., 2019). The transcript abundance of each sample was calculated by StringTie Quantify (Pertea et al., 2015). The transcripts per kilobase of exon model per million mapped reads (TPM) was applied to describe the gene transcript level. The TPM value was transformed into log 2 (TPM + 1) in the presented heatmaps. DESeq2 R package was used to identify the differential expression genes (DEGs, |log 2 (fold change)|

Measurement of anthocyanin content
The total anthocyanin content of fig fruit peel was measured according to Li, An & Wang (2020). After approximately 0.5 g of fig peel was incubated in 5 ml extraction solution (1% (v/v) HCl-methanol) for 24 h at room temperature under shading, the supernatant was centrifuged at 8,000 g for 15 min, and the absorbance at 530, 620, and 650 nm was measured. The anthocyanin content was calculated according to the following formula: × 10 6 , V represents the volume of extraction solution; M represents the weight of each sample, and the absorbency index of total anthocyanin is 4.62 × 10 4 . Three independent replicates were used for mean value calculation.

Quantitative real-time PCR (qRT-PCR) analysis
As for qRT-PCR assay, total RNA extraction from the fig fruit peel of three biological replicates per sample was carried out using the CTAB method. We used the One-Step gDNA Removal and cDNA Synthesis Supermix kit (Transgen Biotech, Peking, China) to reverse transcribe the total RNA. We added 20-µL of reaction mixture to each well, including 1 µL cDNA, 0.5 µL each of 5 and 3 specific primer, 10 µL chamQ SYBR qPCR Master Mix (Vazyme, Shanghai, China), and 8 µL DEPC water. The ABI 7300 Real-Time PCR system was employed to perform the qRT-PCR amplification program with the following settings: 95 • C for 5 min, 40 cycles of: 95 • C for 5 s, 60 • C for 35 s. The relative expression of target genes was calculated with the 2 − CT method (Livak & Schmitten, 2001) and normalized by FcActin (Vangelisti et al., 2019). Each sample was quantified with three replicates. All qRT-PCR primers are listed in Data S2.

Data analysis
Statistical analysis was performed with SPSS 26.0 software. Significance (p < 0.05) was evaluated by analysis of variance (ANOVA) and Duncan's multiple range test. Graphs were generated in Excel 2016. To evaluate the evolutionary relationships among the AtGSTs, OsGSTs, and FcGSTs, a phylogenetic tree was constructed using the NJ method of MEGA X. These proteins were classified into eight distinct subfamilies. Tau, lambda, zeta, TCHQD, DHAR, EF1Bγ , theta, and phi ( Fig. 1). Among these subfamilies, the tau and phi groups were found to be the most abundant in Arabidopsis, rice, apple, peach, melon, pepper, and The gene length of the FcGST genes varied from 291 bp (FcGSTU19) to 6,368 bp (FcGSTT1) ( Table 1). The lengths of the coding DNA sequence (CDS) of 34 FcGSTs (64.15%) were between 603 (FcGSTU33) and 696 bp (FcGSTU21). The polypeptide of most FcGSTs were similar, with approximately 43 members (81.12%) ranging from 201 aa to 269 aa. The average molecular weight (MW) of FcGST protein was 25.73 kDa. The change of putative theoretical isoelectric point (pI) varied from 4.71 (FcGSTU34) to 10.38 (FcGSTU19), and 40 members (75.47%) were acidic proteins. Two online soft tools, PSORT and BUSCA, were used to predict FcGST subcellular localization. Most of the FcGST proteins were predicted to be localized in the cytoplasm (43), followed by microbody (26), mitochondria (15), chloroplast (nine), nucleus (seven), plasma membrane (four), and vacuole (one). Fourteen members (26.42%) were densely mapped on chromosome 9. Chromosome 6 contained eight FcGST genes, followed by chromosome 10, 11, and 3 with six, five, and four members, respectively. Chromosomes 2, 5, and 8 had three genes each, while chromosome 7 and 13 had two genes each. Chromosomes 1, 4, and 12 contained only one FcGST gene each. We also found that members of the same subfamily were mostly grouped together. Thirty-nine GSTU genes were present in eight clusters, including three two-gene clusters, two three-gene clusters and one five-, eight-and 14-gene cluster. Four GSTF genes were located in a single one-gene cluster. Gene family duplication investigation showed that the fig genome contained a total of 13 duplicated FcGST gene pairs ( Fig. 2A). A maximum of five duplicated GST genes were located in chromosome 6, followed by three in chromosome 9 and chromosome 10, and one in chromosome 11. Out of the 13 gene pairs, 12 pairs were found as tandemly duplicated, and one pair appeared segmentally duplicated. Of the tandem duplications, ten pairs were in the tau subfamily and two pairs were in the phi subfamily, indicating that the tandem events contributed to these two subfamilies' expansion. The average sequence similarity of 13 gene pairs was 66.27% (Fig. 2B), and gene pairs of the phi subfamily (FcGSTF3/4, and FcGSTF4/5) showed high protein sequence similarity, implying that the potential roles of these gene pairs might be similar. A synteny map was constructed to further estimate the homologous evolutionary relationship among GSTs of fig, the model species (Arabidopsis and rice), as well as three horticulture crops (apple, peach, and grape) (

Gene structure and conserved motif analysis of FcGSTs
As depicted in Fig. 4B, FcGST genes contained one (FcGSTU19 and FcGSTU39) to 11 (FcGSTT1) exons and showed similar exon-intron organization within the same subfamily, which was consistent with the phylogenetic analysis and data from other species (Figs. 4A To search the conserved amino acid motifs among FcGST proteins, the MEME online tool was utilized. Ten distinct motifs were identified in 53 FcGST proteins, ranging from 11 (motif 8) to 61 (motif 10) amino acids (Data S6). Except for FcGSTU10/34/31, FcEF1Bγ 1, and FcGSTF5, the remaining members all contained motif 3 (Fig. 4C). The rest of the motifs exhibited similar type and order within the same subfamily. For example, all six FcGSTF members contained motifs 3, 9, and 10, and motifs 1 and 3 were only present in FcGSTLs without the other motif. In the N-terminal of most of the FcGSTU members, the motif shared common conserved configuration. In the C-terminal, these genes could be separated into two groups according to motif 5, and motif 2 and 7, which was consistent with phylogenetic clade separation in the tau subfamily. We also observed that motifs 2, 4, 5, 6, 7, and 8 were specific to the tau subfamily.

Protein structure prediction of GST families in fig
The secondary structure of 53 FcGST proteins were predicted using NPS@: SOPMA. The results showed that the presence of α-helix occupied a large position (51.36%), followed by random coils (29.86%), extended β strands (13.49%), and β turn (5.30%) (Data S7). We observed that the members within the same FcGST subfamily possessed similar protein sequences and showed similar secondary structure. Compared with phi family members, more α-helix (54.00%) and fewer random coils (28.13%) and β strands (12.92%) presented in the tau family composed the secondary structure with the exception of FcGSTU19/24/39 (Data S7). Members of lambda showed fewer β strands (10.45%) and β turns (3.78%) than FcGSTF proteins.
Tertiary structure models of five FcGST proteins were predicted by using Swiss-Model (Fig. 5). Each FcGST protein shared over 45% sequence identity with its template (Data S8). VERIFY3D showed that over 86% of the residues of these proteins had an average 3D-1D score ≥ 0.2 (Data S8). The ERRAT score showed an overall quality factor of over 92 (Data S8). Ramachandran plot results showed that the residues in the most favored regions of five proteins were more than 90% (Data S9), indicating that high quality models were built. Corresponding with secondary structure prediction, the α-helix was the major protein structure in 3-D models of all these five FcGSTs, followed by the random coil, β sheet, and turn (Fig. 5). The α-helix mainly presented in the C-terminal domain, which facilitated formation of substrate binding site, whereas the N-terminal domain composed of α-helix and β sheets (Fig. 5;Liu et al., 2013;Song et al., 2021a). The GST domains of three GSTU proteins shared the same conformation of the tertiary structure (Fig. 5).

Cis-acting elements in promoter regions of FcGST genes
To understand the transcriptional characteristics of FcGST genes, the putative cis-elements of 53 FcGST genes promoter regions located 2,000 bp from the upstream of the translation start site (ATG) were retrieved and scanned through the PlantCARE database. Forty-four kinds of cis-acting elements were predicted in the promoter region of 53 FcGST genes (Fig. 7). These elements were classified into five groups: light, hormone, stress, development, and transcription factor (Fig. 7A). Detailed information about the cis-acting elements is provided in Data S10. Approximately 31.42% elements responded to the light (Fig. 7B), and these elements were detected in the promoter region of all FcGST genes, ranging from four to 19 (Fig. 7A), suggesting that light stimulation is an important factor in the regulation of FcGST gene expression. The hormone responsive element was the second largest group in the subfamilies of tau (21.41%), phi (29.46%), and lambda (20.17%), but not the EF1Bγ and theta subfamilies (Fig. 7B). Among these hormone-related elements, a total of 144 MeJA (methyl jasmonate) response elements (CGTCA-motif and TGACG-motif) were detected in 42 FcGSTs, followed by 122 abscisic acid response elements (ABRE) in 43 FcGSTs, 54 gibberellin response elements (GARE-motif, P-box, and TATC-box) in 36 FcGSTs, 48 auxin response elements (AuxRR-core and TGA-element) in 26 FcGSTs, and 28 salicylic acid responsive (TCA-element) elements in 20 FcGSTs (Fig. 7A). The upstream region of all 53 FcGSTs possessed stress-related elements, including anaerobic induction (ARE, and GC-motif), low-temperature response elements (LTR), drought-response elements (MBS, and MYC), and defense and stress (TC-rich repeats) (Fig. 7A). By contrast, only 3.11% cis-acting elements related to plant development were identified, including CAT-box (15), circadian (8), GCN4_motif (7), HD-Zip 1 (5), MBSI (4), MSA-like (2), O2-site (13), and RY-element (2) (Figs. 7A and 7B). Additionally, several transcription factorbinding elements were also detected in FcGST promoters, such as MYB, Myb-binding site, MYB-like sequence, MYB recognition site, and W box (WRKY-binding motif), which may be involved in the direct regulation of these FcGST transcriptions (Fig. 7A).

Expression pattern of FcGSTs during fig fruit ripening
To comprehensively illustrate FcGST gene expression characteristics, we investigated these gene expression pattern on three published RNA-seq data, including six stages of female flower and fruit peel during 'Zibao' fig fruit development (Fig. 8;Zhai et al., 2021), the fruit peel of 'Green Peel' and 'Purple Peel' at young and mature stages (Data S11A; Wang et al., 2017), and female flowers after ethephon application at two, four, and six days (Data S11B; Cui et al., 2020). During fig fruit development, most FcGSTU genes showed slight or undetectable expression (Fig. 8, Data S11A), even under ethephon treatment (Data S11B). However, FcGSTU30/38/40 expressed highly in both female flower and fruit peel, and the expression abundance significantly increased in fruit maturation (Fig. 8, Data S11B). However, the difference was not observed in mature fruit peel between 'Green Peel' and 'Purple-Peel' (Data S11A). FcGSTU5/6/7 showed similar expression trends in fruit peel, which slightly declined to the bottom at P4 but rapidly increased in the mature stage (P5 and P6), and was especially highly expressed in the mature fruit peel of 'Purple-Peel' (Fig. 8, Data S11A). The specific expression changes of FcGSTF1 were noticeable. FcGSTF1 was significantly up-regulated until F4 and then suddenly decreased in F5 and F6, whereas its expression was continually upregulated in fruit peel (Fig. 8). Compared with 'Green Peel', FcGSTF1 exhibited high expression abundance in the mature fruit peel of purple peel (Data S11A), which completely matched with the fruit peel phenotype (Wang et al., 2017). Additionally, we observed that several genes also increased in response to ethylene induction, such as FcGSTU5/6/7/38/40 after two d of ethephon treatment, and FcGSTU30 and FcGSTF1 at six d. Previous study results showed that ethephon did not accelerate anthocyanin accumulation but instead inhibited pigmentation in fig flower (Cui et al., 2020), which was highly consistent with the expression changes of anthocyanin biosynthetic genes (F3H, DFR, and LDOX ) (Cui et al., 2020), FcGSTF1, and FcGSTU5/6/7 (Data S11B). These results showed that FcGSTF1 and FcGSTU5/6/7 expression was related to the development of anthocyanins, and we conjecture that these genes may play additional roles in the coloration of fig fruit.

Anthocyanin accumulation-related GSTs gene features and the expression assay in fig fruit peel
Furthermore, we investigated the gene features and expression patterns in the fruit peel of two common fig cultivars during the fruit ripening. Multiple sequence alignment and phylogenetic analyses with other reported anthocyanin-related GST genes were performed. The results showed that FcGSTF1, CsGSTa, AcGST1, CkmGST3, GhGSTF12, VvGSTF4, FvRAP, PpGST1, PcGST57, MdGSTF6, LcGST4, IbGSTF4, PhAN9, AtGSTF12, and LhGST were clustered together (Fig. 9A). FcGSTU5/6/7 were found to be closely associated with anthocyanin-related tau family members, such as MdGSTU12, CsGSTc, and CsGSTb . Among these genes, the expression of FcGSTF1 was significantly higher than FcGSTU5/6/7, even at the fruit color turning stage (60 days after fruit setting) (Fig. 9F).

FcGSTF1 and FcGSTU5/6/7 are associated with anthocyanin accumulation in fig peel
Anthocyanin is not only an important and beneficial pigments, it is also a key appearance trait of fruit coloration. Recently, several studies have characterized the molecular function of a series of anthocyanin biosynthetic genes during fig fruit coloration, including FcPAL, FcC4H, Fc4CL, FcCHS1, FcCHS2, FcCHI, FcF3H, FcF3 H, FcDFR, FcANS, and FcUFGT (Wang et al., 2017;Lama et al., 2020;Li, An & Wang, 2020). However, the function of FcGST is still unknown in fig, while its homologous genes have been verified as the anthocyanin transporter in various plant species (Kitamura, Shikazono & Tanaka, 2003;Sun, Li & Huang, 2012;Luo et al., 2018;Jiang et al., 2019). In our study, GO analysis showed that several FcGST family genes were enriched in the flavonol biosynthetic process (GO:0051555) and anthocyanin-containing compound metabolic process (GO:0046283) (Fig. 6A). During fig fruit development, anthocyanin accumulation initiation was earlier in female flower tissues than in fruit peel (Wang et al., 2021;Cui et al., 2020). The particular expression pattern of FcGSTF1 was highly matched with the obvious space-time differences of the coloration process between female flower and fruit peel ( Fig. 8; Data S11B). Additionally, FcGSTF1 showed high sequence similarity with the GSTs' known function in anthocyanin transport (Figs. 9B and 9D) (Kitamura, Shikazono & Tanaka, 2003;Luo et al., 2018;Jiang et al., 2019). We also found that FcGSTU5/6/7 were significantly induced in the mature stages (stage 5 and 6) of 'Zibao' fruit peel ( Fig. 8; Data S11A), and showed closer relationships as well as similar protein sequences with anthocyanin-related tau family members (Figs. 9A, 9C and 9D) (Liu et al., 2019a;Liu et al., 2019b;Zhao et al., 2021). Moreover, high consistency between anthocyanin accumulation trends and expression level changes of FcGSTF1 and FcGSTU5/6/7 were also detected (Fig. 9). The results reveal that these genes may play important roles in fig fruit peel coloration, although the specific function of these genes in anthocyanin transport requires additional research. In high plants, MYB is the crucial transcriptional factor (TF) of secondary metabolism, and is not only associated with expression level regulation of anthocyanin biosynthesis with other anthocyanin transport-related GSTFs, and FcGSTU5/6/7 with ZmBZ2, CsGSTb, MdGSTU12, and CsGSTc. The different shades of red represent these sequences' similarity degree. (E) The anthocyanin content in the fruit peel of 'Orphan' and 'Bojihong' 30, 60, and 70 days after fruit setting. (F) Expression profile of FcGSTF1 and FcGSTU5/6/7 were detected in the fruit peel of 'Orphan' and 'Bojihong' 30, 60, and 70 days after fruit setting. Gene relative expression was normalized to the expression level of 'Orphan' fruit peel at 30 days after fruit setting. Different lowercase letters indicate significant difference (p < 0.05).

Characteristics of FcGST gene family members in fig genome
GST is a large, ubiquitous, and ancient protein superfamily, that has been identified in various plant species (Sappl et al., 2004;Jain, Ghanashyam & Bhattacharjee, 2010;Dong et al., 2016). In this study, we identified a total of 53 FcGST gene members in fig, and classified them into five distinct subfamilies (Table 1). Among the different subfamilies, tau was the largest, followed by phi (Fig. 1), which supported their dominant distribution in other species, such as Arabidopsis (Sappl et al., 2004), rice (Jain, Ghanashyam & Bhattacharjee, 2010), and apple (Jiang et al., 2019). GST proteins are localized in various subcellular compartments, and they mainly show cytosolic localization (Raza, 2011). Likewise, most FcGST members were also predicted to profusely present in cytoplasm (Table 1). In plants, most cytosolic GSTs typically function as dimers with molecular weights ranging from 23 to 29 kDa (Frova, 2006), which is similar to GrGSTs (26.48 kDa) in G. raimondii, GaGSTs (26.64 kDa) in G.arboreum (Dong et al., 2016), CmGST (23.51 kDa) in Hami melon (Song et al., 2021a), and FcGSTs (26.25 kDa) in fig (Table 1). Tandem, segmental, and whole-genome duplication are common genetic events that are responsible for gene family expansion and diversification in plants (Flagel & Wendel, 2009). Previous studies suggested that segmental gene duplication mainly occurred for the amplification of tau and phi subfamily GST genes in pear (Wang et al., 2018), tomato (Islam et al., 2017), pepper (Islam et al., 2019), chickpea (Ghangal et al., 2020, and M. truncatula (Hasan et al., 2021). However, 20 of the 48 FcGSTU and FcGSTF genes had tandem duplicate events (Fig. 2) subfamilies. It appears that tandem repeat events mainly contributed towards the evolution of tau and phi subfamily GST genes in fig. The conserved exon-intron organization could reveal information about the evolutionary relationship (Sánchez et al., 2003). In the same subfamily, a specific exonintron arrangement pattern and several length-conserved exons were identified in fig and other plant species (Sappl et al., 2004;Jain, Ghanashyam & Bhattacharjee, 2010;Islam et al., 2017;Jiang et al., 2019;Wang et al., 2020;Hasan et al., 2021) (Fig. 4B; Data S5). Most FcGST members also shared a similar motif organization in the same subfamily (Fig.  4C). These results reveal that FcGST family members are evolutionarily conserved across distinct phylogenetic groups. Moreover, the members of eight different GST subfamilies showed high similarities between fig, Arabidopsis, and rice ( Fig. 1), indicating the ancient evolution of these subfamilies before the split of monocots and dicots (Islam et al., 2017;Islam et al., 2019;Wang et al., 2020).

FcGSTs may play pleiotropic roles in fig
GST family genes show versatile functions in diverse plants (Vaish et al., 2020). GO, KEGG enrichment, and cis-acting element analyses reflected their multipurpose roles in fig (Figs.  6 and 7). A considerable amount of information has been gathered on various abiotic stress management roles of GST family genes in plants (Vaish et al., 2020). In general, the complex stress regulation roles of the GST gene family are highly correlated with reactive oxygen species (ROS) scavenging (Estévez & Hernández, 2020). Here, the molecular function of glutathione peroxidase activity (GO:0004602), biological process of ROS stress response (GO:0006979, GO:0000302, GO:0051775), toxic substance (GO:0009636), freezing (GO:0050826) and salt (GO:0009651), as well as the glutathione metabolism pathway (ko00480), were enriched (Figs. 6A and 6B). Furthermore, we also found that at least one stress responsive element was presented in the FcGST promoter region which might be associated with anaerobic stress, drought, and cold (Fig. 7). Therefore, further studies are needed to verify the positive roles of FcGST genes in the improvement of abiotic stress tolerance in fig.
GSTs are known to respond to multiple phytohormones, including auxin, cytokinin, MeJA, ABA, salicylic acid, ethylene, and other hormones (Hasanuzzaman et al., 2017). GO and promoter prediction analysis suggest that FcGSTs may function in various hormoneregulated fig tree growth and development (Figs. 6A and 7).
In plants, tau subfamily members might be involved in sugar signaling. Five tau genes in Arabidopsis and 26 tau GST members in sorghum were up-regulated under sucrose treatment (Chi et al., 2011). In pear, the expression of PpGST2 in the tau subfamily was significantly increased by glucose during fruit development (Shi et al., 2014). In this study, the expression of

CONCLUSION
In this study, a total of 53 GST genes identified in the fig genome database were divided into five distinct subfamilies. The members of FcGST in the same phylogenetic subfamily showed evolutionary conservation with other plant species in the GST family. The FcGST family genes may play multiple functions in fig tree growth and development. Furthermore, our results suggest that FcGSTF1 and FcGSTU5/6/7 might play roles in regulating fruit peel coloration, but the specific molecular function and regulation mechanism still need to be verified by further study.

ADDITIONAL INFORMATION AND DECLARATIONS Funding
This work was financially supported by the Natural Science Foundation of Anhui Province (No. 1908085QC108), and the University Natural Science Research Project of Anhui Province (No. KJ2020A0043 andNo. KJ2021B09). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.