Comparative Transcriptome Profiles of Human HaCaT Cells in Response to Gynostemma pentaphyllum Extracts Obtained Using Three Independent Methods by RNA Sequencing

Gynostemma pentaphyllum (GP) is widely used in herbal medicine. In this study, we developed a method for the large-scale production of GP cells using plant tissue culture techniques combined with bioreactors. Six metabolites (uridine, adenosine, guanosine, tyrosine, phenylalanine, and tryptophan) were identified in GP extracts. Transcriptome analyses of HaCaT cells treated with GP extracts using three independent methods were conducted. Most differentially expressed genes (DEGs) from the GP-all condition (combination of three GP extracts) showed similar gene expression on treatment with the three individual GP extracts. The most significantly upregulated gene was LTBP1. Additionally, 125 and 51 genes were upregulated and downregulated, respectively, in response to the GP extracts. The upregulated genes were associated with the response to growth factors and heart development. Some of these genes encode components of elastic fibers and the extracellular matrix and are associated with many cancers. Genes related to folate biosynthesis and vitamin D metabolism were also upregulated. In contrast, many downregulated genes were associated with cell adhesion. Moreover, many DEGs were targeted to the synaptic and neuronal projections. Our study has revealed the functional mechanisms of GP extracts’ anti-aging and photoprotective effects on the skin using RNA sequencing.


Introduction
Plants have been regarded as a useful source of natural medicines for curing diverse human diseases since ancient times [1]. Plants are rich in phytochemicals that are useful as natural medicines and provide health benefits. Nowadays, natural medicines are more popular because of their beneficial effects and lower prices, compared to high-cost synthetic drug agents, and are widely used for the prevention and treatment of many diseases [2]. Hundreds of plant species have been widely used as natural herbal materials in many industries, such as pharmaceuticals, cosmetics, and nutraceuticals [3].
Gynostemma pentaphyllum (Thunb.) Makino, known as Jiaogulan in China, is a dicotyledonous and perennial plant belonging to the Cucurbitaceae family, which includes cucumber, pumpkin, and melon [4]. G. pentaphyllum (GP) has a dioecious reproductive system, with distinct individual organisms that produce male or female gametes. GP is an herbaceous vine plant distributed in the United States and Asian countries, such as China, Korea, Japan, Malaysia, and India, and has been traditionally consumed as a food and tea in many countries [5].

Treatment of GP Cells with Radiofrequency and Extraction of GP Cells
Three different GP extracts were prepared in this study. The GP extract from the plants, referred to as GP-P, was prepared as follows. GP plants were dried at 60 • C for two days using an agricultural dryer (Kiturami, KED-S12D1, Seoul, Republic of Korea). Subsequently, the dried GP plants were heat extracted with distilled water at 121 • C for 15 min. The extracts were stored at 4 • C. To prepare GP calli (GP-C), the plant cells (callus) from the tissue culture were further cultured for six days. To prepare the GP callus treated with radiofrequency (GP-CR), the plant cells were subjected to radiofrequency at 360 kHz for 20 min per day, and the treatment was repeated for three days. Thereafter, the cells treated with or without radiofrequency were harvested using a non-woven fabric filter. The harvested plant cells were subjected to heat extraction using distilled water at 121 • C for 15 min. For all three extract preparations, 2 g/L of the sample was used for extraction. Solids were removed by filtration through a mesh after heat extraction.

Cultivation and Differentiation of the Human Cells
Human HaCaT keratinocytes (Cat. No. CRL-2404), obtained from the American Type Culture Collection (ATCC) (Manassas, VA, USA), were cultured in Dulbecco's Modified Eagle's Medium (DMEM; Welgene, Gyeongsan-si, Republic of Korea) supplemented with 10% heat-inactivated fetal bovine serum (FBS) (Thermo Fisher Scientific, Waltham, MA, USA) and 100 U/mL penicillin/streptomycin mixture (Thermo Fisher Scientific) at 37 • C in a humidified atmosphere containing 5% CO 2 . The cells were passaged at 80-90% confluence. mycoplasma contamination is very common in cell cultures. To confirm the presence of mycoplasma in the cell cultures, we performed PCR using specific primers (MGSO and GPO-3) for the Mycoplasma species amplifying a PCR product with 270 bp in length as described in the previous study [30]. As a result, we confirmed that our HaCaT cell cultures were free from the mycoplasma contamination. We worked with HaCaT cells carefully under aseptic conditions to prevent contamination. For the treatment of HaCaT cells, 0.1% of the GP extract from either of the three different extraction conditions (GP-P, GP-C, and GP-CR) was added to the cell culture medium and incubated for 24 h; the control samples were treated with an equal volume of distilled water.

Total RNA Isolation, Library Preparation, and RNA Sequencing
For the transcriptome analyses of HaCaT cells in response to the treatments (GP-P, GP-C, and GP-CR), the cells were seeded at a density of 1 × 10 6 cells per well in a six-well plate and incubated for 24 h. They were then treated with 0.1% of either of the three extracts (GP-P, GP-C, and GP-CR) or distilled water (control) for 24 h. Three biological replicates were harvested for each condition. Total RNA was isolated from the cells using TRIzol reagent (Invitrogen, Waltham, MA, USA) according to the manufacturer's protocol, and the quality of the extracted total RNA was measured using Bioanalyzer 2100 (Agilent, Santa Clara, CA, USA). RNA samples were quantified using a Thermo Scientific Multiskan GO microplate spectrophotometer (Fisher Scientific Ltd., Vantaa, Finland), and purity was assessed using the ratio of absorbances at 260 nm and 280 nm (A ratio between 1.8 and 2.0 was considered acceptable). Total RNA with an RNA integrity number (RIN) value ≥ 7 was used for library preparation. RNA sequencing libraries were prepared using the TruSeq Stranded mRNA LT Sample Prep Kit (Illumina, San Diego, CA, USA) according to the manufacturer's instructions. Twelve libraries with their respective indices were paired-end sequenced using the HiSeq X system (Illumina).

Mapping and Identification of DEGs
The raw data were deposited in the NCBI SRA database under the PRJNA903215 accession number. Raw sequence reads were mapped to the human reference genome (GRCh38) (https://www.ncbi.nlm.nih.gov/genome/guide/human/ (accessed on 10 October 2022)) using BBMap ver. 39.00 with the default parameters. HTSeq ver. 2.0.1 was used to quantify the gene expression for each gene using the GRCh38 general transfer format (GTF) (Homo_sapiens.GRCh38.108.chr.gtf) derived from Ensemble [31]. The number of mapped reads in each transcript was used to identify DEGs using DESeq2 ver. 1.38.3 in DEBrowser v1.24.1 [32]. We compared the GP extract-treated HaCaT cells with those treated with distilled water (control). The number of reads was normalized using the MRE normalization method without any correction. By applying a fold change (FC) of more than four times and p-values less than 0.01, we identified DEGs for each comparison. For the GP-all comparison, all nine datasets from the GP extract-treated samples were compared with the three control datasets.

Gene Set Enrichment Analyses
The DEGs identified from GP-all were divided into upregulated and downregulated genes. The gene symbols in each group were used for gene set enrichment analyses using the WebGestalt program [33]. WebGestalt uses an overrepresentation analysis (ORA) against the Gene Ontology (GO) functional database (Release on 1 January 2019), which is further divided into three categories: biological processes, cellular components, and molecular functions. First, we determined the number of genes mapped to GO Slim. Genes annotated to the selected functional categories were used for the enrichment analysis of the human genome. The minimum and maximum number of IDs in the category were set to 5 and 2000, respectively. Bonferroni and the top 100 were used for the FDR method and significance level, respectively. Further, only the identified functional gene sets with a p < 0.05 were selected.

Production of Plant Cells Using Plant Tissue Culture Technology
The prepared GP leaves were sterilized and used for the callus induction ( Figure 1). The induced GP callus was used for suspension cell cultivation in a liquid medium. Finally, GP cells were mass-cultivated using bioreactors.
Three different GP extracts, such as the GP extract from the plants (GP-P), the extract from the GP plant cells (callus) (GP-C), and the GP plant cells treated with radiofrequency (GP-CR) were prepared in this study. Next, we conducted an HPLC analysis of three different extracts ( Figure S1). The chromatograms showed the representative six highest peaks in all three GP extracts ( Figure S1A). Interestingly, the six peaks were much higher in GP-C and GP-CR than in GP-P. We examined the six peaks in each GP extract using MS/MS analysis ( Figure S1B)and identified uridine (peak 1), adenosine (peak 2), guanosine (peak 3), tyrosine (peak 4), phenylalanine (peak 5), and tryptophan (peak 6) by comparison to the MS/MS of the authentic compounds ( Figure S1C).

Production of Plant Cells Using Plant Tissue Culture Technology
The prepared GP leaves were sterilized and used for the callus induction ( Figure 1). The induced GP callus was used for suspension cell cultivation in a liquid medium. Finally, GP cells were mass-cultivated using bioreactors. Three different GP extracts, such as the GP extract from the plants (GP-P), the extract from the GP plant cells (callus) (GP-C), and the GP plant cells treated with radiofrequency (GP-CR) were prepared in this study. Next, we conducted an HPLC analysis of three different extracts ( Figure S1). The chromatograms showed the representative six highest peaks in all three GP extracts ( Figure S1A). Interestingly, the six peaks were much higher in GP-C and GP-CR than in GP-P. We examined the six peaks in each GP extract using MS/MS analysis ( Figure S1B)and identified uridine (peak 1), adenosine (peak 2), guanosine (peak 3), tyrosine (peak 4), phenylalanine (peak 5), and tryptophan (peak 6) by comparison to the MS/MS of the authentic compounds ( Figure S1C).
Total RNA samples extracted from HaCaT cells treated with either extract or distilled water (control) were used to prepare 12 libraries for RNA sequencing ( Table 1). The read numbers ranged from 49,818,356 (GP-C3) to 32,703,588 reads (GP-CR3) ( Table 2). The GC content of the raw data was higher than 50% in all 12 libraries. GP-P1 had the highest GC content (51.19%), and GP-CR1 had the lowest (50.48%). The quality scores of 20 (Q20) and 30 (Q30) in all libraries were higher than 98% and 94%, respectively, indicating that the sequenced raw data was high quality. GP-P GP extract from plants GP-P3 3 Total RNA samples extracted from HaCaT cells treated with either extract or distilled water (control) were used to prepare 12 libraries for RNA sequencing ( Table 1). The read numbers ranged from 49,818,356 (GP-C3) to 32,703,588 reads (GP-CR3) ( Table 2). The GC content of the raw data was higher than 50% in all 12 libraries. GP-P1 had the highest GC content (51.19%), and GP-CR1 had the lowest (50.48%). The quality scores of 20 (Q20) and 30 (Q30) in all libraries were higher than 98% and 94%, respectively, indicating that the sequenced raw data was high quality.

Mapping and Identification of DEGs
Raw data were mapped onto the human reference genome (GRCh38) using the BBMap program with default parameters. The proportion of mapped reads in each library ranged from 90.7% to 89.9% ( Figure 2). Approximately 10% of the sequenced reads in each library were not associated with human transcripts.

Mapping and Identification of DEGs
Raw data were mapped onto the human reference genome (GRCh38) using the BBMap program with default parameters. The proportion of mapped reads in each library ranged from 90.7% to 89.9% ( Figure 2). Approximately 10% of the sequenced reads in each library were not associated with human transcripts. We identified DEGs by comparing the GP extract-treated and control samples. The numbers of DEGs identified under treatments with the three extracts were very similar ( Figure 3). Among the extracts, GP-C (93 DEGs) showed the highest number of DEGs, followed by GP-CR (90 DEGs) and GP-P (87 DEGs). In addition, we combined all GP extract-treated samples, referred to as GP-all, and compared them to the control conditions, resulting in the identification of 176 DEGs. In both GP-P and GP-all conditions, the number of upregulated genes was higher than that of downregulated genes. The number of upregulated genes (125 genes) was almost three times that of the We identified DEGs by comparing the GP extract-treated and control samples. The numbers of DEGs identified under treatments with the three extracts were very similar ( Figure 3). Among the extracts, GP-C (93 DEGs) showed the highest number of DEGs, followed by GP-CR (90 DEGs) and GP-P (87 DEGs). In addition, we combined all GP extract-treated samples, referred to as GP-all, and compared them to the control conditions, resulting in the identification of 176 DEGs. In both GP-P and GP-all conditions, the number of upregulated genes was higher than that of downregulated genes. The number of upregulated genes (125 genes) was almost three times that of the downregulated (51 genes) in GP-all. In contrast, the number of downregulated genes was much higher than that of the upregulated genes in the GP-C and GP-CR conditions.
Life 2023, 13, x FOR PEER REVIEW 7 of 19 downregulated (51 genes) in GP-all. In contrast, the number of downregulated genes was much higher than that of the upregulated genes in the GP-C and GP-CR conditions.

DEGs in Response to GP Extracts
We selected the top ten significant DEGs in each condition based on the Manhattan distance ( Figure 4). Among the 87 DEGs in GP-P (Table S1), the top ten DEGs comprised five significantly downregulated genes (CARD17, LCN12, FUOM, LOC105374145, and LOC107984639) and five significantly upregulated genes (KCND3, LOC101928092, NLRP6, S100A8, and ELDR) ( Figure 4A). In addition, we identified five long intergenic non-coding (LINC) RNAs LINC00539, LINC01103, LINC01765, LINC01792, and LINC02355. downregulated (51 genes) in GP-all. In contrast, the number of downregulated genes was much higher than that of the upregulated genes in the GP-C and GP-CR conditions.
We compared the number of DEGs identified under the four conditions. In total, 168 DEGs were highly upregulated ( Figure 5A), while 174 were highly downregulated ( Figure 5B). Of the 168 upregulated genes, LTBP1 and LINC00539 were commonly identified under all four conditions ( Figure 5A). The expression of LINC02355 was highly downregulated under all four conditions ( Figure 4B). Many DEGs were condition-specific. For example, 69 upregulated and 32 downregulated genes were identified only under the GP-all condition. We examined the expression of the 176 DEGs identified from the GP-all condition in the other three conditions using a heatmap ( Figure 5C). The heatmap showed that the expression of most DEGs identified in the GP-all condition was quite similar to that in the other three conditions.

Functional Classification of DEGs in Response to GP Extracts
Based on the heatmap results, most DEGs identified in GP-all showed similar gene expression patterns. Therefore, we conducted a gene enrichment analysis using only the DEGs from GP-all. DEGs from GP-all were divided into upregulated and downregulated genes. Of the 125 upregulated genes, 121 were assigned GO term annotations. Genes assigned to biological regulation (50 genes) were most abundant, followed by those associated with metabolic process (48 genes) and response to stimulus (40 genes) when annotated by the biological process ( Figure 6A). Genes associated with the membrane (37 genes) were most abundant, followed by those associated with the nucleus (28 genes) and endomembrane system (26 genes) when annotated by cellular components. Genes related to protein binding (48 genes) were most abundant, followed by ion binding (29 genes) and transferase activity (16 genes) when annotated by molecular function.
Of the 51 downregulated genes, a high number were assigned to biological processes (26 genes), responses to stimuli (19 genes), and metabolic processes (19 genes) ( Figure 6B). Many downregulated genes were associated with membranes (19 genes), protein-containing complexes (13 genes), and cell projections (nine genes). Two GO terms, endosome (four genes) and chromosome (one gene), were identified only among the upregulated genes. Compared to the upregulated genes assigned to 16 GO terms according to molecular functions, the downregulated genes were assigned to only 13 GO terms. Three GO terms, antioxidant activity (one gene), oxygen binding (one gene), and enzyme regulator activity (eight genes), were identified only in the upregulated genes.
To identify the enriched functions for the 176 DEGs identified from GP-all, we conducted gene set analyses using the WebGestalt program, which contains several databases. For this, 125 upregulated and 51 downregulated genes were subjected to gene set analysis (Tables S5 and S6). GO enrichment analyses showed a similar number of GO terms for up-(52 terms) and downregulated genes (53 terms) according to the biological processes ( Figure 7). However, based on cellular components and molecular functions, the number of GO terms for upregulated genes was higher than that for downregulated genes.

Functional Classification of DEGs in Response to GP Extracts
Based on the heatmap results, most DEGs identified in GP-all showed similar gene expression patterns. Therefore, we conducted a gene enrichment analysis using only the Of the 51 downregulated genes, a high number were assigned to biological processes (26 genes), responses to stimuli (19 genes), and metabolic processes (19 genes) ( Figure 6B). Many downregulated genes were associated with membranes (19 genes), protein-containing complexes (13 genes), and cell projections (nine genes). Two GO terms, endosome (four genes) and chromosome (one gene), were identified only among the upregulated genes. Compared to the upregulated genes assigned to 16 GO terms according to molecular functions, the downregulated genes were assigned to only 13 GO terms. Three GO terms, antioxidant activity (one gene), oxygen binding (one gene), and enzyme regulator activity (eight genes), were identified only in the upregulated genes.
To identify the enriched functions for the 176 DEGs identified from GP-all, we conducted gene set analyses using the WebGestalt program, which contains several databases. For this, 125 upregulated and 51 downregulated genes were subjected to gene set Of the 51 downregulated genes, a GO enrichment analysis revealed that nine genes associated with cell adhesion (GO:0007155) were strongly downregulated (Table S8 and Figure 9A). In addition, two genes related to the positive regulation of viral genome replication (GO:0045070) were strongly downregulated. We found that many downregulated genes were localized to diverse cellular components, such as the synapse (seven genes), neuron part (eight genes), apical part of the cell (four genes), an intrinsic component of the plasma membrane (nine genes), chloride channel complex (two genes), and exocyst (one gene) ( Figure 9B). analysis (Tables S5 and S6). GO enrichment analyses showed a similar number of GO terms for up-(52 terms) and downregulated genes (53 terms) according to the biological processes ( Figure 7). However, based on cellular components and molecular functions, the number of GO terms for upregulated genes was higher than that for downregulated genes. For upregulated genes, two GO terms, response to growth factor (GO:0070848) and heart development (GO:0007507), were the most significant functional GO terms assigned to the biological processes ( Figure 8A). The 11 genes assigned to growth factors were LTBP1, MYOCD, VTN, SOX5, ERBB4, TH, CLDN5, FGF1, GRB10, TMEM204, and SOX2 (Table S7). The nine genes related to heart development were LTBP1, SCN5A, DNAH5, MYOCD, ERBB4, TH, CLDN5, CSRP3, and GLI1. Moreover, monooxygenase activity (GO:0004497)-associated GO term was the most significant molecular function, and four genes, TPH1, CYP24A1, CYP2F1, and TH, were involved (Table S7 and Figure 8B). Many upregulated genes were targeted to several cellular components, such as synapse (nine genes), supramolecular fiber (eight genes), neuron projection (10 genes), and integral components of the plasma membrane (12 genes), Golgi apparatus (11 genes), endoplasmic reticulum (15 genes), cytoplasmic region (six genes), and NLRP3 inflammasome complex (one gene) (Table S7). Of the 51 downregulated genes, a GO enrichment analysis revealed that nine genes associated with cell adhesion (GO:0007155) were strongly downregulated (Table S8 and Figure 9A). In addition, two genes related to the positive regulation of viral genome replication (GO:0045070) were strongly downregulated. We found that many downregulated genes were localized to diverse cellular components, such as the synapse (seven genes), neuron part (eight genes), apical part of the cell (four genes), an intrinsic component of   (Table S8 and Figure 9A). In addition, two genes related to the positive regulation of viral genome replication (GO:0045070) were strongly downregulated. We found that many downregulated genes were localized to diverse cellular components, such as the synapse (seven genes), neuron part (eight genes), apical part of the cell (four genes), an intrinsic component of the plasma membrane (nine genes), chloride channel complex (two genes), and exocyst (one gene) ( Figure 9B). An enrichment analysis using the 125 upregulated genes in the KEGG database revealed that folate biosynthesis (hsa00790) was a significant metabolic pathway (Table S5). In contrast, among the 51 downregulated genes, five metabolic pathways, including neuroactive ligand-receptor interaction (hsa04080), renin-angiotensin system (hsa04614), maturity-onset diabetes of the young (hsa04950), nicotinate and nicotinamide metabolism (hsa00760), and mucin-type O-glycan biosynthesis (hsa00512) were significantly enriched (Table S6).
Enrichment analyses using the PANTHER database revealed that the circadian clock system (P00015) and vitamin D metabolism and activity pathways (P04396), were significantly enriched in the upregulated genes, while the plasminogen activating cascade (P00050) and axon guidance mediated by Slit/Robo (P00008) pathways were highly enriched in the downregulated genes.

Discussion
In this study, we established a mass production system for GP cells based on tissue culture and bioreactors and identified six metabolites in GP extracts. Three of the identified metabolites (uridine, adenosine, and guanosine) were nucleosides. A previous study suggested that guanosine and uridine have therapeutic roles in asthma and exert anti-inflammatory effects [34]. Furthermore, we identified three aromatic amino acids, tyrosine, phenylalanine, and tryptophan, in the extracts. Tryptophan and phenylalanine play important roles in the human diet. For example, tryptophan metabolism is associated with serotonin production, whereas phenylalanine is necessary for tyrosine production [35]. Notably, GP-C and GP-CR contained much higher amounts of the six metabolites than GP-P. This result suggested the enrichment of major metabolites in the extract from the cultured GP cells.
The effects of radiofrequency on plant growth have not yet been well-studied. A recent study showed that a weak radiofrequency magnetic field changed the cryptochromedependent plant growth response and gene expression in Arabidopsis [36]. We tested the possible applications of radiofrequency in plant cell culture. Similarly, we assumed that treatment with weak radiofrequency might alter the metabolites of GP cells; hence, treatment with the extract from such radiofrequency-treated GP cells may alter gene expression in HaCaT cells. However, we did not observe any significant differences in the gene expression profiles between GP-C and GP-CR, although the expression of some genes was condition-specific. We supposed that the radiofrequency treatment in our study might not have been high enough to change the metabolites of GP cells.
Based on the metabolites identified from GP extracts, it might be of interest to examine the possible effects of GP extracts on human skin cells. To this end, we compared the gene expression profiles of HaCaT cells after treatment with the three different GP extracts. The gene expression profiles in response to the treatments were similar. For example, genes highly upregulated in the GP-P condition were also upregulated in the other three conditions and vice versa. Therefore, the metabolites and their concentrations in the three GP extracts may be very similar, with minor differences. Based on these results, we concluded that the extracts from plant cell tissue culture and bioreactors could replace those derived from whole plants.
Most DEGs identified in GP-all showed similar expressions to the other three groups. Interestingly, we found that not only genes coding for proteins but also several lncRNAs were identified as DEGs after GP extract treatment. LncRNAs are non-coding RNA transcripts longer than 200 nucleotides in size, and at least 15,000 lncRNAs have been reported in genomic regions that do not include protein-coding genes in humans [37]. They also contain promoter-or enhancer-related RNA close to the gene region [37], and therefore, participate in a wide range of biological processes, including gene expression. In this study, we identified two upregulated (LINC00052 and LINC00539) and one downregulated (LINC02355) lncRNAs. Several studies have demonstrated that LINC00052 is involved in tumorigenesis, progression, and metastasis in several types of human cancers [38]. A previous study using gene expression data showed that LINC00539 might be involved in the immune response against lung adenocarcinoma [39]. These results suggest that GP extract might be associated with the immune response to tumors by inducing the expression of LINC00052 and LINC00539. However, the functional role of LINC02355 has not been reported. Human tissue-specific transcriptome data showed that LINC02355 was highly expressed in testis and thyroid tissues [40].
The most significantly upregulated gene under the four conditions was LTBP1. LTBP1 is a member of the latent transforming growth factor β (TGFβ)-binding protein (LTBPs) family, which is composed of four proteins (LTBP1-LTBP4). LTBPs associate with microfibrils to anchor TGF-β in the extracellular matrix (ECM) and are required for the correct assembly of ECM components [41]. LTBP1 plays an important role in skin and bone ECM assembly and homeostasis [42].
GP extract induced the expression of 11 genes involved in the response to growth factors (GO:0070848). Of these, the MYOCD gene encoding myocardin plays an important role as a master regulator of smooth muscle gene expression [43]. VTN gene encodes vitronectin protein, a glycoprotein in the hemopexin family that is abundantly present in the serum of the ECM and bone [44]. In addition, two genes (SOX2 and SOX5) are transcription factors in the SRY-box transcription factor (SOX) family, and their expression was highly upregulated by GP extract treatment. SOX family transcription factors are involved in diverse cellular developmental processes with functional redundancy [45]. ERBB4 encodes the epidermal growth factor receptor-4 protein and is required for normal tissue development, such as the heart and nervous system, and tumor suppression [46]. CLDN5 is a member of the claudin family and is an integral membrane protein specific to tight junctions that provide a protective paracellular barrier [47]. FGF1, which encodes fibroblast growth factor 1, plays a vital role in embryonic development, wound healing, neurogenesis, angiogenesis, and the control of type 2 diabetes mellitus [48]. GRB10 is a member of the growth factor receptor-bound protein family and a negative regulator of insulin signaling and action [49]. TMEM204 is a member of the transmembrane protein (TMEM) family and is required for angiogenesis and tumorigenesis [50].
We found that nine genes upregulated by the GP extract were involved in heart development (GO:0007507). Of these nine genes, five, namely LTBP1, MYOCD, ERBB4, TH, and CLDN5, were required for the response to growth factors (GO:0070848). Therefore, these five genes specifically function as growth factors in heart development. In addition to the other four genes associated with heart development, SCN5A, which encodes the alpha subunit of the main cardiac sodium channel, is associated with multiple cardiac disorders, such as cardiac conduction system dysfunction and dilated cardiomyopathy [51]. A previous study showed that mutations in the DNAH5 gene encoding dynein cause primary ciliary dyskinesia with outer dynein arm defects [52]. Moreover, mutations in CSRP3, a member of the CSRP family containing the LIM domain, cause hypertrophic cardiomyopathy, a genetic cardiac disease, and the most common cause of sudden cardiac death [53]. Glioma-associated oncogene homolog 1 (GLI1) is a zinc finger transcription factor involved in the hedgehog (HH) signaling pathway, required for normal cell growth and differentiation [54]. In our study, a Wikipathway analysis revealed that the HH signaling pathway (WP47) was enriched in the genes upregulated by the GP extract. Our results are consistent with those of a previous study in which the aqueous extract of GP showed positive cardiovascular effects [55]. However, aberrant HH signaling due to truncated GLI1 can cause several human cancers [56].
GP extract has been widely used as an anticancer agent, and several studies have shown the inhibitory activity of GP extracts against cancer cells in vivo and in vitro [57]. In our study, several GP extract-upregulated genes, such as ERBB4 and GL1, are associated with diverse cancers. The main compounds in GP extract were gypenosides, followed by sterols, flavonoids, and polysaccharides. Therefore, we hypothesized that GP gypenosides could play important roles as anticancer agents.
Reactome pathway analysis revealed that several genes associated with elastic fibers (R-HSA-2129379), collagen chain trimerization (R-HSA-8948216), and elastic fiber formation (R-HSA-1566948) were upregulated. Elastic fibers are important connective tissue components that provide elasticity and resilience to the skin, lungs, and blood vessels [58]. They consist of elastin and many elastin-associated microfibrillar proteins that assemble into a complex fiber network [58].
Aging is easily visible as changes in the skin. In normal human skin, the loss of elastic fibers with age results from the loss of skin thickness [59]. Moreover, the loss of skin elasticity with age is caused by the loss of fibroblast cells owing to a reduction in biosynthetic activity and modification of ECM macromolecules [59]. Many genes involved in the response to growth factors, such as LTBP4 have also been associated with elastic fibers. The genes highly upregulated by the GP extract encode components of elastic fibers and the ECM. Treatment of human skin cells with GP extract may promote the expression of genes associated with growth factors. Moreover, an increase in growth factor proteins may prevent the aging of human skin.
Among the known metabolic pathways, two upregulated genes, TH and TPH1 were associated with folate biosynthesis (hsa00790). TH encodes tyrosine hydroxylase, whereas TPH1 encodes tryptophan hydroxylase 1. Moreover, CYP24A1 encoding cytochrome P450 family 24 subfamily A member 1, is involved in vitamin D metabolism. Several studies have suggested that vitamin D and folate metabolism may be correlated with vascular health and human skin pigmentation [60]. The high expression of these three genes on treatment with the GP extract could cause the accumulation of vitamin D in the skin, and the active metabolites of vitamin D have a wide range of anti-aging and photoprotective effects on the skin [61].
In contrast, the most significant function of the downregulated genes was cell adhesion (GO:0007155). Cell adhesion molecules allow cells to adhere to each other and to the ECM [62]. Therefore, they play important roles in the formation and structural integrity of the epithelium of human skin [63]. Of the identified genes associated with cell adhesion, DDR2 encoding discoidin domain receptor 2, a collagen receptor, regulates cell proliferation, and the absence of DDR2 results in dwarfism in mouse mutants [64]. Thus, the downregulation of DDR2 by the GP extract may have negative effects on cell adhesion. Another identified gene downregulated by GP extract was ADAMTS18, which is a member of the disintegrin metalloproteinase with thrombospondin domains (ADAMTS) family composed of 19 genes. ADAMTS genes participate in the assembly and degradation of ECM during development, morphogenesis, tissue repair, and remodeling [65]. Loss of ADAMTS18 function promotes the growth, migration, and metastasis of melanomas [66]. Thus, the downregulation of ADAMTS18 gene expression by GP extract may have a positive effect on cell adhesion.
Enriched cellular components in both up-and downregulated genes were associated with the synapse (GO:0045202), neuron projection (GO:0043005), neuron part (GO:0097458), and an intrinsic component of the plasma membrane (GO:0031226). As expected, the GP extract affected the expression of genes associated with the plasma membrane. Notably, many DEGs were associated with neurons. Consistent with our results, previous studies have shown that GP has positive effects on injured spinal motor neurons [67] and neuroprotective effects in a rat model of cardiopulmonary resuscitation [68].
GP extract also induced the expression of genes involved in several metabolic pathways, such as amino acid metabolism, tryptophan metabolism, and heme biosynthesis, which are essential processes in human organs. For example, tryptophan metabolism is highly associated with inflammation, energy homeostasis, and brain function [69].
Taken together, our study provides comprehensive expression profiles of HaCaT cells in response to GP extracts and reveals many crucial biological processes that may be affected by GP extracts.

Conclusions
A wide range of effects of GP extracts has been reported in numerous studies. However, changes in the gene expression profile in human cells in response to GP extract have not yet been investigated. In this study, we identified six metabolites (uridine, adenosine, guanosine, tyrosine, phenylalanine, and tryptophan) in GP extracts, and compared the transcriptome profiles of HaCaT cells in response to treatment with GP extracts obtained through three independent methods from plants, calli, and calli after radiofrequency treatment. Some genes showed condition-specific expression; however, most DEGs identified in the GP-all condition exhibited a similar gene expression in the three GP extracts. The most significantly upregulated gene under the four different conditions was LTBP1, which is required for skin and bone ECM assembly as well as homeostasis. Moreover, we identified two upregulated (LINC00052 and LINC00539) and one downregulated (LINC02355) lncRNAs. On combining the data, we identified 125 upregulated and 51 downregulated genes in response to all extracts. The enriched functions of the 125 upregulated genes were associated with response to growth factors and heart development. Moreover, genes encoding components of elastic fibers and ECM were upregulated by the GP extract. Some upregulated genes, such as ERBB4 and GL1, are also associated with diverse cancers. Furthermore, the upregulation of genes related to folate biosynthesis and vitamin D metabolism suggests the functional roles of the extracts in anti-aging and photoprotective effects on the skin. In contrast, the most enriched function among the downregulated genes was cell adhesion. Moreover, many DEGs were targeted to synapse and neuron projections, indicating their possible roles in neuroprotection. In summary, we provide a comprehensive overview of gene expression changes in HaCaT cells after treatment with GP extracts and reveal their functional roles in anti-aging and photoprotective effects on the skin using RNA sequencing.
Supplementary Materials: The following supporting information can be downloaded at https: //www.mdpi.com/article/10.3390/life13020423/s1: Table S1: List of differentially expressed genes by the GP extract from the plant; Table S2: List of differentially expressed genes by the GP extract from the callus; Table S3: List of differentially expressed genes by the GP extract from the callus treated with radiofrequency; Table S4: List of differentially expressed genes by the three GP extracts; Table S5: Enriched functions in 125 upregulated genes; Table S6: Enriched functions in 51 downregulated genes; Table S7: Representative functions enriched in 125 upregulated genes and related gene list; Table S8: Representative functions enriched in 51 downregulated genes and related gene list. Figure S1: HPLC chromatograms of Gynostemma pentaphyllum and tandem mass spectrometry of each peak fraction.

Data Availability Statement:
The raw data were deposited in the NCBI SRA database with the following accession numbers: SRR22342813-SRR22342824.

Conflicts of Interest:
The authors declare no conflict of interest. Seung Hye Paek, Soo-Yun Kim, Sung Joo Jang, Sak Lee, Jeong Hun Lee, and Sang Hyun Moh were employed by the company BIO-FD&C Co., Ltd. The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest. The company had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.