Integrating Network Pharmacology with Molecular Docking to Unravel the Active Compounds and Potential Mechanism of Simiao Pill Treating Rheumatoid Arthritis

Objective To explore the main components and unravel the potential mechanism of simiao pill (SM) on rheumatoid arthritis (RA) based on network pharmacological analysis and molecular docking. Methods Related compounds were obtained from TCMSP and BATMAN-TCM database. Oral bioavailability and drug-likeness were then screened by using absorption, distribution, metabolism, and excretion (ADME) criteria. Additionally, target genes related to RA were acquired from GeneCards and OMIM database. Correlations about SM-RA, compounds-targets, and pathways-targets-compounds were visualized through Cytoscape 3.7.1. The protein-protein interaction (PPI) network was constructed by STRING. Gene Ontology (GO) analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis were performed via R packages. Molecular docking analysis was constructed by the Molecular Operating Environment (MOE). Results A total of 72 potential compounds and 77 associated targets of SM were identified. The compounds-targets network analysis indicated that the 6 compounds, including quercetin, kaempferol, baicalein, wogonin, beta-sitosterol, and eugenol, were linked to ≥10 target genes, and the 10 target genes (PTGS1, ESR1, AR, PGR, CHRM3, PPARG, CHRM2, BCL2, CASP3, and RELA) were core target genes in the network. Enrichment analysis indicated that PI3K-Akt, TNF, and IL-17 signaling pathway may be a critical signaling pathway in the network pharmacology. Molecular docking showed that quercetin, kaempferol, baicalein, and wogonin have good binding activity with IL6, VEGFA, EGFR, and NFKBIA targets. Conclusion The integrative investigation based on bioinformatics/network topology strategy may elaborate on the multicomponent synergy mechanisms of SM against RA and provide the way out to develop new combination medicines for RA.


Introduction
Rheumatoid arthritis (RA) is a chronic polyarticular symmetric disease. It is characterized by chronic inflammation of the synovial membrane, which can destroy articular cartilage and juxta-articular bone [1]. RA affects 0.3%-1% of the population worldwide [2]. If insufficiently treated, it usually leads to persistent joint inflammation, progressive joint destruction, continuing functional decline, extra-articular manifestations, disability, and increased mortality [3,4]. Although current available therapeutic approaches against RA, including nonsteroidal anti-inflammatory drugs (NSAIDs), disease-modifying antirheumatic drugs (DMARDs), and corticosteroid, allow for excellent disease control, novel therapies are needed because RA remains incurable [5]. Furthermore, the long-term use of these drugs may cause multiple side effects and lead to limited therapeutic responses. erefore, novel treatments are in urgent demand.
Traditional Chinese medicine (TCM) has been extensively applied for the treatment of RA for centuries in Asia and has been gradually accepted for worldwide clinical applications [6,7]. Numerous studies have indicated that TCM can be served as complementary and alternative RA drugs for therapeutic effects and with fewer side effects [8,9].
Nowadays, network pharmacology integrates network biology and polypharmacology based on existing databases, providing a novel approach for exploring the mechanisms and synergistic effect of TCM formulas as disease treatments [14][15][16]. Combining network science with ancient TCM formulas to investigate multiple molecular mechanisms has achieved successful attempts in the previous researches [17][18][19][20]. erefore, in this study, a network pharmacology-based study was conducted to predict bioactive compounds and elucidate the comprehensive pharmacological mechanisms about the antirheumatic effect of SM. In addition, molecular docking analysis was performed to validate in silico to predict molecular interactions between compounds and targets.

Materials and Methods
Network pharmacology-based prediction of SM treating RA was constructed by the following (Figure 1): (1) data collection and preparation, including retrieving the ingredients list of SM formula, screening for candidate compounds, identifying SM and RA targets, and intersecting the identified targets of compounds and disease; (2) topological analysis of network and protein-protein interaction (PPI) network construction; (3) enrichment analysis; and (4) molecular docking analysis. e in silico integrative ADME (absorption, distribution, metabolism, and excretion) model administrated by TCMSP is employed for pharmaceutical research. As an oral drug, two related-ADME models, oral bioavailability (OB), and drug-likeness (DL) are applied to identify the potential bioactive compounds in this study. Only the compounds with OB ≥ 30 and DL ≥ 0.18 that satisfied the criteria suggested by the TCMSP database (removed the duplicated) are retained as the candidate compounds for further study [21]. In addition, among the compounds with OB ＜30 or DL ＜0.18, which are searched with "compound (name)" and "rheumatoid arthritis" [all fields] in PubMed databases to find relevant researches, the compounds in purified form focused on anti-RA mechanisms are also considered to be bioactive compounds (removed the duplicated) and included for further study.

Predictions of Target Genes Related to the Identified
Compounds. All the potential compounds were input into TCMSP to capture the relationships between drugs and targets. Since the obtained targets include various biological species, all target names were also put into UniProt databases (http://www.uniprot.org/) to search for target gene names selected by human species.

Venn Analysis.
All target genes of identified compounds and RA are put into Bioinformatics and Evolutionary Genomics system (bioinformatics.psb.ugent.be/ webtools/Venn/), respectively, to produce a Venn diagram, which indicates the intersection of identified targets of drug and disease. e above 77 target genes acquired from the Venn diagram intersection were imported into STRING (https://string-db.org/, version 11.0) to construct a PPI network for understanding protein interaction systematically. e PPI network is constructed by setting the organism as "human sapiens", setting the minimum required interaction score to "medium confidence (0.40)", and excluding the disconnected protein nodes. In addition, statistics of protein interactions are figured out according to the PPI network, and a related bar plot diagram is constructed with R 3.6.0 subsequently.

Enrichment
Analysis. R 3.6.0 and related R packages (colorspace, stringi, DOSE, clusterProfiler, and pathview) are applied to carry out Gene Ontology (GO) enrichment analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of intersection target genes of SM-RA. P values ＜0.05 and q values ＜0.05 are considered statistically significant based on Fisher's test.

Molecular Docking Analysis.
e 3D structures of candidate targets were obtained from the PDB database (http://www.rcsb.org/) in PDB format by setting the organism to "Homo sapiens only". e 3D conformers of candidate compounds are acquired from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/) with SDF format. Subsequently, they were imported to the Molecular Operating Environment (MOE) to get the docking score. e greater the absolute value of the docking score, the better.

Identified Disease Target Genes.
e target genes related to RA were searched in GeneCards and OMIM databases, which include 3768 genes in GeneCards and 1 gene in OMIM, with no overlapping target gene.

Intersection of Identified Targets of Compounds and
Disease. In the Venn diagram intersection of identified targets about identified compounds and of RA ( Figure 2), a total of 77 target genes are acquired.

PPI Network.
e PPI network is established by setting the confidence level of more than 0.40 and hiding the independent target protein nodes. e PPI network nodes represent proteins and edges represent protein-protein interactions. e network has 75 nodes and 1604 edges ( Figure 5). In addition, we analyzed the importance prioritization (adjacent nodes count of each protein) of proteins according to the network, and the leading 30 genes with higher connection were visualized by constructing a bar plot diagram ( Figure 6), which indicates the 30 genes or proteins that may play a bridge role in connecting other nodes in the PPI network. ese 30 genes or proteins include inflammation-associated genes (IL6 [36], NFKBIA [37]), cell proliferation-, differentiation-, and transformation-related genes (FOS [38], EGFR [39], MAPK8 [40], NR3C1 [41], RHOA [42], and PARP1 [43]), cell apoptosis-related genes (CASP3, CASP8 [44], CASP9 [45], MYC [46], CYCS [47], HIF1A [48], MCL1 [49], and GSK3B [50]), cell cycle-related gene (CCND1 [51]), hormone-related genes (INS [52], ESR1 [53], AR [54], and PGR [55]), angiogenesis-related gene (VEGFA [56]), and transcription factor (RELA [57]). Figure 7, the top 20 enrichment terms are visualized by the bar plot diagram. e results demonstrated that numerous targets are involved in various BPs associated with immune response and inflammation, such as the response to a steroid hormone, response to oxidative stress, and regulation of the apoptotic signaling pathway, which confirmed strongly the correlation with the pathogenesis in RA. e CC results showed that most of the targets are localized to the cellular membrane and nuclear chromatin part. e MF results indicated that many targets are associated with nuclear receptor activity and transcription factor activity.

KEGG Enrichment Analysis.
e KEGG pathways are applied to examine the function and signaling pathways of the identified target genes, with the top 20 of the potential pathways (P < 0.05 and q < 0.05) shown by a bar plot diagram ( Figure 8) and visualized with the pathways-targetscompounds network (Figure 9). e results showed that numerous targets are associated with certain virus infections (such as Epstein-Barr virus infection) and cancer, which are associated with the onset and prognosis of RA.

Molecular Docking Analysis.
e selected targets, including IL6, VEGFA, EGFR, and NFKBIA, play a significant role in the SM-RA network. e candidate compounds, including quercetin, kaempferol, baicalein, and wogonin, are the top 4 compounds (ranking by related target genes count) in the SM-RA network. ese 4 target genes and 4 compounds are imported into MOE for molecular docking verification. e docking scores are shown in Table 3. e action mode of NFKBIA and quercetin, kaempferol, baicalein, and wogonin and the action mode of wogonin and IL6, VEGFA, EGFR, and NFKBIA are shown in Figure 10.

Discussion
In the present network pharmacological analysis, a total of 479 compounds were identified in the four herbal medicines of SM, and 72 compounds were yielded by ADME criteria screening. A total of 386 targets related to potential compounds and 3769 targets associated with RA were identified, and 77 target genes were obtained from the interaction of targets about SM identified compounds and RA. SM-RA network analysis visualized the interaction of multicomponents and multitargets about SM on RA. e compoundstargets network analysis indicated that the 6 compounds, including quercetin, kaempferol, baicalein, wogonin, betasitosterol, and eugenol, were linked to ≥10 target genes, and the 10 target genes (PTGS1, ESR1, AR, PGR, CHRM3, PPARG, CHRM2, BCL2, CASP3, and RELA) were core target genes in the network. GO enrichment analysis indicated that numerous targets are involved in response to a steroid hormone, oxidative stress, and regulation of the About 72 identified compounds, particularly the 6 compounds, including quercetin, kaempferol, baicalein, wogonin, beta-sitosterol, and eugenol, were linked to more than 10 targets, indicating that these compounds might play a vital role in the process of RA treatment. Furthermore, certain compounds have exhibited the potential antirheumatic therapeutic activities except for wogonin (Table 4). For instance, quercetin has been reported to decrease levels of tumor necrosis factor-α (TNF-α), interleukin-1β (IL-1β), interleukin-17 (IL-17), and monocyte chemotactic protein-1 (MCP-1) [58] and significantly reduced damage to interchondral joints, infiltration of inflammatory cells, and pannus formation [59]. Besides, kaempferol suppresses the proliferation and migration of RAFLS and the release of activated T-cell-mediated inflammatory cytokines and reduces osteoclast differentiation through targeting on the fibroblast growth factor receptor 3-(FGFR3-) ribosomal S6 kinase 2 (RSK2) signaling axis [60]. In addition, baicalein inhibits human rheumatoid arthritis fibroblast-like synoviocytes (RAFLS) proliferation involving suppression of nuclear factor kappa B (NF-κB) transcriptional activity and recombinant macrophage migration inhibitory factor-(MIF-) mediated signaling [61]. What's more, β-Sitosterol could modulate the functions of macrophages and attenuates rheumatoid inflammation in CIA mice [62]. For eugenol, it is reported to be effective in ameliorating oxidative stress and inflammation in arthritic rats [25,26]. Moreover, among the other 66 compounds, some articles previously reported the antirheumatic effect. For example, ferulic acid is reported to suppress osteoclast differentiation and bone erosion via the inhibition of receptor activator of nuclear factor lB ligand-(RANKL-) dependent NF-κB signaling pathway [63], and berberine could attenuate adjuvant-induced arthritic fibroblast-like synoviocytes (AA-FLS) proliferation and regulate the 17/Treg imbalance [64]. Collectively, these active components exhibit antirheumatic effects from various aspects, including anti-inflammatory, immunoregulatory, reducing bone erosion and destruction, and attenuating oxidative stress. erefore, these might indicate the collective effectiveness and diversity of constituents in SM for treating RA.
Of the leading 30 target genes with a higher connection in the PPI network, IL6, VEGFA, EGFR, and NFKBIA play a critical role in the development of RA, which has been aforementioned. Besides, in the visualized pathways-targets network, IL6, VEGFA, EGFR, and NFKBIA are involved in numerous pathways, indicating that SM may exert anti-RA effects through multipathways and multitargets combined interaction. Furthermore, the molecular docking analysis was constructed to investigate the interaction of some candidate compounds and targets. For example, the absolute value of docking scores about NFKBIA and quercetin, kaempferol, baicalein, and wogonin is the highest in each group, indicating that NFKBIA has a higher binding affinity than other target genes. For wogonin, although there have been no relevant studies about the effect in RA, the docking results indicated that wogonin performed good binding activity with IL6, VEGFA, EGFR, and NFKBIA. In brief, the high binding affinities of these active components indicated that the therapeutic effects of SM treating RA were probably through the modulation of several related targets.
As shown, the anti-RA effect of identified compounds (quercetin, kaempferol, baicalein, beta-sitosterol, and eugenol) is partially associated with the potential target genes, including NFKBIA, IL6, and MAPK, and potential signals, including PI3K-Akt, TNF, and IL-17, indicating the interaction between multicomponents, multitargets, and multisignaling of SM treating RA. Figure 9: e pathways-targets-compounds network. e green diamonds represent pathways, the blue-purple ellipses represent genes, and the pink V's represent compounds.

Limitation
is study has some limitations. It provides only a predictive overview of the pharmacological mechanisms of SM against RA based on the existing database, and further experiment verification in vivo and in vitro is necessary to ensure the reliability and reasonability of predicted results. First, posttranscriptional processing, translation regulation, and posttranslational processing and regulation play a critical role in gene expression regulation, and most of the research on mechanisms about SM treating RA is gene level in this study; therefore, an in-depth study needs to explore the related mechanism. Second, the key proteins and KEGG pathways need to be verified. ird, the anti-RA effect needs to be further verified in the animal model.
Besides, for clinical applying of SM treating RA, owing to the ethnic, genetic, and possible etiological differences, potential mechanisms about the related therapeutic module coinciding with clinical applications are worthy of further experimental investigation. In addition, importantly, dosage exploration, oral bioavailability, water-solubility, pharmacokinetics, and potential side effects of SM will also need a thorough exploration.

Conclusion
In summary, a bioinformatics/topology-based strategy, including ADEM screening, bioinformatics, network topology, enrichment analysis, and molecular analysis, was applied for identification of the molecular mechanisms of SM against RA. e integrated strategy might make the decipherment of biological mechanisms more accurate and efficient. e SM-RA network, compounds-targets network, and pathways-targets network analysis visualized the interaction of multicomponents and multitargets about SM treating RA. In particular, quercetin, kaempferol, baicalein, wogonin, beta-sitosterol, and eugenol might be the candidate therapeutic agents, and PTGS1, ESR1, AR, PGR, CHRM3, PPARG, CHRM2, BCL2, CASP3, and RELA were identified as potential drug targets. e enrichment and PPI analysis revealed the biological functions of the grouping networks related to the pathogenesis of RA. e multicomponent cosynergism of the herbal combinations about SM was elaborated. e study also revealed the multifunctional synergetic mechanisms of SM, including certain virus infection and cancer, and PI3K-Akt, TNF, and IL-17 signaling pathway.
Data Availability e figures and tables used to support the findings of this study are included within the article, and the original data are available from the first author or corresponding author upon request.

Disclosure
is research did not receive any specific funding.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.