The multi-target mechanism of action of Selaginella doederleinii Hieron in the treatment of nasopharyngeal carcinoma: a network pharmacology and multi-omics analysis

Abstract

Nasopharyngeal carcinoma (NPC) presents significant treatment challenges due to its complex etiology and late-stage diagnosis. The traditional Chinese medicine Selaginella doederleinii Hieron (S. doederleinii) has shown potentiality in NPC treatment due to its multi-target, multi-pathway anti-cancer mechanisms. First, we identified NPC related target genes from databases like GeneCards, OMIM, and DisGeNET, and performed WGCNA analysis on the GSE53819 dataset to identify several important gene modules related to NPC. Active components and their targets in S. doederleinii were screened from the TCMSP and other databases, identifying 32 overlapping genes. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis revealed that these genes are primarily involved in critical biological processes like protein phosphorylation and cell cycle regulation. A protein-protein interaction network was constructed, and cytoHubba identified six key genes (BCL2, MAPK14, ABCB1, PLK1, ATM, HMOX1). Kaplan-Meier analysis and immune infiltration analysis further showed that these key genes are closely related to the prognosis and immune microenvironment of NPC patients. Single-cell RNA sequencing analysis revealed the expression distribution of key genes across different immune cell types and explored their roles in the differentiation process of malignant cells through pseudotime trajectory analysis. Molecular docking and dynamics simulation results indicated that the Berberine-MAPK14 and Matairesinol-PLK1 complexes have high binding affinity and stability. Binding free energy calculations confirmed the stability of these complexes. Based on our comprehensive multi-level analysis, the active components of S. doederleinii may play a significant role in the treatment of NPC through multi-pathway and multi-target synergistic effects.

Keywords: Selaginella doederleinii Hieron; MAPK14; Molecular dynamics simulation; Nasopharyngeal carcinoma; Network pharmacology; PLK1.

Fig. 1

Weighted gene co-expression network analysis. (A) Selection of soft thresholds for network topology analysis. The left panel shows the scale-free fit index (R²) under different soft thresholds, and the right panel shows the mean connectivity under different soft thresholds. (B) Clustering dendrogram of gene co-expression modules. The color bars represent different gene modules. (C) Heatmap of module-trait relationships. Colors indicate the strength and direction of correlations, with red representing positive correlations and blue representing negative correlations.

Fig. 2

Elucidation of core targets in S. doederleinii treatment of NPC. (A) Venn diagram of genes related to disease, prognosis, and S. doederleinii. (B) Maximum connectivity subnetwork of intersecting genes. Node colors indicate gene importance, with red being high and yellow being low. (C) UpSet plot showing the overlap of the top 6 targets from six algorithms. (D) Correlation heatmap of key genes. Colors indicate the strength and direction of correlations, with red representing positive correlations and blue representing negative correlations. Numbers indicate correlation coefficients.

Fig. 3

GO enrichment and KEGG pathway analysis. (A) BP. (B) CC. (C) MF. (D) KEGG pathway enrichment analysis. The size of the dots represents the number of genes, and the color indicates the -log10 (Q-value).

Fig. 4

Kaplan–Meier survival analysis. The impact of high and low expression of key genes BCL2 (A), PLK1 (B), ATM (C), HMOX1 (D), MAPK14 (E), and ABCB1 (F) on survival rates. The shaded areas represent the 95% confidence intervals, and the table at the bottom shows the number of patients at risk at each time point.

Fig. 5

Immune infiltration analysis. (A) Comparison of immune cell infiltration levels between tumor and normal tissues. (B) Correlation heatmap between key genes and immune cell infiltration. Colors indicate the strength and direction of correlations, with red representing positive correlations and blue representing negative correlations.

Fig. 6

Single-Cell RNA Sequencing Analysis. (A) t-SNE plot showing the clustering of different cell types. (B) Expression of cell markers in different cell types. (C) Expression distribution of PLK1 across different cell types. (D) Expression levels of key genes in different cell types. The size of the dots represents the percentage of cells expressing the gene, and the color indicates the average expression level.

Fig. 7

Molecular docking results showing the binding modes of four compounds with four key gene proteins. (A) Docking model of ABCB1 with Gallic acid. (B) Docking model of BCL2 with Nobiletin. (C) Docking model of MAPK14 with Berberine. (D) Docking model of PLK1 with Matairesinol.

Fig. 8

Molecular dynamics simulation results demonstrating the stability and conformational changes of ligand–protein complexes over a 100 ns simulation. (A) RMSD curve showing the stability of the Berberine-MAPK14 complex. (B) RMSD curve showing the stability of the Matairesinol-PLK1 complex. (C) RMSF curve for the Berberine-MAPK14 complex. (D) RMSF curve for the Matairesinol-PLK1 complex. (E) Number of hydrogen bonds in the Berberine-MAPK14 complex. (F) Number of hydrogen bonds in the Matairesinol-PLK1 complex.

Fig. 9

Rg and SASA Analysis. (A) Changes in Rg of the Berberine-MAPK14 complex. (B) Changes in Rg of the Matairesinol-PLK1 complex. (C) Changes in SASA of the Berberine-MAPK14 complex. (D) Changes in SASA of the Matairesinol-PLK1 complex.