The potential mechanism of Huangqin for treatment of systemic lupus erythematosus based on network pharmacology, molecular docking and molecular dynamics simulation

Identification of active components and targets of Huangqin

The Traditional Chinese Medicine Systems Pharmacology Database and Analysis Platform (TCMSP, tcmsp-e.com) was utilized to retrieve information on the physicochemical properties and target action information of the active components of Scutellaria baicalensis (Ru et al., 2014). Criteria of effective ingredients including oral bioavailability (OB) and drug-likeness (DL). Based on the criteria provided in existing literature, the selection criteria were set to OB > 30%, DL > 0.18, and a total of 36 effective components and 512 targets were selected. All target information was standardized through the UniProt database (https://www.uniprot.org/) and a network of active components and target genes was constructed using Cytoscape 3.7.2.

Disease target screening and acquisition

The GeneCards database (https://www.genecards.org/) was used with the search query “[disorders](systemic lupus erythematosus)” to identify related genes, and the KEGG database (https://www.kegg.jp/) was searched for “Systemic lupus erythematosus” to supplement disease genes.

The GSE database (GSE51997) was utilized to obtain gene expression profiles from normal individuals and lupus patients (Kyogoku et al., 2013; Menssen et al., 2009). The dataset includes four healthy individuals and six lupus patients. Data was imported into R 4.3.1 and differential gene analysis was conducted using the limma R package with a set condition of P < 0.05, |logFc| > 1. The ggplot2 and ComplexHeatmap R packages were used to create volcano plots and gene heatmaps.

Construction and analysis of the PPI network

Add the intersection targets of Huangqin and SLE in the STRING database (https://cn.string-db.org/) (Szklarczyk et al., 2021), with “Homo Sapiens” set as the species, and all other parameters set to default, to obtain the Protein-Protein interaction (PPI) network. The network file was imported into Cytoscape v3.7.2 for processing, and the protein-protein interaction network diagram was redrawn. The MCC algorithm in the cytoHubba plugin of Cytoscape was used to identify the top five hub genes (Chin et al., 2014). The MCODE program in Cytoscape was used for protein clustering, and the cluster was further subjected to BP enrichment analysis.

GO enrichment analysis and KEGG enrichment analysis

The Database for Annotation, Visualization and Integrated Discovery (DAVID, https://david.ncifcrf.gov) is used to perform GO enrichment analysis and KEGG enrichment analysis on the key targets. GO analysis was conducted in three aspects: biological process (BP), molecular function (MF), and cellular component (CC). Results with an FDP ≤ 0.05 were selected, and visualization was performed using the “ggplot2” R software package. The top 50 KEGG pathways were categorized and summarized based on the six classifications in the KEGG PATHWAY database. Non-SLE-related pathways were excluded from the KEGG pathway enrichment analysis. A targeted pathway network was constructed using Cytoscape 3.7.2, the targeted areas and color tones in the network were arranged according to the Degree value.

Molecular docking validation

Molecular docking between the core components and the core targets was performed based on the above analysis. The 3D crystal structures of the target proteins were obtained in the Protein Data Bank (PDB) (Drew et al., 1981). Heir PDB format files were downloaded, dehydrated and hydrogenated in AutoDockTools (version 1.5.7), then selected as receptors and saved as pdbqt files. The mol2 files of the active components were obtained from the TCMSP database, hydrogenated in AutoDockTools, selected as ligands, and exported as pdbqt files. Molecular docking was performed using AutoDock Vina (version 1.2.3) (Eberhardt et al., 2021; Trott & Olson, 2010), and the affinity of each component to the target was obtained. Components and targets with good binding activity were selected based on affinity and visualized using Discovery Studio 2019.

Molecular dynamic simulation

Molecular dynamics (MD) simulations were performed using the Gromacs2022 program, with the Generalized Amber Force Field (GAFF) for small molecules, the AMBER14SB force field for proteins, and the TIP3P water model, with the system temperature set at 298 K and the simulation time at 100 nanoseconds. After the simulation was completed, the simulation trajectories were analyzed using VMD and PyMOL, and the g_mmpbsa program was used to perform the MMPBSA binding free energy analysis between the protein and the small molecule ligand.

Effective components and targets of Huangqin

A total of 36 active components and 512 drug targets of Scutellaria baicalensis were identified through the TCMSP database (tcmsp-e.com). The top five components with the highest degree values analyzed by NetworkAnalyst in Cytoscape were: baicalin, acacetin, oroxylin a, neobaicalin, and stigmasterol (Table 1, Table S1). Hexagons represent the active components of Huangqin and rhombuses represent the targets of the active components. The degree value represents the number of edges connected to the node (Fig. 2).

Identification of targets related to Scutellaria baicalensis intervention in SLE

Using the GeneCards database (https://www.genecards.org/), 5,947 related genes were filtered out, and the top 1,493 were selected based on the Score value. The KEGG database (https://www.kegg.jp/) provided 23 disease genes for “Systemic lupus erythematosus.”

Based on the GSE database GSE51997 microarray dataset, which includes four healthy individuals and six lupus patients, a total of 160 up-regulated genes and 35 down-regulated genes were obtained (Table S2). Volcano maps and heat maps of DEGs were shown in Figs. 3A and 3B.

After merging the genes obtained from the three databases and removing duplicates, a total of 1,639 disease-related genes were obtained. A Venn diagram was used to identify the intersecting targets between the predicted Huangqin targets and SLE targets. The analysis identified 28 targets as potential candidates for the anti-SLE action of Huangqin (Fig. 4A).

Construction and analysis of PPI network

A PPI network containing 28 nodes and 334 edges was constructed by adding the intersection targets to the STRING database (https://cn.string-db.org/) with “Homo Sapiens” specified as the species. The network file was imported into Cytoscape v3.7.2 for processing (Fig. 4B, Table S3). The darker the color of the nodes in the diagram, the greater the impact on SLE. Protein clustering of the targets was performed using the MCODE plugin in Cytoscape, and BP enrichment analysis was further conducted on the cluster, with the seed protein of the cluster being CYCS (Fig. 5A), mainly involved in the signal transduction of intracellular steroid hormone receptors, positive regulation of transcription, and response to estradiol (Fig. 5B). The top five hub genes identified by the MCC algorithm in the cytoHubba plugin of Cytoscape are IL6, CASP3, BCL2, ESR1, and FOS (Fig. 5C).

GO and KEGG enrichment analysis

Gene ontology (GO) analysis showed that the main enriched biological processes (BP) targeted by Scutellaria baicalensis in SLE include positive regulation of transcription by RNA polymerase II promoter, positive regulation of transcription, response to estradiol, response to external stimulus, mammary gland alveolar development, and response to hypoxia; the main cellular components (CC) include macromolecular complexes, chromatin, caspase complexes, nucleoplasm, and transcription factor complexes; the main molecular functions enriched include RNA polymerase II transcription factor activity, transcription coactivator binding, enzyme binding, and sequence-specific DNA binding (Fig. 6A, Table S4).

The main Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways involved in SLE by Scutellaria baicalensis include the Pathways in cancer, Lipid and atherosclerosis, Kaposi sacroma, measles, and so on (Fig. 6B, Table S5).

Figures 7A and 7B showed the top 50 enriched pathways, which can be divided into four categories: cellular processes (cell growth and death, cell communication), organism systems (endocrine system, immune system, nervous system), human diseases (tumors, cardiovascular diseases, infectious diseases, and endocrine diseases), and signal transduction (HIF-1 signaling pathway, PI3K-Akt signaling pathway, TNF signaling pathway, MAPK signaling pathway, and VEGF signaling pathway). Pathways unrelated to SLE were removed from the KEGG enrichment results, and a total of 28 enrichment results were obtained (Table 2). A targeted pathway network was constructed using Cytoscape 3.7.2, and the targeted areas and color tones in the network were arranged according to the Degree value. Among them, RELA, FOS, BCL2, PRCKA, CASP3, and CASP8 target genes rank at the top (Fig. 8). Each pathway interacts with common targets, indicating that Scutellaria baicalensis can treat SLE by coordinating multiple pathways.

Molecular docking validation

Based on the above analysis, molecular docking was conducted between 5 core components (baicalin, acacetin, oroxylin A, neobaicalein, and stigmasterol) and 5 core targets (IL6, CASP3, BCL2, ESR1, and FOS). Affinity ≤ −4.25 kcal/mol indicates binding activity between the ligand and the target; affinity ≤ −5.0 kcal/mol indicates good binding activity; affinity ≤ −7.0 kcal/mol indicates strong docking activity. The best activities were selected in BCL2-stigmasterol (−8.1, Fig. 9A), IL-6-baicalein (−7.6, Fig. 9B), ESR1-baicalein (−7.6, Fig. 9C), and ESR1-acacetin (−7.5, Fig. 9D).

Molecular dynamics simulation

To explore the stability of the protein-ligand interactions in greater depth, with a binding energy of ≤−7.6 kcal/mol as the limit, we selected the first four molecular docking systems for MD simulations: BCL2-stigmastero, IL-6-baicalein, ESR1-baicalein, and ESR1-acacetin, using GROMACS 2020.6 software. MD results such as RMSD, Rg value, and SASA provide an important basis for measuring the stability of complex systems of compounds and proteins, and the stability of protein tertiary structures after combining small molecules. As shown in Fig. 10A, the RMSD of all the complexes and proteins gradually stabilizes as the simulation progresses, indicating that the complexes are gradually stabilizing. The radius of gyration (Rg) can be used to describe the changes in the overall structure and is indicative of the compactness of protein structures, a larger change in Rg suggests that the system is more expanded. As depicted in Fig. 10B, the Rg values of the complexes remained stable throughout the simulation. The buried surface area (Buried SASA) of a small molecule embedded within a protein can reflect the size of the binding interface between the small molecule and the protein. As shown in Fig. 10C, the Buried SASA remains essentially stable, indicating that the contact area between the small molecule and the protein remains constant. The number of hydrogen bonds reflects the strength of protein-ligand binding, with ESR1-baicalein showing the highest hydrogen bond density and strength (Fig. 10D).

The binding free energies (∆Gbind) of the four protein-ligand complexes were further calculated using the MM/GBSA method, three replicates per experiment, and the mean values + SD (E) are shown (Table 3). ΔEele represents the electrostatic interactions between the small molecule and the protein, ΔEvdw represents the van der Waals interactions, ΔEpol represents the polar solvation energy, which can indicate the electrostatic potential energy, and ΔEnonpol represents the nonpolar solvation energy, which can indicate the hydrophobic interactions. ΔEMMPBSA is the sum of ΔEele, ΔEvdw, ΔEpol, and ΔEnonpol, and the Gibbs Binding Energy (binding energy) ΔGbind is the sum of ΔEMMPBSA and −TΔS. As shown in Table 3, IL-6-Baicalein shows the lowest ∆Gbind value, indicating the strongest binding. Decompose the ΔEMMPBSA to obtain the contribution of each amino acid to the whole binding energy, the important residues were shown in Fig. S1A. We selected the top two residues and analysed the evolution of the distances between ligand features and the top two residues (Fig. S1B).

Analysis of potential active components

In this study, we preliminarily screened the effective components of Scutellaria baicalensis for treating SLE using network pharmacology, including baicalin, acacetin, quercetin, flavonoids of baicalin, and stigmasterol, which are also confirmed by high-performance liquid chromatography (HPLC) (Luo et al., 2024; Wang et al., 2025). Baicalin, a flavonoid compound extracted from the root of Scutellaria baicalensis, is the main active component of Scutellaria baicalensis and has been found in recent years to have good anti-tumor effects, including inhibiting cell proliferation, inducing cell apoptosis, autophagic cell death, and inhibiting cell adhesion, migration, and invasion. Animal experiments have shown that baicalin can improve renal function and pathological manifestations in lupus model mice by downregulating the activation of NLRP3 inflammasomes and the phosphorylation level of NF-κB (Li et al., 2019). Acacetin exhibits multiple effects, including antitumor, anti-inflammatory, and antiviral activities. Acacetin can inhibit the expression of Toll-like receptor 4 (TLR 4), tumor necrosis factor-related apoptosis-inducing ligand (TRAIL-R1), interleukin-6 (IL-6), tumor necrosis factor-alpha (TNF-α), and interleukin-5 (IL-5), thereby showing anti-inflammatory potential (Singh et al., 2020). In terms of anti-rheumatic disease, acacetin can inhibit the expression of matrix metalloproteinases in fibroblast-like synoviocytes through the MAPK signaling pathway, thereby playing a role in anti-rheumatoid arthritis (Chen et al., 2015).

Neobaicalein is a flavonoid compound isolated from Scutellaria baicalensis with significant anti-inflammatory effects. T helper cells 17 (Th17) and regulatory T (Treg) cells play an important role in the pathological process of autoimmune diseases. The imbalance of Th17/Treg cells is related to the occurrence and development of organ inflammation in systemic lupus erythematosus (Shan, Jin & Xu, 2020). In vitro experiments have confirmed that neobaicalin can inhibit the activation of STAT3 signaling by blocking the IL-17 signaling pathway, thereby improving the ratio of Th17 cells to Treg cells (Chen et al., 2022). It has been shown that stigmasterol can activate the NF-κB pathway through AMPK, thereby playing an anti-inflammatory role.

Targets analysis

This study indicates that the main targets regulated by Scutellaria baicalensis in SLE are IL-6, CASP3, BCL 2, ESR1, FOS, RELA, and PRCKA. IL-6 is related to T cell-dependent B cell activation and antibody production (Dienz et al., 2009), reduced natural killer (NK) cell cytotoxicity, and the amplification of pro-inflammatory cytokine cascades (Crayne et al., 2019). It has been reported that there is an increased B cell response to IL-6 in lupus patients (Kitani et al., 1992), and the inhibition of IL-6 levels in identified lupus patients helps to restore the balance of these naive B cells and T cells (Shirota et al., 2013). Anti-IL-6 treatment studies for SLE have been completed, in which the anti-IL-6 monoclonal antibody PF-04236921 was given to patients with active SLE and a placebo. The SLE response index-4 (SRI-4) and BICLA response rates were significantly different between the 10 mg treatment group and the placebo group in patients with higher disease activity (Wallace et al., 2017). Caspase-3 has a proteolytic role in programmed cell death (Nagata, 2018). A characteristic of T cells in SLE patients is an abnormal TCR-mediated signal transduction response, and it has been shown that casp3 can restore normal T cell function by limiting the key molecule CD3ξ of T cell receptor (Krishnan et al., 2005). FOS, as a universal immediate-early nuclear transcription factor, is involved in the regulation of cell proliferation, differentiation, transformation, and apoptosis. In animal experiments, FOS transgenic mice were more resistant to experimental autoimmune encephalomyelitis by inhibiting the production of inflammatory cytokines in dendritic cells (Wang et al., 2021b). The Fosl 1-JunB complex is related to the synthesis of IL-17 and the differentiation of Th17 cells. Their overexpression can lead to severe inflammatory symptoms and increased IL-17 expression. The RELA gene encodes the transcription factor p65, a subunit of NF-κB. The NF-κB signal is a key regulator of the immune system and inflammatory response. It has been shown that miR1013p can reduce the inflammatory response of peripheral blood mononuclear cells in SLE by inhibiting the expression of MAPK1 and blocking the NF-κB pathway (Zhao et al., 2021). PRCKA encodes a member of the serine/threonine-specific protein kinase family, protein kinase C (Coussens et al., 1986), which is involved in mediating cell growth and inflammatory responses (Cheng et al., 2020). PRKCA inhibits the release of pro-inflammatory cytokines through macrophages and the MAPK signaling pathway, thereby alleviating acute lung injury in lipopolysaccharide-induced mice (Wang et al., 2021c).

Biological enrichment analysis

GO analysis found that Scutellaria baicalensis plays a role in treating SLE by regulating oxidative stress, transcription-related, and apoptosis. KEGG signal pathway analysis indicates that the main signal pathways targeted by Scutellaria baicalensis intervention in SLE include the HIF-1 signaling pathway, PI3K-Akt signaling pathway, MAPK signaling pathway, IL-17 signaling pathway, which involves inhibiting inflammatory factors, affecting the expression and activation of immune cells, and regulating cell autophagy.

HIF-1 plays an important role in cellular hypoxia, inflammation, and tumorigenesis. In T cells, HIF-1 enhances the production of inflammatory cytokines and cytotoxic molecules. HIF-1α directly activates the expression of RORγt and inhibits the transcriptional activity of FoxP 3, regulating the development of Treg and Th17 cells (Dang et al., 2011). Overexpression of HIF-1α has been detected in human lupus CD4+ T cell subsets, indicating that abnormal expression of HIF-1α may be involved in the immune dysregulation of SLE (Garchow, Maque Acosta & Kiriakidou, 2021). Animal experiments have shown that increased expression of HIF-1 in CLE skin-infiltrating T cells promotes skin tissue damage, and inhibiting HIF-1 can eliminate disease progression (Little et al., 2023). The same changes have been observed in lupus nephritis mouse models (Chen et al., 2020). The PI3K-Akt signaling pathway can regulate a wide range of cellular processes, including cell proliferation, growth, and metabolism, especially in the nervous system. Excessive activation of the PI3K/AKT pathway in gliomas promotes the proliferation, invasion, and migration of glioma cells (Gao et al., 2022). The mitogen-activated protein kinase (MAPK) cascade signaling pathway is a key signaling pathway that regulates various cellular processes, including proliferation, differentiation, apoptosis, and stress responses. Animal experiments suggest that 1,25-dihydroxyvitamin D3 can significantly inhibit the NF-κB and MAPK signaling pathways to improve lupus activity, providing a potential target for the treatment of LN (Li et al., 2022).

Interleukin-17 (IL-17) is the signature cytokine of Th17 cells. The dysregulation of IL-17 expression is associated with a variety of inflammatory diseases, including inflammatory bowel disease (IBD), multiple sclerosis (MS), psoriasis, systemic lupus erythematosus (SLE), rheumatoid arthritis (RA), and asthma (Patel & Kuchroo, 2015). A study report showed that IL-17-dependent signaling in BXD 2 mouse cells leads to the production of autoantibodies and an increase in the number of germinal centers. In this study, IL-17-induced genes encode G protein signaling regulators 13 and 16 (Rgs 13 and Rgs 16), which inhibit the movement of B cells in the lymph node and promote the spontaneous formation of autoreactive germinal centers, thereby increasing the production of autoantibodies (Hsu et al., 2008). Significantly elevated levels of IL-17 can be detected in the plasma and serum of SLE patients (Crispin et al., 2008).

Molecular docking and dynamic simulation

Computational simulations were used to explore the binding affinity and mechanisms of Huangqin with key SLE-related targets. The molecular docking results show that the docking scores were all less than −5 kcal/mol. The best activities were selected in BCL2-stigmasterol, IL-6-baicalein, ESR1-baicalein, and ESR1-acacetin, primarily driven by hydrogen bonding and hydrophobic interactions. Subsequently, molecular dynamics simulations were conducted to analyze the dynamic behavior and stability of the complexes over time. The values of RMSD, Rg, and SASA and the number of hydrogen bonds indicated that the four complexes possessed reliable structural stability and compactness. Additionally, the binding free energies (∆Gbind) of the four protein–ligand complexes were analyzed by the MM/PBSA method, with the IL-6-baicalein showing the strongest binding affinities, which might inspire targeted drug delivery systems (Pingping, Nan & Yong, 2025) or derivative compounds optimized for efficacy and bioavailability (Gao et al., 2025).