This article is mentioned in:

Introduction

Uterine cervical cancer (CC) is the fourth most prevalent cancer type among women globally, underscoring the critical need for ongoing advancements in its management (1). Despite progress in the treatment of CC, the existing therapeutic options have persistent limitations, including side effects and the emergence of drug resistance (2). The histopathological classification of CC identifies two main subtypes: Squamous cell carcinoma (SCC) and adenocarcinoma (AC). While the treatment approach for both subtypes is generally similar, AC is associated with a poorer prognosis than SCC (3). Therefore, addressing the challenges posed by AC remains an essential area for improving overall therapeutic outcomes in CC.

Microbiota are considered to serve a pivotal role in carcinogenesis (4,5). In the context of CC, beyond the well-established risk factor of human papillomavirus, studies have highlighted the contribution of the microbiome to CC development and progression (6,7). It is widely acknowledged that the depletion of Lactobacillus can lead to dysbiosis, characterized by an increase in pathogenic microbial diversity, such as Atopobium, Gardnerella vaginalis, Prevotella, Sneathia, Streptococcus and Megasphaera (7–9). This microbial imbalance has been associated with damage to the mucus and cytoskeleton structures, alteration of antimicrobial peptides, promotion of proinflammatory cytokine production, and initiation of cervical intraepithelial neoplasia (CIN) (8,10).

Streptococcus, a prominent genus of gram-positive bacteria, has garnered attention in the context of patients with CC (11). Mulato-Briones et al (12) identified Streptococcus as a predominant bacterium in patients with CC through the culture of cervical exocervix samples. Streptococcus is also a potential biomarker for distinguishing between invasive CC and CIN (13). Notably, Streptococcus agalactiae [group B Streptococcus (GBS)] often colonizes asymptomatically in the female genital tract and can subsequently become a pathogen, leading to severe infections in neonates and adults under certain conditions (14). Previous studies have indicated that GBS stimulates proinflammatory cytokines, including IL-8, IL-1b, IL-6, IL-17 and TNF-α, during its interaction with host epithelial cells (15–18). GBS exhibits a notable capacity to invade cervical epithelial cells and maintain colonization, particularly when compared with vaginal epithelial cells (17). Zhang et al (19) revealed that GBS exerts an indirect effect by mediating human papillomavirus infection in patients with CIN. Additionally, GBS has a higher abundance in CIN2 than in CIN1, and has been suggested as a crucial biomarker for predicting the severity of CIN (20). However, to the best of our knowledge, the underlying mechanism of the association between GBS and CC progression, especially AC progression, remains unclear.

Transcriptomics, through differential gene expression and gene ontology analysis, provides a valuable tool for evaluating mechanistic changes during the interaction between host cells and bacteria (21–23). Molecular alterations in cancer-related genes, coupled with the corresponding changes in signaling pathways, guide the development of novel precision medicine treatments for cancer (24). In alignment with this, leveraging transcriptomic analyses also facilitates the exploration of drug repurposing strategies, allowing the identification of existing drugs with potential anticancer properties (25).

In the present study, transcriptomic analyses were performed to elucidate the gene and pathway transcriptome profile of the host CC after GBS infection, particularly for the AC type. Building upon observed RNA sequencing (RNA-seq) data changes, the present approach aimed to connect the findings from the transcriptomic analysis with the potential for drug repurposing. By identifying specific gene expression patterns and altered pathways, this strategy provides a targeted and efficient route for identifying existing drugs that could be repurposed as anticancer therapies, paving the way for further investigation and clinical development.

Materials and methods

Public data acquisition

Expression and clinical data of 304 samples from The Cancer Genome Atlas (TCGA)-CESC were acquired from TCGA Genomic Data Commons portal (https://portal.gdc.cancer.gov/). Metagenomic profiles of the same cohort were sourced from TCGA metagenomic microbiome study (https://ftp.microbio.me/pub/cancer_microbiome_analysis/TCGA). The detailed pipeline for obtaining metagenomic profiles was outlined in a previous study (26). Relative abundance normalization was applied to the metagenomic data, and Streptococcus absence was defined as a relative abundance equal to 0.

DESeq2 was employed to investigate differentially expressed genes (DEGs) between Streptococcus-present and -absent groups. Significance criteria included an adjusted P<0.05 and an absolute log2-fold change ≥0.5 (27).

Weighted gene coexpression network analysis (WGCNA) is a computational method that organizes genes into clusters or modules based on their coexpression patterns across different samples. This approach reveals the complex relationships between genes, providing insights into their roles in biological processes and their connections to disease phenotypes (28,29). In the present study, WGCNA was applied to identify significant genes associated with Streptococcus in samples from patients with AC. A soft threshold of β=7 and a scale-free topology fitting index (R2) of 0.95 were used for matrix transformation according to the scale-free topological criteria. The modules exhibiting the highest correlation with Streptococcus presence were selected. Key genes related to Streptococcus in the AC samples were then identified if they met the following criteria: Gene significance (GS) >0.5, module membership (MM) >0.8 and P<0.05.

Pathway enrichment analysis was conducted using two web-based tools, including g:Profiler (https://biit.cs.ut.ee/gprofiler/gost; version 2023) (30) and Enrichr (https://maayanlab.cloud/Enrichr/; version June 8, 2023) (31), employing the Kyoto Encyclopedia of Genes and Genomes pathway database (32) and Molecular Signatures Database Hallmark 2020 (33). Statistical significance was established with an adjusted P<0.1 and a minimum of four genes per pathway. Significant enriched pathways were visualized using ggplot2 version 3.4.2 (34).

HeLa cell culture

HeLa cells, also known as human AC cell lines (KCLB number 10002; lot no. 59726; passage no. 98; Korean Cell Line Bank; Korean Cell Line Research Foundation), were thawed, cultured and passaged within T75 cell culture flasks using Minimum Essential Medium (Gibco; Thermo Fisher Scientific, Inc.) supplemented with 10% fetal bovine serum (Gibco; Thermo Fisher Scientific, Inc.) in a humidified atmosphere with 5% CO2 at 37°C.

Two subculturing instances were conducted to ensure cellular stability for experimentation. Initially, cells were seeded at 5–6×105 cells per flask. Media exchange occurred every 2 days or as needed based on cellular conditions, with harvesting conducted once cells achieved 90–95% confluence.

Streptococcus agalactiae exposure experiment using HeLa cells

GBS or Streptococcus agalactiae strain NCTC 818 [G19] (13813; American Type Culture Collection) was aerobically grown at 37°C in Tryptic Soy Broth (Becton, Dickinson and Company) according to the manufacturer's protocol.

HeLa cells were infected with GBS for 0, 2, 6 and 24 h at a multiplicity of infection equal to 100 (GBS:cell=100:1). Right before the experiment, GBS in the bacterial medium was precipitated and resuspended in human medium (Minimum Essential Medium supplemented with 10% fetal bovine serum), ensuring the concentration remained consistent with the 100:1 ratio. For the control group, human medium without bacteria was administered. Following this, the flasks were incubated for 0 (baseline), 2, 6 and 24 h in a humidified atmosphere with 5% CO2 at 37°C.

All experiments were conducted in triplicate.

RNA isolation and sequencing

Total RNA was isolated from HeLa cells cocultured with GBS using the RNeasy Mini Kit (74104; Qiagen GmbH) according to the manufacturer's protocol. Ribosomal RNA depletion was performed using the MGIEasy rRNA Depletion Kit (1000005953; MGI Tech Co., Ltd) and library preparation was performed using the MGIEasy RNA Directional Library Prep Set (1000006386; MGI Tech Co., Ltd.). The final library concentration was 113 nM. Sample quantification was performed using the Qubit ssDNA Assay Kit (Q10212; Invitrogen; Thermo Fisher Scientific, Inc.) and the Invitrogen Qubit Fluorometer (Q33216; Thermo Fisher Scientific, Inc.). Sample quality was assessed using an Agilent 2100 Bioanalyzer (G2939AA; Agilent Technologies, Inc.). Sequencing was conducted on the DNBSEQ-G400 sequencer (MGI Tech Co., Ltd.) using the DNBSEQ-G400RS-High throughput sequencing FCL PE100 kit (1000016949; MGI Tech Co., Ltd.), generating 100-bp paired-end reads.

FASTQ files from the sequencer underwent quality assessment using FASTQC v0.12.1 and Multiqc version 1.24.1 (35,36). Cutadapt version 4.9 was employed to eliminate low-quality reads and sequencing adapters (37). The alignment and quantification of RNA-seq data were performed using STAR version 2.7.11b (38). Detailed information on sequenced and mapped samples can be found in Table SI.

Time-series transcriptomic analysis of GBS-HeLa cell interactions

TimeSeriesAnalysis (TiSA) is a tool including analysis and visualization packages for RNA-seq and microarrays. It facilitates the extraction of significant genes from time series transcriptomic data by assessing differential gene expression along both condition and temporal axes (39).

In the present study, TiSA was used to identify significant genes meeting the criteria of an absolute log2-fold change ≥1 and an adjusted P<0.05 from longitudinal transcriptomic data, elucidating the interaction between HeLa cells and GBS.

A principal component analysis (PCA) plot was generated using the first and second principal components to visualize the differences between control and GBS-exposed HeLa cell groups at various time points. This was conducted using the plot_PCA_TS() function from the TiSA package. The partitioning algorithm based on recursive thresholding (PART) from the clusterGenomics R package, version 1.0 (40), was applied for clustering. DEGs and clusters were visualized using heatmaps, created with the ComplexHeatmap package version 2.8.0 (http://www.bioconductor.org/packages/devel/bioc/html/ComplexHeatmap.html). To elucidate the biological meaning of the clusters, Enrichr and g:Profiler were used, with the criteria as aforementioned.

Reverse transcription-quantitative PCR (RT-qPCR)

Total RNA was extracted from GBS-exposed HeLa cells after 24 h to validate the mRNA expression of the genes belonging to the MAPK pathway [fibroblast growth factor 21 (FGF21), nerve growth factor (NGF), IL1A and IL1B] and mTORC1 pathway [DNA damage inducible transcript 3 (DDIT3), cystathionine γ-lyase (CTH), asparagine synthetase (glutamine-hydrolyzing) (ASNS) and nuclear protein 1 (NUPR1) using the RNeasy Mini Kit (74104; Qiagen GmbH). RT-qPCR was performed using a one-step TOPreal™ SYBR Green RT-qPCR Kit (cat. no. RT432S; Enzynomics Co., Ltd.) on a Thermal Cycler Dice Real Time System TP800 (Takara Bio, Inc.). This one-step kit combines high-yield TOPscript reverse transcriptase (cat. no. RT002; Enzynomics Co., Ltd.) with chemically modified Taq polymerases. According to the manufacturer's instructions, reverse transcription was performed at 50°C for 30 min, and the thermocycling conditions were as follows: Initial denaturation at 95°C for 10 min, followed by 45 cycles of denaturation at 95°C for 5 sec, and annealing and elongation at 60°C for 30 sec. The primer sequences utilized are presented in Table SII.

All reactions were performed in triplicate. The relative expression levels of the target genes were normalized to the reference gene (GAPDH) using the ΔCq method. Fold changes in gene expression were calculated using the 2−ΔΔCq method, with the control sample serving as the calibrator (41).

Western blotting

To validate the activation of significantly enriched pathways in HeLa cells after 24 h of exposure to GBS, western blotting was performed. HeLa cells were co-cultured with GBS for 24 h, washed with cold PBS and lysed using cold M-PER™ Mammalian Protein Extraction Reagent containing 0.2% Halt™ Protease Inhibitor Cocktail (Thermo Fisher Scientific, Inc.). Tumor cell lysates (20 µg/lane; protein concentration quantified using a Bradford assay) were loaded onto a 4–20% Mini-PROTEAN TGX precast gel (Bio-Rad Laboratories, Inc.) and electrophoresed for 10 min at 80 V, followed by 60–90 min at 120 V. The proteins were transferred to PVDF membranes for 45 min at 12 V using the semi-dry blotting system (ATTO Corporation). After blocking with 5% skim milk for 1 h at room temperature, the membranes were incubated overnight at 4°C with the following primary antibodies: ERK1/2, phosphorylated ERK1/2 (p-ERK1/2), phosphorylated mTOR at Ser2448 (p-mTOR), p38, phosphorylated p38 (p-p38), panAKT (each at a dilution of 1:1,000), mTOR (dilution, 1:2,000), phosphorylated AKT at Ser473 (p-AKT) and Ki67 (each at a dilution of 1:5,000). Detailed information about the antibodies is provided in Table SIII.

Following incubation, the membranes were treated for 1 h at room temperature with anti-rabbit IgG HRP-linked antibody and anti-mouse IgG HRP for β-actin antibody at a dilution of 1:5,000. The membranes were washed three times with 1X TBS with 0.1% Tween-20 for 10 min. The PVDF membranes were then processed with a PicoEPD Western Blot Detection kit (EBP1073; ELPIS-Biotech, Inc.). The chemiluminescence images were captured using the Amersham Imager 600 (GE Healthcare), and the resulting bands were semi-quantified to determine relative protein levels using ImageJ software version 1.54 (National Institutes of Health). The experiments were conducted in triplicate for both control and treatment groups.

Therapeutic agents

Candidate therapeutic agents were selected through a drug repurposing approach. One valuable tool in drug repurposing is the Connectivity Map (CMap) tool (42). CMap is a powerful platform for unraveling connections among drugs, genes and diseases. The fundamental principle involves comparing the gene expression profiles induced by a specific drug with those associated with particular diseases (42).

Genes from TiSA clusters of interest were uploaded to the CMap web-based tool version 1.1.1.43 (https://clue.io/). Potential anticancer drugs were identified based on their tau score, which is converted from normalized connectivity scores, comparing them with values of that disease for all the drugs in the reference database (43,44). Tau scores exceeding 90 or below-90 between two signatures indicate strong connectivity. A negative connectivity signifies the reversal of the disease signature by the drug (43).

Binimetinib (MEK162) was purchased from Selleck Chemicals (cat. no. S7007). A 10 mM stock solution was generated in DMSO and stored at −20°C.

Ridaforolimus (deforolimus; cat. no. HY-50908; MedChemExpress) was prepared as 10 mM stock solutions with DMSO and stored at −20°C.

Specially tailored working solutions of these drugs were prepared for immediate use. In the case of monotherapy, binimetinib and ridaforolimus were serially diluted to concentrations of 0.625, 1.25, 2.5, 5, 10, 15 and 20 µM. For combination therapy of binimetinib and ridaforolimus, a 1:1 ratio of the drugs was applied, with serial dilutions of 1.25, 2.5, 5, 10, 15 and 20 µM used. For all treatments the cells were incubated at 37°C with 5% CO2. Each experiment was conducted in triplicate.

Trypan blue assay

A trypan blue assay was utilized to evaluate the viability of HeLa cells following exposure to GBS or therapeutic agents. After washing with PBS, adherent cells in a 90-mm plate were dissociated using 2 ml of 0.25% trypsin-EDTA for 2 min at room temperature. The cells were then resuspended in 8 ml cell culture medium. A portion of this cell suspension was mixed with an equal volume of 0.4% trypan blue (Gibco; Thermo Fisher Scientific, Inc.) for 2 min at room temperature. The number of viable and non-viable cells was counted in triplicate using a hemocytometer.

MTS assay

After reaching confluence in a T75 flask, HeLa cells were subcultured into 96-well plates at a concentration of 2×104 cells/100 µl. Following a 24-h incubation period, binimetinib and/or ridaforolimus were added into each well using serial dilutions of 1.25, 2.5, 5, 10, 15 and 20 µM.

Subsequently, the cells were incubated at 37°C with 5% CO2 for either 24 or 48 h. Afterwards, 20 µl CellTiter 96® AQueous One Solution (G3582; Promega Corporation) reagent was added to each well. The plate underwent an additional 4-h incubation at 37°C with 5% CO2, and the absorbance at 490 nm was measured using a microplate reader (Multiskan™ FC; Thermo Fisher Scientific, Inc.).

The absorbance intensity at 490 nm for both the control and treatment groups was normalized to the absorbance of the background. Subsequently, cell viability percentages for different treatment groups were calculated using the following formula: % viability=adjusted optical density (OD)sample/mean (adjusted ODcontrol), where adjusted ODsample=ODsample-mean (ODbackground), and adjusted ODcontrol=ODcontrol-mean (ODbackground).

Cell apoptosis assay

To measure the induction of apoptosis, 5×105 HeLa cells were plated in 90-mm plates and allowed to adhere and proliferate for 48 h. Cells were then treated with either 15 µM binimetinib, 20 µM ridaforolimus, a combination of both (5 µM binimetinib + 5 µM ridaforolimus) or DMSO (vehicle treatment) as a control at 37°C with 5% CO2.

After 48 h, cells were harvested and transferred into tubes at a concentration of 5×105 cells/ml. Cell death was assessed using an annexin V-FITC apoptosis detection kit (ab14085; Abcam) with 5 µl annexin V-FITC and 5 µl PI, followed by a 5-min incubation at room temperature in the dark. The apoptosis assay data were obtained using a BD FACSCalibur flow cytometer (BD Biosciences). The graph showing the percentages of apoptotic and necrotic cells was created using the Floreada tool (https://floreada.io; WASM version SIMD).

Dose-response curve and synergy scores

To estimate the efficacy of anticancer drugs, the dose-response curve R package version 3.0.1 was used to determine the IC50 (45). The IC50 values represent the concentrations of each drug required to achieve a 50% inhibition of cancer cell viability.

To predict the combination (binimetinib and ridaforolimus) effect, the synergy score was calculated using the SynergyFinder R package version 3.8.2 (46). If the synergy score is between −10 and 10, the effect of the combination of the two drugs is additive, if the score is >10, the effect is synergistic, and if the score is <-10, the combination is considered antagonistic (46,47).

Effect of therapeutic agents on GBS-exposed HeLa cells

To evaluate the effect of therapeutic agents on GBS-exposed HeLa cells and to compare these effects with those on HeLa cells not exposed to GBS, GBS and HeLa cells were co-cultured at a ratio of 100:1 in a humidified atmosphere with 5% CO2 at 37°C for 6 h. Following this exposure, the GBS-exposed HeLa cells were treated with 15 µM binimetinib, 20 µM ridaforolimus, a combination of both (5 µM binimetinib + 5 µM ridaforolimus) or DMSO as a control at 37°C with 5% CO2. The number of dead and viable cells after 24 h of treatment was measured using a trypan blue assay in triplicate as aforementioned.

Statistical analysis

To compare two groups, an unpaired Student's t-test was performed. For analyses involving multiple groups, one-way ANOVA followed by post hoc testing using the Holm-Bonferroni correction method was performed (48). All statistical analyses were performed using the stats package version 4.3.3 (https://search.r-project.org/R/refmans/stats/html/00Index.html) in R version 4.3.3 (https://www.R-project.org/). The results were visualized using the ggplot2 package version 3.4.2 (34) and ggpubr package version 0.6.0 (https://github.com/kassambara/ggpubr). Analyses were based on three experimental repeats. In the plots, error bars were included to represent the dispersion of the variable (standard deviation). P<0.05 was considered to indicate a statistically significant difference.

Results

Gene expression and pathway profiles of Streptococcus-associated CC using TCGA-CESC public data

To investigate the relationship between Streptococcus and CC, a cohort of 304 TCGA CC cases, including transcriptomic data and Streptococcus abundance profiles, were analyzed (Table SIV). As shown in Fig. 1A, 31 cases (10.2%) exhibited the presence of Streptococcus.

Figure 1. Gene expression and pathway profiles of Streptococcus-associated CC using The Cancer Genome Atlas-CESC public data. (A) Distribution of Streptococcus-present and -absent cases. (B) DEG profiles between Streptococcus-present and -absent CC groups, with red and blue indicating upregulated and downregulated genes, and grey representing non-significantly DEGs. (C) Pathway profiles of upregulated DEGs. CC, cervical cancer; DEG, differentially expressed gene; padj, adjusted P-value.

DEG profiles (Fig. 1B) revealed 128 upregulated and 209 downregulated DEGs identified using DESeq2 (Table SV). Pathway enrichment analysis, conducted using the web-based tools g:Profiler and Enrichr, of the list of 128 upregulated DEGs, highlighted statistically and biologically significant pathways associated with cancer, including ‘MAPK signaling pathway’, ‘mTORC1 signaling’, ‘Ras signaling pathway’ and ‘Wnt signaling pathway’, as well as metabolism pathways such as ‘metabolism of xenobiotics by cytochrome P450’, ‘fatty acid metabolism’, ‘glutathione metabolism’ and ‘arachidonic acid metabolism’ (Fig. 1C; Table SVI).

To further explore the specificity of the relationship between Streptococcus and AC, 31 patients with AC from a cohort of 304 patients with CC in TCGA were analyzed. Among them, 3 patients with AC exhibited the presence of Streptococcus (Table SIV). Using the WGCNA tool, two significant key modules with high correlation values, r=0.77 (P=5×10−7) for the black module and r=0.59 (P=4×10−4) for the brown module, were identified (Fig. 2A-C; Tables SVII and SVIII).

Figure 2. Streptococcus-related transcriptomic and pathway profiles in 31 patients with cervical adenocarcinoma from The Cancer Genome Atlas analyzed using weighted gene co-expression network analysis. (A) Heatmap of the correlation between Streptococcus abundance and module eigengenes. Each row corresponds to a module, and the color gradient indicates correlation coefficients, with blue representing negative correlations and red representing positive correlations. Corresponding P-values are displayed in brackets. (B) Scatter plot of the correlation between MM of the brown module and GS for Streptococcus abundance. (C) Scatter plot of correlation between MM of the black module and GS for Streptococcus abundance. (D) Bar plot showing results of pathway enrichment analysis for significant genes from the brown and black modules with GS >0.5, MM >0.8 and P<0.05. AMPK, AMP-activated protein kinase; GS, gene significance; ME, module eigengenes; MM, module membership; PPAR, peroxisome proliferator-activated receptor.

A total of 500 key genes related to Streptococcus were selected from the two key modules based on the following criteria: GS >0.5, MM >0.8 and P<0.05. Pathway enrichment analysis showed that these genes were involved multiple oncogenic pathways, including ‘MAPK signaling pathway’, ‘mTOR signaling pathway’, ‘PI3K-AKT signaling pathway’, ‘TGF-beta signaling pathway’ and ‘Ras signaling pathway’ (Fig. 2D; Table SIX).

Time-series transcriptomic analysis of GBS-HeLa cell interaction

To further enhance the insights gained from TCGA-CESC data, a time-series experiment was performed to investigate HeLa cell-GBS interactions. Briefly, HeLa cells were infected with GBS for 0, 2, 6 and 24 h, with corresponding control groups treated solely with bacteria-free human medium. Subsequently, the transcriptomic data generated from this experiment were analyzed using the TiSA package.

PCA plot, and DEGs between conditional and temporal groups

The time series PCA plot demonstrated a clear separation between the GBS-exposed HeLa cell groups at 6 and 24 h and the other groups. PC1, representing variability over time, and PC2, representing differences between groups, were utilized for this analysis. As shown in Fig. 3, the dominance of PC1, accounting for 45% of the variation, while PC2 accounted for 14%, indicated that the variation over time was more pronounced than the variation between groups, which was especially evident after 6 h. Notably, the groups associated with GBS stimulation at 6 and 24 h exhibited a distinct separation from the other groups.

Figure 3. PCA plot. The PCA plot illustrates each sample from two groups: The GBS group and the control group. In the plot, HeLa cells exposed to GBS are represented in red for the time points 0, 2, 6 and 24 h. The control group is denoted by the blue color. Additionally, each time point is distinguished by specific shapes: Circles for 0 h, triangles for 2 h, crosses for 6 h and multiplication signs for 24 h. GBS, group B Streptococcus; PCA, principal component analysis.

The conditional differential gene expression analysis for four-time points is presented in Fig. 4. The number of DEGs markedly increased after 6 h of exposure. Specifically, there were 16, 11, 460 and 1,494 DEGs after 0, 2, 6 and 24 h of exposure, respectively.

Figure 4. Heatmap of conditional and temporal differentially expressed genes. The rows indicate the condition, and the columns indicate the genes. The conditions are divided into exposure to GBS in purple and the control group in orange. The genes are grouped according to the time points: 0 h in red, 2 h in green, 6 h in blue and 24 h in cyan. FC, fold change; GBS, group B Streptococcus; TP, timepoint.

Clustering analysis

All DEGs identified from the conditional differential gene expression analysis, totaling 1,981 genes, were utilized for PART clustering, resulting in six clusters, as shown in Fig. 5A and Table SX. Genes in clusters 1, 2, 3 and 4 showed an upward trajectory between 6 and 24 h, while genes in cluster 6 decreased during the same period. Cluster 5 exhibited mixed expression patterns (Fig. 5B).

Figure 5. Clustering results. (A) Heatmap of partitioning algorithm based on recursive thresholding clustering including 6 clusters with different colors, following exposure of HeLa cells to GBS across different time points. Condition groups are represented at the top, with the control group indicated in blue and the GBS group indicated in red. The yellow and orange bar colors represent the time points below the condition groups. (B) Gene expression trajectory plot for six clusters identified from time series analysis comparing between GBS exposure and control groups. Each plot shows the average expression trend for all genes in each cluster across the time points, with the x-axis representing time points (h after exposure) and the y-axis showing normalized gene expression levels. Red curves represent the gene expression trajectories in GBS-exposed HeLa cells, and blue curves represent the control group. GBS, group B Streptococcus.

Pathway enrichment analysis

Functional analysis of each cluster was performed using Enrichr and g:Profiler. Cluster 2 exhibited the most significant pathways, encompassing both oncogenic pathways (‘MAPK signaling pathway’, ‘mTORC1 signaling’ and ‘p53 pathway’) and immune response pathways (‘TNF signaling pathway’, ‘NF-kappa B signaling pathway’, ‘IL-17 signaling pathway’ and ‘NOD-like receptor signaling pathway’) (Table SXI).

The pathway enrichment analyses revealed the upregulation of ‘MAPK signaling pathway’ and ‘mTORC1 signaling’ in the HeLa cells-GBS exposure experiment (Fig. 6) and the TCGA-CESC public data (Fig. 1C). Additionally, activation of the ‘MAPK signaling pathway’ and ‘mTOR signaling pathway’ was observed in the group of 3 patients with AC with Streptococcus presence (Fig. 2D). These findings suggested that the MAPK and mTORC1 pathways may serve significant roles in AC associated with GBS.

Figure 6. Enriched pathways of cluster 2. The ‘MAPK signaling pathway’ and ‘mTORC1 signaling’ (shown in red) were identified in both the time-series analysis and the public data analysis. NOD, nucleotide-binding oligomerization domain.

Validation of MAPK and mTORC1 signaling pathway activation after GBS exposure

To validate the activation of the MAPK and mTORC1 pathways, the upregulated expression of genes associated with the MAPK pathway (FGF21, NGF, IL1A and IL1B) and the mTORC1 pathway (DDIT3, CTH, ASNS and NUPR1) from cluster 2 was first confirmed using RT-qPCR. As shown in Fig. 7C, most of these genes exhibited higher expression levels in HeLa cells exposed to GBS for 24 h compared with control cells, with the exception of NGF. This discrepancy may arise from differences in sensitivity or specificity between RNA-seq and RT-qPCR methods. Overall, the RT-qPCR results support the activation of the MAPK and mTORC1 pathways in response to GBS exposure.

Figure 7. Effects of GBS on HeLa cells. (A) Line graph showing HeLa cell proliferation after 6 and 24 h of GBS exposure, with comparisons with control cells and between two time points, as measured by the trypan blue assay. (B) Western blot analysis and corresponding bar plot showing Ki67 expression in HeLa cells at 6 and 24 h post-GBS exposure. (C) Bar plot showing the logarithm fold changes in the expression of eight targeted genes as determined by reverse transcription-quantitative PCR. (D) Western blot analysis and corresponding bar plots showing the levels of proteins in HeLa cells 24 h after GBS exposure. Blue, orange and red indicate the control group, the 6 and 24-h GBS-treated groups, respectively. *P<0.05, **P<0.01, ***P<0.001. GBS, group B Streptococcus; ns, not significant; p-, phosphorylated; t-, total.

Western blot analysis revealed a significant increase in the ratios of p-ERK1/2 to total (t-)ERK1/2, p-AKT to t-AKT and p-mTOR to t-mTOR in HeLa cells exposed to GBS for 24 h compared with control cells. Although the ratio of p-p38 to t-p38 was higher in the treated HeLa cells than in the control group, the difference was not statistically significant (Fig. 7D). The elevated phosphorylation of ERK1/2 indicated enhanced MAPK pathway activity, while the upregulation of p-AKT and p-mTOR suggested activation of the mTORC1 pathway, further supporting the involvement of these signaling cascades in the cellular response to GBS exposure.

Additionally, the trypan blue assay demonstrated a decrease in viable HeLa cell numbers at 6 h compared with control cells, followed by a significant increase in HeLa cell numbers, reaching levels similar to control cells by 24 h (Fig. 7A). This cell proliferation phenomenon was confirmed by examining the expression of Ki67, a well-known cell proliferation marker. The highest Ki67 expression was observed after 6 h of exposure and the expression was still increased at 24 h compared with that in the control group (Fig. 7B). This pattern may reflect cancer cell proliferation influenced by MAPK and mTORC1 pathway activity after 6 h of exposure to GBS.

Overall, transcriptomic analysis using public data, time-series co-culture experiments, and validation through RT-qPCR and western blotting assays underscored the pivotal role of the MAPK and mTORC1 signaling pathways. These promising results led to the hypothesis that inhibiting these pathways could be effective in treating AC.

Therapeutic agent selection

To further investigate candidate drugs targeting the MAPK and mTORC1 signaling pathways, the list of genes from cluster 2 was uploaded to the CMap web-based tool. Ridaforolimus (or deforolimus) was chosen as the mTOR inhibitor due to its highest negative range score (Table SXII). Notably, to the best of our knowledge, it has not been extensively studied in the context of AC or HeLa cells.

While CMap suggested several anticancer drugs to inhibit the MAPK pathway, these drugs were previously investigated in HeLa cells (49–51). Therefore, binimetinib, a MEK inhibitor that has shown efficacy in other cancer types (52,53) but has not yet been explored in AC, was selected for further evaluation in the present study.

Evaluating the effects of therapeutic agents on HeLa cells

Although the viability of HeLa cells only notably decreased after 24 h of exposure to doses >20 µM compared with the control group, a significant reduction in cell viability was observed across all binimetinib concentrations (1.25–20 µM) following 48 h of treatment (Figs. 8A and S1). The IC50 was determined to be 14.40±2.55 µM (Fig. 8A).

Figure 8. MTS assay results for monotherapies in HeLa cells. (A) Binimetinib and (B) ridaforolimus. Blue and dark red represent drug exposure durations of 24 and 48 h, respectively.

For the apoptosis assay, binimetinib was administered at a concentration of 15 µM, corresponding to the IC50 value observed in the MTS assay. The apoptotic effect after 48 h of binimetinib treatment was higher than that observed in the control group; nevertheless, no significant difference was detected between the two groups (Figs. 9 and S2B).

Figure 9. Apoptosis assay results for therapeutic agents in HeLa cells. Green, orange, red and blue represent the binimetinib, ridaforolimus, combination (binimetinib + ridaforolimus) and control groups, respectively. **P<0.01. ns, not significant.

Ridaforolimus

The MTS assay revealed minimal inhibitory effects of ridaforolimus on HeLa cell proliferation. No significant reduction in cell viability was observed after treatment with ridaforolimus for 24 and 48 h (Figs. 8B and S3). The IC50 of 59.61±10.01 µM further supported this observation (Fig. 8B).

Based on this finding, the highest concentration of ridaforolimus (20 µM) was used for the apoptosis assay. This dose exhibited a weak cytotoxic effect (Figs. 9 and S2C).

Combination effect of binimetinib and ridaforolimus on HeLa cells

A 1:1 ratio of binimetinib and ridaforolimus was administered to assess the impact of combination therapy. SynergyFinder was used to predict potential synergy, and the results suggested that this combination may exhibit an additive effect, as indicated by a Loewe synergy score of ~4.31, with statistical significance at P<0.05 (Fig. 10). Compared with ridaforolimus alone, the combination of binimetinib and ridaforolimus showed significantly greater inhibition of cell viability in the MTS assay. However, significant differences between the combination and binimetinib alone were observed only at concentrations of 10 and 20 µM (Fig. 11).

Figure 10. Contour plot. Loewe synergy score of the binimetinib and ridaforolimus combination in HeLa cells after 48 h.

Figure 11. Comparison of MTS assay results between monotherapies and combination therapy in HeLa cells after 48 h. Green, orange and red represent the cell viability percentages for the binimetinib, ridaforolimus and combination (binimetinib + ridaforolimus) groups, respectively. **P<0.01, ***P<0.001. ns, not significant.

For the apoptosis assay, 5 µM binimetinib and 5 µM ridaforolimus were used, as this combination effectively reduced cell viability by 50% as indicated by the MTS assay. Although the predictive tool suggested an additive effect, the combination therapy induced notable cytotoxic effects, as evidenced by the highest percentage of annexin V observed 48 h after treatment in the apoptosis assay (Figs. 9 and S2D).

Furthermore, the IC50 of this combination, calculated using the dose-response curve package, was determined to be 6.67±0.45 µM for each drug (Fig. S4). In practical terms, this translates to a dose reduction of 2.16 times for binimetinib and 8.94 times for ridaforolimus to reach the IC50 when combined, emphasizing the potent impact of this combination.

Evaluating the effects of the therapeutic agents on HeLa cells with GBS exposure

Based on the promising therapeutic effects observed, we hypothesized that these anticancer agents (binimetinib, ridaforolimus and their combination) might affect HeLa cells exposed to GBS more than the group of cells without exposure to GBS. Experiments were conducted to confirm this hypothesis.

After 6 h of exposure of HeLa cells to GBS, the cells were treated with therapeutic agents (monotherapy and combination therapy). As shown in Fig. S5A, both monotherapies and combination therapies exhibited enhanced effects on HeLa cells exposed to GBS compared with the control group after 24 h of treatment, with the combination therapy demonstrating a more pronounced effect than monotherapy.

In addition, the data specifically indicated that combination therapy was more effective in HeLa cells exposed to GBS than in those not exposed to GBS (Fig. S5B).

Discussion

Numerous studies have suggested that Streptococcus may be involved in cancer initiation, progression and metastasis. For instance, Streptococcus anginosus has been associated with activating gastric tumorigenesis (54). Streptococcus gallolyticus subsp. Gallolyticus is strongly linked to promoting human colon cancer cell proliferation by increasing β-catenin signaling (55). Similarly, Streptococcus pneumoniae is related to upregulated mTOR2/AKT signaling pathways for invasion and migration in lung cancer cells (56). Despite these observations, the understanding of the interaction between Streptococcus and CC remains limited. Therefore, the present study aimed to elucidate the underlying mechanism of this interaction, with a specific focus on the relationship between GBS and AC. The selection of GBS for the present study was based on its prevalence as a Gram-positive coccus, often carried asymptomatically by numerous women (57). Furthermore, it has been reported to be associated with various types of CIN in CC (19,20). Regarding the CC subtype selection, AC was chosen due to its poor prognosis compared with SCC (3). Notably, treatment strategies for AC lack specificity in clinical settings (58,59), adding significance to the present investigation.

The comprehensive transcriptomic analysis using in vivo public data (TCGA-CESC) and in vitro time series data revealed that the ‘MAPK signaling pathway’ and ‘mTORC1 signaling’ emerged as crucial elements in the host-GBS interaction.

Elevated MAPK pathway activation is implicated in cancer occurrence and progression, making it a pertinent target for investigation (60,61). Research also indicates that pathogens often target host signaling pathways to regulate a number of cell processes in cancer, and the most common pathway observed is MAPK signaling (62). Among the MAPK cascades, ERK1/2 is particularly noted for its role in regulating cell proliferation, survival, metabolism, migration and differentiation (63,64), while p38 is more closely associated with stress responses and inflammation (65). Nevertheless, the present study revealed that the interaction between GBS and HeLa cells involved MAPK signaling via ERK1/2, serving a more dominant role than the p38 pathway.

The activation of mTORC1 pathways was confirmed by the upregulation of p-AKT and p-mTOR in HeLa cells exposed to GBS for 24 h, as demonstrated by western blot analysis. The activated p-AKT Ser473 can phosphorylate and inhibit TSC complex subunit 2 (TSC2), a negative regulator of mTORC1. The inhibition of TSC2 leads to the activation of the Rheb GTPase, which in turn directly activates the mTORC1 signaling pathway (66,67). This pathway serves a crucial role in various cellular processes, including cell viability, proliferation and metabolism (68,69). Dysregulation of mTORC1 signaling has been implicated in the development and progression of cancer (66–68).

The present study also validated the expression of genes related to the MAPK pathway (FGF21, NGF, IL1A and IL1B) and the mTORC1 pathway (DDIT3, CTH, ASNS and NUPR1) from cluster 2 using RT-qPCR. Most of these genes were upregulated in HeLa cells exposed to GBS for 24 h compared with control cells at 24 h. Previous studies have shown that these genes are associated with enhancement of MAPK and mTORC1 pathway activation (70–77), as well as cancer cell progression and proliferation (78–84).

To evaluate the association between cell proliferation and GBS, HeLa cells were examined after exposure to GBS for 6–24 h. A significant decrease in HeLa cell numbers was observed after 6 h of GBS exposure, indicating an initial cytotoxic effect. This phenomenon was similar to the results of a study by Tyrrell et al (85), which examined cell death at 6 h of direct contact between HeLa cells and GBS. However, after 6 h of exposure, the number of viable cells increased and reached levels similar to the control cells at 24 h. This phenomenon was confirmed by elevated Ki67 expression, a common cell proliferation marker, present during the active phases of the cell cycle (G1, S, G2 and M) and absent in resting cells (G0) (86,87). Ki67 expression was significantly higher in GBS-treated HeLa cells at 6 and 24 h compared with that in control cells, indicating that GBS treatment stimulated cell proliferation. The peak in Ki67 expression at 6 h may reflect a rapid initiation of cell cycle activity, with numerous cells entering or actively progressing through the cycle. Since Ki67 levels typically peak at G2 or during the mitosis phase (88), the subsequent decrease at 24 h could suggest that fewer cells were entering the cell cycle or that cells that initially proliferated had completed division, resulting in a lower Ki67 signal at this time point.

According to Burnham et al (89), GBS could protect HeLa cells from caspase-3-induced apoptosis during early exposure by triggering the PI3K/AKT signaling pathway. Since mTORC1 signaling is a subsequent cascade of PI3K/AKT, the present results contributed to understanding the complex mechanisms underlying the interaction between HeLa cells and GBS more comprehensively than before.

In general, while GBS may initially cause cell stress or death, it subsequently activates survival and proliferative pathways, specifically the MAPK and mTORC1 pathways, leading to the recovery and proliferation of HeLa cells (Fig. 12).

Figure 12. Schematic diagram of GBS-induced HeLa cell proliferation via MAPK and mTOR pathway activation. ASNS, asparagine synthetase (glutamine-hydrolyzing); CTH, cystathionine γ-lyase; DDIT3, DNA damage inducible transcript 3; FGF21, fibroblast growth factor 21; GBS, group B Streptococcus; NGF, nerve growth factor; NUPR1, nuclear protein 1; p, phosphorylated.

While several signaling pathways, such as the ERK/MAPK, PI3K/AKT, EGFR/VEGF and Wnt signaling pathways, are activated in CC development, evidence for specific pathways in AC is limited (24,90,91). By highlighting the importance of MAPK and mTORC1 signaling pathways in the interplay between AC cells and GBS, the present study provided indirect evidence of their crucial role in AC development mechanisms. Furthermore, building upon this observation, a promising therapeutic approach for AC may be targeting these signaling pathways through a drug repositioning strategy.

Drug repurposing refers to the process of identifying novel uses for existing drugs that were initially developed for a different therapeutic indication (92,93). This can potentially reduce the time and costs associated with drug development. In alignment with this strategy, the present study investigated the MEK inhibitor binimetinib and the mTOR inhibitor ridaforolimus, which were initially developed for the treatment of other cancer types (94–97) but have not been previously studied in the context of AC therapy.

While binimetinib showed a modest cytotoxic effect, it demonstrated a significant cytostatic impact on HeLa cells after 48 h, as indicated by a decrease in the absorbance at 490 nm compared with the control group in the MTS assay. The cytostatic effect refers to its ability to inhibit or slow down the viability and proliferation of cancer cells, providing a means of controlling cancer progression (98). Conversely, ridaforolimus exhibited limited efficacy in achieving the desired IC50 level after 48 h, accompanied by mild cytotoxic effects. These findings prompted further investigation into strategies to enhance the effectiveness of both binimetinib and ridaforolimus.

Evidence suggests the efficacy of mTOR inhibitors, such as rapamycin and its analogs, in cancer treatment (66,99). However, the development of resistance poses a challenge to their long-term effectiveness (100). Targeted inhibition of mTOR can induce MAPK reactivation and lead to resistance to single mTOR inhibition (24,101). To address this challenge, a combinatorial approach targeting mTOR and MAPK signaling should be considered for a more effective response to therapies (102,103). Accordingly, a combination study involving binimetinib and ridaforolimus was conducted. Despite the synergy score suggesting an additive effect, the MTS and apoptosis assays revealed that the combination outperformed using binimetinib or ridaforolimus alone by significantly inhibiting cell proliferation at the concentrations of 10 and 20 µM and exhibiting a robust cytotoxic effect. This outcome underscored the importance of exploring combination therapies, as they may yield more effective results than monotherapy treatments in AC. In particular, the combination (binimetinib and ridaforolimus) was more effective in the group of HeLa cells with GBS exposure than in the group without exposure. This indicated the important role of this combination in treating patients with AC with GBS dominant in their microbiome profile.

A limitation of the present study is the absence of in vivo experiments, which restricts the ability to fully understand the physiological relevance and systemic impact of MAPK and mTORC1 pathway activation in response to GBS exposure. Future in vivo studies are needed to validate the findings and assess their implications in a more complex biological context, as well as to assess the potential efficacy and safety of these drugs for AC management.

In conclusion, the present transcriptomic analysis revealed the complex molecular interaction between GBS and AC, highlighting the crucial involvement of the ‘MAPK signaling pathway’ and ‘mTORC1 signaling’. The combination of two anticancer drugs targeting these pathways, binimetinib and ridaforolimus, might become a potential therapy for the treatment of AC, offering a promising direction for future research and clinical investigation.

Supplementary Material

Supporting Data

Supporting Data

Supporting Data

Supporting Data

Supporting Data

Supporting Data

Supporting Data

Supporting Data

Supporting Data

Supporting Data

Supporting Data

Supporting Data

Acknowledgements

Not applicable.

Funding

This research was supported by the Brain Pool program funded by the Ministry of Science and ICT through the National Research Foundation of Korea (grant no. RS-2023-00283791).

Availability of data and materials

The RNA sequencing data have been deposited in the National Center for Biotechnology Information Gene Expression Omnibus database under accession number GSE264492 or at the following URL: https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE264492. All other data generated in the present study may be requested from the corresponding author.

Authors' contributions

HSH and JHS were responsible for study conception. HSH, JHS, HDTN and TML were responsible for study design. HDTN, TML, YC, DL, DRJ and YJ performed the experiments. HDTN, TML, DL and HSH analyzed the data. HDTN wrote the original draft. IS, JC, NJYP, GOC, YC, OEL, YJ and DRJ contributed to data interpretation and revised the manuscript. HDTN and TML confirm the authenticity of all the raw data. All authors have read and approved the final manuscript.

Ethics approval and consent to participate

Not applicable.

Patient consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

References