Deciphering the Dose-Dependent Effects of Thymoquinone on Transcriptomic Changes and Cellular Proliferation in Glioblastoma

1. Introduction

Approximately 1.9 million cases of cancer are estimated to be diagnosed in 2023, with over 600,000 individuals expected to die from the disease. Brain and other nervous system tumors account for 3% of all cancer-related deaths in both sexes. Glioblastoma multiforme, the most common primary malignant brain tumor in adults, has an incidence rate of 3.26 per 100,000 populations. The term 'glioblastoma' is reserved for the IDH-wildtype diffuse gliomas grade 4, as defined by Louis et al. in 2021 [1]. These tumors are characterized by rapid proliferation and poor differentiation. Genomic mutational studies, including those by The Cancer Genome Atlas (TCGA) Research Network in 2008, show that 80–90% of glioblastoma tumors harbor mutations affecting the PI3K, AKT, RAS, MAPK, p53, and/or RB signaling pathways. Glioblastomas are also highly angiogenic cancers with markedly increased VEGF expression, as noted by Reardon and Wen in 2015 [2]. However, the therapeutic outcome following antiangiogenic therapy is highly variable in improving overall survival, according to studies by Chinot et al. in 2014 and Gilbert et al. in 2014 [3,4]. The current standard of care for glioblastoma patients, as outlined by Batash et al. in 2017 [5] and the National Comprehensive Cancer Network (NCCN) in 2015, includes surgical resection, radiotherapy, and adjuvant chemotherapy with temozolomide. Despite recent advances in diagnosis and treatment, the overall survival for patients with glioblastomas is very poor, with a median overall survival of 13.5 months, as reported by Marenco-Hillembrand et al. in 2020 [6] and Ostrom et al. in 2022 [7]. Treatment outcomes in GBM patients are hindered by several factors, including high rates of de novo and acquired resistance, incomplete surgical resection due to the tumor's highly infiltrative nature, low selective cytotoxicity, serious side effects, an immunosuppressive microenvironment, redundancy in aberrantly activated cell signaling pathways, and the challenge of effectively delivering therapeutics due to the blood–brain barrier. These challenges are highlighted in studies by Davis in 2016 [8] and Reardon & Wen in 2015 [2].

In recent years, there has been a surge of interest in exploring the therapeutic potential of medicinal plants for anticancer agents within the scientific community. Many conventional anticancer drugs are derived from bioactive constituents, including Vinca alkaloids (vinblastine and vincristine), cytotoxic podophyllotoxins, and Taxol [9]. Nigella sativa Linn., commonly known as Blackseed or Black cumin, is an annual flowering plant native to South and Southwest Asia, and is prevalent in Northern Africa, the Middle East, and Southern Europe [10]. It is one of the most promising and extensively studied medicinal plants for anticancer agents. Phytochemical studies have identified several bioactive compounds in Nigella sativa, including sterols, saponins, phenolic compounds, alkaloids, novel lipid constituents, fatty acids, and volatile oils [11]. The most pharmacologically active constituent of Nigella Sativa's volatile oil is Thymoquinone (TQ) (2-methyl-5-isopropyl-1,4-benzoquinone) [12,13]. TQ is extracted from the essential oil of Nigella Sativa seeds using supercritical fluid extraction [14] and has been found to possess antimicrobial, antihistaminic, antidiabetic, anti-inflammatory, antioxidant, hypolipidemic, and anticancer properties [15,16,17]. Several in vitro and in vivo studies have illustrated the cytotoxic effects of TQ on various cancer cell lines, including breast adenocarcinoma, leukemia, lung cancer, colorectal carcinoma, pancreatic cancer, osteosarcoma, prostate cancer, and glioblastoma [17]. TQ affects multiple pathways and processes relevant to cancer, exhibiting anticancer activity through numerous mechanisms, specifically by showing selective antioxidant and oxidant activity, interfering with DNA structure, affecting carcinogenic signaling molecules/pathways, and immunomodulation [18]. In vivo studies on mice have suggested that TQ can cross the blood-brain barrier due to its size and lipophilicity [19]. Although some underlying mechanisms for the anticancer activities of TQ in glioblastoma cells have been identified, many potential pathways and targets remain to be fully elucidated. Research on the exact molecular pathways of its anticancer activities in glioblastoma cells is still in its early stages. Most research on TQ in glioblastoma cells has been conducted using the U-87 MG cell lines, with other glioma cell lines yet to be extensively studied. The aim of this study was to use RNA sequencing to examine changes in the expression of genes in A172 glioma cells in response to treatment with TQ at concentrations where significant cell death occurs. We then explored the pathways associated with those differentially expressed genes.

2. Materials and Methods

2.1. Cell Culture, Treatment and Cell Viability Assays

The glioblastoma cell line, A172, was purchased from the American Type Culture Collection (ATCC #CRL-1620) (Manassas, VA, USA). Cells were maintained in Dulbecco's modified Eagle medium (DMEM) (ATCC, Manassas, VA) containing 4,500 mg/L D-glucose, L-glutamine, and 110 mg\L sodium pyruvates in 75 or 175 cm2 cell culture flasks. The media was supplemented with 10% (v/v) Fetal bovine serum (FBS) (Atlas Biologicals, Fort Collins, CO). Additionally, an antibiotic-antimycotic mixture: of penicillin (10 000 U/mL), streptomycin (10mg/ mL), and amphotericin (25µg/mL) (Millipore Sigma Life Sciences, Burlington, MA) was also supplemented to the media. All cell lines were maintained in a humidified incubator at 37^oC with 5% atmospheric CO2. For passaging the cells, trypsin-EDTA (ATCC, Manassas, VA) was used. 50mM stock solution of thymoquinone (Sigma Aldrich, MA) was prepared in filter sterilized DMSO (Millipore Sigma). Using the complete DMEM medium, suitable working concentrations were prepared from the stock solution. The final concentration of DMSO was 0.1% for all treatments. Cells were treated with different concentrations of thymoquinone ranging from 0 to 50 µM for 48 hours. After overnight attachment of cells, the media was removed and replaced with the corresponding volume of TQ-containing and untreated media. To test the effect of different concentrations of thymoquinone on cell viability, Toluidine blue staining Assay and Cell Counting Kit-8 was used.

2.2. Toluidine Blue Staining

A172 cells were seeded into 96 well plates and incubated overnight. 24 hours after incubation the media was removed, 100µl of thymoquinone-treated and untreated media was added to each well. Cells were treated with 0, 10, 25, and 50µM concentrations of TQ for 48 hours. Following the TQ treatments, the media was removed from each well. Cells were washed once with phosphate-buffered saline (PBS) (ATCC). 100µl of toluidine blue staining solution (1% toluidine blue (LabChem Inc, Zelienople, PA) and 1% borax (LabChem)) was added to each well and incubated for 20 minutes at room temperature. After 20-25 minutes, the staining solution was directly discarded to the sink. The remaining unbound solution was washed by immersing the plate in water. The plate was inverted and left to dry. The next day, the cells absorbance at 600nm was measured using SpectraMax iD3 Multi-Mode Microplate Reader (Molecular Devices, San Jose, CA).

2.3. Cell Counting Kit-8

Cell Counting Kit-8 (CCK-8) is a highly sensitive colorimetric assay widely used for the determination of a number of viable cells in cell viability and cytotoxicity assays. CCK-8 utilizes Dojindo's highly water-soluble tetrazolium salt, WST-8. The WST-8 assay was performed according to the manufacturer's protocol (Dojindo Laboratories, Gaithersburg, MD). A172 cells were seeded into 96 well plates at a density of not more than 5000 cells per well. Cells were incubated overnight for reattachment to the culture plate. The next day, the media was removed and replaced with fresh thymoquinone-treated and untreated media. Control, 10µM, 25µM, and 50µM of thymoquinone were used as treatment. After 48 hours of treatment, 10 µL of CCK-8 solution (Dojindo Laboratories) was added to each 96 wells and incubated for 3 hours at 37°C and 5% CO2. After incubation, the absorbance was measured at 450nm using SpectraMax iD3 Multi-Mode Microplate Reader (Molecular Devices).

2.4. RNA Extraction, Library Preparation, and RNA-Sequencing

A172 cells were seeded in 60mm culture dishes and incubated overnight. After 24 hours of incubation, the cells were treated with 25 and 50 µM concentrations of thymoquinone for 48 hours. Total RNA was extracted using Qiagen RNeasy Total RNA Isolation Kit (QIAGEN, Valencia, CA). The concentration and integrity of RNA samples were evaluated using a Nanodrop spectrophotometer (Thermo Fisher Scientific) and 0.7% % agarose gel electrophoresis. RNA was quantified using Qubit 3.0 (Thermo Fisher Scientific, USA). Then, the next-generation sequencing libraries were constructed as per the manufacturer's protocol ((NEBNext®Ultra™ RNA Library Prep Kit for Illumina®; Illumina, New England Biolabs Inc., Ipswich, MA, USA, Cat #E7770S/L, #E7775S/L). 1 µg of total RNA from each of the three biological replicates were used for the library preparation. The Poly(A) mRNA isolation was performed using NEBNext Poly(A) mRNA Magnetic Isolation Module (NEB #E7490). Similarly, NEBNext First-Strand Synthesis Reaction Buffer and NEBNext Random Primers mix were used for the fragmentation and priming of mRNA. First-strand cDNA was synthesized using NEBNext First Strand Synthesis Enzyme Mix, and the second-strand cDNA was synthesized using NEBNext Second Strand Synthesis Enzyme Mix. The double-stranded cDNA was purified using SPRIselect Beads. Thus, purified double-stranded DNA was treated with an End Prep Enzyme Mix to repair any broken or damaged ends and add a single A-base to the 3' end of the DNA. This was followed by the adapter ligation on both ends of cDNA. The size selection of the adapter-ligated DNA was performed using SPRIselect Beads to remove any DNA fragments that are too short or too long for efficient sequencing. PCR enrichment of size-selected DNA was performed using i5 and i7 primers. Purification of PCR products was performed using SPRIselect beads. The final library quality and insert size were determined using a bioanalyzer (Invitrogen, USA), and the library was quantified using a Qubit fluorometer (Invitrogen, USA). The library was diluted to 4nM concentration and sequenced using Illumina's NextSeq 500 platform with paired-end sequencing chemistry. The resulted image files in the BCL format were converted to FASTQ with 2x150 bp reads using the bcl2fastq tool (Illumina, USA).

2.5. Differential Gene Expression and Pathway Analysis

The raw data was obtained in FASTQ format. The quality of the raw data was assessed using the tool FastQC. For quality control, these raw reads were processed by Trimmomatic (v0.30) to remove low-quality bases (Quality Score <30) and those containing the adapter sequences (Bolger et al., 2014). Clean data obtained after quality trim was used for downstream analysis. The reference human genome sequences (GRCh38.p13) and the gene annotation files were then downloaded from gencode website (https://www.gencodegenes.org/human/). STAR (version 2.7.3a) [20] was used to index the reference human genome sequences in order to create the necessary data structures for efficient alignment. After creating the reference genome indexed file, the clean data were aligned to them via STAR aligner [20]. In the next step, the aligned transcriptomic data in SAM format was converted into BAM format using SAMTools. The aligned BAM file was indexed using the SAMtools. The read count table was generated from the BAM alignment file and general feature format (GFF) of genome annotation using the HTSeq R package [21]. The Differentially expressed genes (DEGs) among different experimental pair-wise combinations were identified using the DESeq2 [22]. The DEGs were filtered based on the minimum Log2 Fold-Change (Log2FC) greater than 1 or less than -1 with a false discovery rate (FDR)⩽0.05. GO and Pathway Enrichment analyses were performed using ShinyGO [23]. Similarly, heatmap was generated using statistical package pheatmap in RStudio, venn-diagram using Venny 2.1, and volcano plot was generated using online tool Galaxy Bioinformatics.

2.6. Real-time reverse transcription polymerase chain reaction (qRT-PCR)

Briefly, cells were treated with control, 10, 25, and 50µM concentrations of TQ for 48 hours. Total RNA was isolated and purified. Then, 500ng of total RNA was reverse transcribed into cDNA using random primers in a ImProm-II Reverse Transcription System (Promega, Madison, WI) following the manufacturer protocol. The RT reaction was carried out in a total volume of 20 µl. The quality and quantity of cDNA were assessed using a Nanodrop spectrophotometer (Thermo Fisher Scientific). cDNA from 27ng of RNA template and 100μM of each target genes primers was used. The mRNA expression level of the target gene was normalized to 18S rRNA mRNA expression. Fold change in the mRNA expression relative to the control was calculated using the RQ Study component of the ABI SDS v1.2.3 software and the comparative method based on 2-ΔΔCt values [24].

2.7. Statistical Analysis

To test for significance, a one-way analysis of variance (ANOVA) was used. This test was used to determine whether there were statistically significant differences among the means of control and treatment groups. Data were presented as means + SEM. The confidence interval was 95%, and the P<0.05 was considered significant. For the multiple comparisons of the means, Tukey HSD post hoc test was conducted on R Studio.

3. Results

3.1. TQ treatment for 24 hours had no significant inhibitory effect on A172 cell viability.

Since, the gene expression fold changes were not statistically significant at lower concentrations (2.5 µM and 5 µM) of TQ, we increased the TQ treatment concentration as well as the treatment time. We treated the A172 cells with increasing concentrations of thymoquinone (10, 25, and 50 µM) and treated to 24 hours as well as 48 hours. We then used WST-8 and Toluidine Blue assay to determine the effect of TQ on A172 cell viability. As determined by the WST-8 assay in Figure 1, no significant growth inhibitory effect of TQ on A172 cells were observed across control, 10 µM, 25 µM, and 50 µM TQ treatment for 24 hours. If anything, TQ was observed to increase the numbers of viable cells in all the treatments, although any differences were small. Among the different concentrations, the proliferative response of A172 cells was higher at 25 µM. The 25 µM treatment significantly increased the viability of cells when used for 24 hours. In contrast, the 10 and 50 µM treatments did not significantly increase the proliferation or number of A172 cells.

3.2. TQ treatment for 48 hours inhibited the viability of A172 cells in dose dependent manner.

We further investigated the cytotoxic effects of 10, 25, and 50 µM TQ on A172 cells for 48 hours using WST-8 and Toluidine Blue assay. As shown in Figure 1, thymoquinone treatment for 48 hours inhibited and increased numbers of viable cells as analyzed by WST-8 assay. The addition of 10 µM TQ numerically, but not significantly, increased the viability of A172 cells. In contrast, 25 and 50 µM treatment for 48 hours resulted in dose-dependent inhibition of cell viability, with this being statistically significant at 50 µM. The colorimetric assay using Cell Counting Kit-8 is not cell permeable, so it is much less cytotoxic to the cells than other cytotoxic assays. Therefore, it is possible to continue further experiments with the same cells (Dojindo Molecular Technologies, Inc., 2020). Hence, after performing the WST-8 assay, we used the same cell and proceeded directly to the toluidine blue cell viability assay. We simultaneously fixed and stained the cells with 100µl of toluidine blue staining solution and measured the absorbance at 600 nm. As shown in Figure 1, TQ decreased the number of viable cells in a dose-dependent manner. The growth inhibition at all concentration (10 µM, 25 µM, and 50 µM) was statistically significant. The greatest decrease in the viability of A172 cells was observed at the TQ concentration of 50 µM.

Overall, the results obtained from both assays were similar except for the 10 µM treatment. The results of the toluidine blue assay were more dramatic, although other cytotoxic assays like MTT and WST-8 are generally recognized as more sensitive. We took the result obtained from our WST-8 assay as a reference for our further experiment. One limitation of the toluidine blue assay is there are chances of experimenter error, like rigorous washing of plates with water. This might have washed out some living cells attached to the cell plate and affected our readings perhaps due to changes in cell adhesion. Our results indicated that TQ inhibits cell viability and/or proliferation of A172 cells at higher concentrations for a longer exposure (48 hours).

Figure 1. Viability of A172 cells after treatment with control, 10, 25, and 50 µM concentrations of TQ for 24 hours (A) and 48 hours (B and C). Cell viability was measured using WST-8 assay (A and B) and Toluidine Blue assay(C). The data are presented as the mean ±Standard Error of mean. 0.1% DMSO was used for each treatment including control. (n=3 for each treatment, *p<0.05, **p<0.01, ***p<0.001).

Figure 1. Viability of A172 cells after treatment with control, 10, 25, and 50 µM concentrations of TQ for 24 hours (A) and 48 hours (B and C). Cell viability was measured using WST-8 assay (A and B) and Toluidine Blue assay(C). The data are presented as the mean ±Standard Error of mean. 0.1% DMSO was used for each treatment including control. (n=3 for each treatment, *p<0.05, **p<0.01, ***p<0.001).

3.3. 2.5 µM and 5 µM TQ for 24 hours did not cause differential gene expression in A172 cells.

Initially we treated A172 cells with control, 2.5 µM and 5 µM TQ for 24 hours. We then extracted the RNA, inspected quality and quantity using Agilent 2100 Bioanalyzer and Qubit 4 Fluorometer respectively, prepared library, and performed RNA-seq. In this case, we wanted to check for gene expression changes at concentrations of TQ that did not cause significant differences in the numbers of viable cells. Altogether 352854966 raw reads were obtained. After removing the low-quality reads (Phred score QV<30) using a read trimming tool Trimmomatic, quality-filtered reads were mapped to the Human reference genome (GRCh38.p13) https://www.gencodegenes.org/human [25]. Using the DESeq2 in R, we quantified the gene expression changes. We found that 2.5 µM, and 5 µM TQ for 24 hours didn't cause any significant changes in the expression of genes in A172 cells. More details of this RNA-seq data are presented in Table 1.

3.4. TQ treatment for 48 hours induced differentially expressed genes in A172 cells.

A172 cells were treated with control (0.1% DMSO), 25 µM, and 50 µM TQ concentration for 48 hours, RNA was extracted, and RNA-sequencing was performed as described in materials and method section. A total of 332,999,712 raw reads were primarily generated from RNA-seq. After removing the low-quality reads (Phred score QScore<30) and those containing sequencing adapters, 321,194,991 clean reads were remained. We then mapped the remaining clean reads to the reference human genome https://www.gencodegenes.org/human. The details of data quality, filtration, mapping percentage are presented in Table 2. Using DESeq2 in RStudio version 3.0., we determined the genes with significant differential gene expression profiles between the control and different treatment conditions (Control vs. 25µM and control vs. 50 µM). Differential gene expressions were determined as the Log2 Fold-Change (Log2FC) greater than 1 or less than -1 with a false discovery rate (FDR)⩽0.05. Compared to the control, 1548 genes were differentially expressed in 25 µM TQ for 48 hours in which 684 genes were upregulated and 864 were downregulated. These 1548 DEGs are visually represented in a volcano plot in Figure 2. In 50 µM TQ for 48 hours, 2797 genes were differentially expressed, of which 1528 were upregulated and 1269 were downregulated. Volcano plots for 2797 DEGs in 25 µM treatment are represented in Figure 2. We also determined the common significant DEGs in 25 and 50 µM treatment for 48 hours. In total, we found 1202 common differentially expressed genes. This includes 580 common significantly upregulated genes and 622 common significantly downregulated genes. Similarly, 104 genes were found to be significantly upregulated only on 25 µM for 48 hours. Also, 948 genes displayed significant upregulation solely on 50 µM treatment for 48 hours. Furthermore, 242 and 647 genes exhibited significant downregulation exclusively on 25 µM and 50 µM treatment for 48 hours, respectively. Venn diagram representing common and unique significantly upregulated genes and common and unique significantly downregulated genes on two treatment conditions are presented in Figure 3. Heat map of significant DEGs in 25 µM and 50 µM TQ treatment is presented in Figure 4.

3.5. KEGG Pathway Enrichment Analysis

We performed KEGG pathway enrichment analysis for common significant DEGs in 25 µM and 50 µM TQ Treatment for 48 hours. We determined the important enriched biochemical, metabolic, and signal transduction pathways related to the activity of thymoquinone on human glioblastoma cells (A172). Using p<0.05 and a pathway size minimum of 15, we found that the common upregulated genes were primarily enriched in 85 pathways. The top 20 highly enriched pathways are shown in Figure 5. In summary, those common upregulated genes were enriched in several pathways known to be involved in tumorigenesis and the progression of glioblastoma. Some of these pathways are P53 signaling pathway, ECM-receptor interactions, HIF-1 signaling pathway, TNF signaling pathway, WNT signaling pathway, MAPK signaling pathway, PI3-AKT signaling pathway, pathways in cancer, focal-adhesions, and cytokine-cytokine interactions.

Similarly, KEGG pathway analysis revealed that the common downregulated genes were enriched in the pathways related to ECM-receptor interaction, AGE-RAGE signaling pathway in diabetic complications, protein digestion and absorption, oxidative phosphorylation, calcium signaling pathway, focal adhesion, chemical carcinogenesis, and growth hormone synthesis, secretion, and action. Top 20 significantly enriched pathways are presented in Figure 5. The details of the enriched pathways, genes and fold changes are presented in Table 3 and Table 4. The KEGG pathways with the identified genes from the current study are presented in the Supplementary Figures 1–5.

3.6. Gene Ontology Enrichment

Gene ontology enrichment analysis of common significant DEGs in 25 µM and 50 µM TQ Treatment for 48 hours was conducted in order to explore the functions of those DEGs. The gene functions were categorized into three classes: biological process, cellular components, and molecular functions. The common upregulated genes in the category of the biological process are the genes related to MAPK cascade, blood vessel development, cell migration, negative regulation of signal transduction, regulation of programmed cell death, and regulation of protein phosphorylation. Similarly, some cellular components highly enriched in those common downregulated genes include the intrinsic component of synaptic vesicle membrane, anchored component of membrane, extracellular matrix, and external encapsulating structure. Furthermore, the predominantly enriched molecular function terms were DNA binding, bending, insulin-like growth factor binding, extracellular matrix binding, growth factor activity, and cytokine activity.

Similarly, the common downregulated genes were significantly enriched in the biological process related to an extracellular matrix organization, extracellular structure organization, external encapsulating structure organization, cell adhesion, and focal adhesion. Similarly, the markedly enriched cellular components were fibrillar collagen trimer, collagen-containing extracellular matrix, actin cytoskeleton, extracellular matrix, external encapsulating structure, and integral component of postsynaptic density membrane. Similarly, the common downregulated genes are enriched in the molecular functions related to collagen binding, peptide hormone binding, and extracellular matrix structural constituents. The 10 most significantly enriched biological processes, cellular components, and molecular functions are presented in Figure 6.

3.7. qRT-PCR Verification

To further validate the gene expression results from Illumina RNA-Sequencing data, qRT-PCR was performed to investigate the mRNA expression of the gene WNT7B, CHAC1, DUSP5, and CD300A in TQ treated A172 cells. Briefly, the same total RNA used for the RNA-seq was used for qRT-PCR. The details of the primers and sequences are presented in the Supplementary Table S1. All these genes are shown to be significantly upregulated by qRT-PCR result which is consistent with the RNA-seq findings. Gene expression fold change in NGS and qRT-PCR is presented in Table 5. qRT-PCR results are shown graphically in Figure 7. These findings demonstrate that the expression pattern of all genes is consistent with Illumina RNA-Sequencing data, confirming the validity of our experimental results.

4. Discussion

In our study investigating the impact of thymoquinone (TQ) on A172 glioma cell viability, we observed a dose-dependent inhibition of cell growth at higher concentrations, particularly after 48 hours of treatment. These results are consistent with existing literature, such as Ballout and Gali-Muhtasib, 2020 [26], who reported similar time and concentration-dependent effects of TQ in neuroblastoma cells, and Racoma et al. (2013) [27], who noted TQ's dose-dependent apoptotic induction in cancer cells. The dual role of TQ in modulating cell proliferation, as seen in our findings and echoed by Fatfat et al. (2021) [28], highlights its complex biological activity. Our reliance on the WST-8 assay for its reliability, as also suggested by Eid et al. (2023) [29], helped mitigate potential experimental errors common in assays like Toluidine Blue, as discussed by Adilovic et al. (2020) [30]. Additionally, our observation of TQ's selective cytotoxicity aligns with Mokashi 2004 [31], who emphasized its potential as a cancer therapeutic due to its ability to selectively target cancer cells.

In this study, TQ was found to inhibit the extracellular matrix-receptor (ECM-receptor) interaction pathway in A172 glioblastoma cells, affecting 8 differentially expressed genes. The ECM, which forms the cellular microenvironment and consists of glycoproteins, proteoglycans, and glycosaminoglycans [32], is often overexpressed in glioma cells. Components such as hyaluronic acid, brevican, tenascin-C, fibronectin, thrombospondin, and specific integrins and receptors are known to promote cell adhesion and migration [33]. These components, particularly HA and collagens, also act as barriers that limit drug diffusion and penetration into tumors [34]. Additionally, collagen can inhibit T cell migration at tumor sites, impairing immune surveillance [35]. TQ treatment's downregulation of several ECM components suggests potential therapeutic effects, such as decreased cell adhesion, impaired cancer cell migration, and hindered invasion through the basement membrane.

Further we noted, TQ inhibited the calcium signaling pathway in A172 glioblastoma cells, as evidenced by the downregulation of genes associated with this pathway following treatment with 25 and 50 µM TQ for 48 hours. Calcium signaling, crucial in controlling diverse cellular functions and linked with the progression of GBMs and other cancers [36,37], saw notable gene expression changes. CAMK2A, a key kinase in Ca2+-induced signaling and cell cycle progression, was significantly downregulated, aligning with findings by Colomer & Means (2007) [38] and Takemoto-Kimura et al. (2017) [39] on its role in proliferation. Studies by Wang et al. (2022) [40] and Yu et al. (2021) [41] further emphasize its impact on glioma cell functions. The TACR1 gene, encoding a receptor linked to various cancers [42], also showed decreased expression, highlighting TQ's potential in cancer cell apoptosis. Moreover, the reduced expression of GRM1, involved in cell cycle arrest and apoptosis in glioma [43,44], and the downregulation of PLCG2, known for its overexpression in glioma and role in intracellular Ca2+ signaling [45,46], further underscore TQ's multifaceted impact on glioblastoma cell signaling and tumor progression.

In our study, KEGG pathway analysis indicated that higher concentrations of Thymoquinone (TQ) (50 µM) significantly downregulated genes in the PI3K-Akt signaling pathway in A172 glioblastoma cells, involving 18 differentially expressed genes like LPAR2, FGF10, COL9A3, PDGFB, IL2RB, IL7, PDGFRB, ITGB6, PPP2R2B, COL1A2, THBS3, PDGFD, GNG7, ERBB4, COL4A5, LAMA2, COL6A6, ITGA1. This pathway, crucial in processes like proliferation, differentiation, migration, metabolism, and survival, is often dysregulated in GBM due to RTK gene mutations and PI3K pathway activation. Genes like PDGFB, ERBB4, and IL-7, which play roles in activating PI3K-Akt signaling [47,48,49,50], were also downregulated. Interestingly, we observed the upregulation of RASD1, a gene that inhibits glioma cell migration and invasion by deactivating the AKT/mTOR signaling pathway [51]. This might contribute to the observed inhibition of the AKT pathway. Additionally, PPP2R2B, often downregulated and methylated in GBM, was also downregulated in TQ-treated cells, correlating with shorter survival in GBM patients [52]. These findings suggest TQ's potential as a therapeutic target in GBM by influencing key pathways and gene expressions.

In our study, we observed that 25 and 50 µM TQ treatments for 48 hours upregulated the P53 signaling pathway in A172 glioblastoma cells, without altering the level of p53 itself. This suggests TQ modulates the p53 pathway by influencing the expression of its downstream target genes. Specifically, TQ upregulated genes like PMAIP1 (NOXA), GADD45A, BBC3 (PUMA), and DUSP5. PMAIP1 and BBC3, transcription targets of p53 involved in DNA damage-induced apoptosis [53,54,55,56], are pro-apoptotic proteins of the Bcl-2 family, which, when activated, bind and inhibit anti-apoptotic proteins. The downregulation of the anti-apoptotic protein Bcl-XL was also noted. GADD45A, rapidly induced in response to DNA damage [57,58], and studies have shown its importance in apoptosis induction by anticancer agents [59,60,61]. DUSP5, a dual-specificity phosphatase, deactivates protein kinases and its overexpression has been linked to suppressed growth in various cancer cells [62]. The increase in proapoptotic proteins and decrease in anti-apoptotic proteins in TQ-treated cells suggest TQ may enhance mitochondrial membrane permeability, potentially activating the intrinsic pathway of apoptosis in glioblastoma cells.

In our study, TQ treatment of A172 glioblastoma cells with 25 and 50 µM concentrations for 48 hours significantly upregulated two candidate tumor suppressor genes: Sprouty RTJ Signaling Antagonist 4 (SPRY4) and Brain Expressed X-Linked 2 (BEX2). SPRY4, upregulated by 2.5-fold, is part of the Spry protein family, known to be repressed in GBM and often missing or deleted in gliomas [63]. As a negative regulator of MAPK activation, its ectopic expression has been shown to inhibit proliferation and migration in GBM cells, suggesting its role as a tumor suppressor [64]. BEX2, on the other hand, showed an upregulation by 3.2-fold in 50 µM TQ treatment and 2.5-fold in 25 µM treatment. Part of the BEX gene family, it was found to be silenced in U-87 and primary glioma cell lines, with its re-expression increasing sensitivity to chemotherapy-induced apoptosis and demonstrating tumor suppressor effects [65]. However, some studies show contrasting results about BEX2, such as its high expression in glioma tissue and its role in inhibiting glioma cell migration and invasion by affecting β-catenin levels [66]. These findings suggest that TQ may play a role in modulating tumor suppressor gene expression, contributing to its potential therapeutic effects in glioblastoma treatment.

In our study on A172 glioblastoma cells, TQ regulated the expression of key components and regulators of the Wnt signaling pathway. TQ increased the expression of Wnt Family Member 7B (Wnt7B) by 2.2-fold and decreased WNT6 expression by 1.6-fold. The Wnt pathway, known for its activation in glioblastoma, plays a crucial role in processes like cell proliferation, apoptosis, migration, and invasion. Studies have indicated varied expression levels of Wnt pathway components in gliomas, such as elevated WNT3A and 5A and decreased WNT7B [67], and the oncogenic role of overexpressed WNT6 [68]. In our study, Phophodiesterase 2A (PDE2A) was upregulated, a gene known to suppress Wnt/β-catenin signaling in glioma stem-like cells by modulating cAMP accumulation and GSK-3β phosphorylation [69]. Additionally, TQ treatment led to downregulation of Sema3C, an overexpressed gene in most GBMs that activates the canonical Wnt pathway [70]. Our findings suggest that TQ modulates Wnt signaling by altering the expression of its core components (WNT7B and WNT6) and regulators (PDE2A and Sema3C), thereby potentially downregulating this oncogenic pathway and contributing to its anticancer mechanism in vitro.

In our study on A172 glioblastoma cells, TQ treatment led to the upregulation of important genes like CHAC1 and DNER, known for their roles in GBM development. CHAC1, a cytosolic protein, typically downregulated in GBM cell lines, was identified as a key gene upregulated by TMZ treatment, enhancing glioma apoptotic death and inhibiting Notch3-mediated pathways [71]. DNER, a noncanonical Notch ligand, was found to suppress glioma growth by inhibiting the oncogene TOR4A [72] and hindered the growth and induced differentiation of GBM-derived neurospheres [73]. Conversely, TQ treatment resulted in the downregulation of several genes overexpressed in GBM cells, including potential oncogenes like AEBP1, MIAT, GHR, LMO1, ELF3 [74,75,76,77,78], and genes involved in tumor proliferation and migration such as EPHA4, COL3A1, PCDH10, ROBO1, ADAMTS5, PCDH18, ST8SIA1 [79,80,81,82,83,84]. Additionally, TQ significantly downregulated genes associated with poor GBM prognosis, including PCSK5, KCNC1, MXRA5, SEMA3C, MFAP2, MTERF2, KDM2B, FOXP2 [85,86,87,88,89,90,91], indicating TQ's potential in targeting key molecular pathways in GBM.

5. Conclusions

The primary aim of our research was to determine the effective concentration of Thymoquinone (TQ) on glioblastoma multiforme (GBM) cell viability and to utilize RNA sequencing (RNA-seq) at this concentration to identify differentially expressed genes and associated pathways. We discovered that low TQ concentrations (below 10 µM) did not significantly alter gene expression in glioblastoma cells within the exposure times tested. However, a decrease in viable cell numbers was observed at higher concentrations (25 and 50 µM). RNA-seq analysis revealed numerous differentially expressed genes, and pathway enrichment analysis showed that TQ exerts its anticancer effects by modulating pathways such as PI3K-AKT, calcium signaling, focal adhesion, ECM-receptor interaction, and by downregulating potential oncogenes and upregulating the P53 signaling pathway and tumor suppressor genes. TQ also influenced the expression of core components and regulators of Wnt signaling in GBM cells. Our study provides valuable insights into TQ's anticancer mechanisms and supports its potential as a therapeutic agent in glioblastoma treatment. However, these findings are specific to A172 cells and might not directly translate to other cell lines or clinical settings. Future research should include in vitro experiments with different GBM cell lines and in vivo studies to validate these results. Investigating the functional consequences of downregulated genes and identifying specific molecular targets of TQ within the affected pathways will enhance understanding. While we initially hypothesized that TQ would inhibit pathways like HIF-1, TNF, and MAPK, our findings showed upregulation, warranting further investigation. The study also suggests the need to consider factors like the interaction of TQ with serum proteins and the potential influence of albumin in FBS on our results, indicating the value of conducting future experiments in serum-free media.