Uncovering Mechanisms of Zanthoxylum piperitum Fruits for the Alleviation of Rheumatoid Arthritis Based on Network Pharmacology

Simple Summary The aim of the study is to investigate the bioactives of Zanthoxylum piperitum fruits on rheumatoid arthritis. The methodology to identify the relationship between signaling pathways, targets, and bioactives is based on network pharmacology. The results show that Zanthoxylum piperitum fruits might alleviate inflammatory symptoms of rheumatoid arthritis. Thus, we suggest that Zanthoxylum piperitum fruit is a promising herbal plant to reduce the level of cytokines against rheumatoid arthritis. Abstract Zanthoxylum piperitum fruits (ZPFs) have been demonstrated favorable clinical efficacy on rheumatoid arthritis (RA), but its compounds and mechanisms against RA have not been elucidated. This study was to investigate the compounds and mechanisms of ZPFs to alleviate RA via network pharmacology. The compounds from ZPFs were detected by gas chromatography–mass spectrometry (GC-MS) and screened to select drug-likeness compounds through SwissADME. Targets associated with bioactive compounds or RA were identified utilizing bioinformatics databases. The signaling pathways related to RA were constructed; interactions among targets; and signaling pathways-targets-compounds (STC) were analyzed by RPackage. Finally, a molecular docking test (MDT) was performed to validate affinity between targets and compounds on key signaling pathway(s). GC-MS detected a total of 85 compounds from ZPFs, and drug-likeness properties accepted all compounds. A total of 216 targets associated with compounds 3377 RA targets and 101 targets between them were finally identified. Then, a bubble chart exhibited that inactivation of MAPK (mitogen-activated protein kinase) and activation of PPAR (peroxisome proliferator-activated receptor) signaling pathway might be key pathways against RA. Overall, this work suggests that seven compounds from ZPFs and eight targets might be multiple targets on RA and provide integrated pharmacological evidence to support the clinical efficacy of ZPFs on RA.


Introduction
Rheumatoid arthritis (RA) is a long-term systemic autoimmune disorder that deteriorates the synovial joints and is associated with gradual disability [1]. RA is a progressive inflammation caused by joint damage and its functional loss around the articular [2,3]. RA can present irrespective of age, diagnosed in around 1% of the population, brings huge social-economic burden [4]. The main factor causing RA is uncontrollable cytokine secretion due to bone damage; however, the etiology of RA is unknown [5,6]. Commonly, the anti-RA drugs administered are disease-modifying arthritis drugs (DMARDs) and nonsteroidal anti-inflammatory drugs (NSAIDs) in most countries [7]. At present, prolonged administration of these drugs is involved in severe side effects such as upset stomach, nephrotoxicity, and electrolyte imbalance [8,9]. In contrast, an animal test demonstrated

Plant Material Collection and Classification
The Zanthoxylum piperitum fruits (ZPFs) were collected from (latitude: 37.628975, longitude: 126.742978), Gyeonggi-do, Korea, in August 2020, and the plant was identified by Dr. Dong Ha Cho, Plant biologist and Professor, Department of Bio-Health Convergence, College of Biomedical Science, Kangwon National University. A voucher number (UUC 270) has been stored at Kenaf Corporation in the Department of Bio-Health Convergence, and the material can be used only for research purposes.

Plant Preparation, Extraction
The ZPFs were dried in a shady area at room temperature (20-22 • C) for 21 days, and dried ZPFs made powder using an electric blender. Approximately 50 g of ZPFs powder was soaked in 800 mL of 100% methanol (Daejung, Siheung city, Gyeonggi-do, Korea) for 10 days and repeated three times to achieve a high yield rate. The solvent extract was collected, filtered, and evaporated using a vacuum evaporator (IKA-RV8, Staufen city, Germany). The evaporated sample was dried under a boiling water bath (IKA-HB10, Staufen city, Germany) at 40 • C for around 8 h to obtain solid extraction.

GC-MS Analysis Condition
Agilent 7890A was used to carry out GC-MS analysis. GC was equipped with a DB-5 (30 m × 0.25 mm × 0.25 µm) capillary column. Initially, the instrument was maintained at a temperature of 100 • C for 2.1 min. The temperature was rose to 300 • C at the rate of 25 • C/min and maintained for 20 min. Injection port temperature and helium flow rate were ensured as 250 • C and 1.5 mL/min, respectively. The ionization voltage was 70 eV. The samples injected in split mode at 10:1. The MS scan range was set at 35-900 (m/z). The fragmentation patterns of mass spectra were compared with those stored in the using W8N05ST Library MS database (analyzed 28 April 2021). The percentage of each compound was calculated from the relative peak area of each compound in the chromatogram.
The concept of integration was used in the ChemStation integrator algorithms (analyzed 28 April 2021) [20].

PPI Networks and Bubble Chart
STRING (https://string-db.org/) (accessed on 19 May 2021) [23] was utilized to analyze the PPI network with final overlapping targets. The RPackage was utilized to identify the degree of value. Then, signaling pathways associated with the occurrence and development of RA were visualized on a bubble chart by RPackage; two key signaling pathways with the highest and the lowest rich factor were selected to analyze the relationships against RA.

Construction of STC Networks
The STC networks were utilized to construct a size plot based on the degree of values. In this size map, green rectangles (nodes) represented signaling pathways; gold triangles (nodes) stood for target proteins, and red circles (nodes) stood for compounds; its size represented degree value. The size of gold triangles represented the amount of connectivity with signaling pathways; the size of red circles represented the amount of connectivity with target proteins. The combined networks were constructed by using RPackage (analyzed 20 May 2021).

Preparation of Compounds from ZPFs for MDT
The ligand molecules were converted .sdf from PubChem into .pdb format using Pymol (accessed on 21 May 2021), and the ligand molecules were converted into .pdbqt format through Autodock (accessed on 21 May 2021).

Ligand-Protein Docking
The ligand molecules were docked with target proteins utilizing autodock4 by setting up four energy ranges and eight exhaustiveness as default to obtain 10 different poses of ligand molecules [24]. The 2D binding interactions were used with LigPlot+ v.2.2 (https://www.ebi.ac.uk/thornton-srv/software/LigPlus/) (accessed on 21 May 2021). After docking, ligands of the lowest binding energy (highest affinity) were selected to visualize the ligand-protein interaction in Pymol (Schrödinger, New York, NY, USA).

Toxicological Properties Prediction by AdmetSAR
Toxicological properties of the key bioactive were established using the admetSAR web-service tool (http://lmmd.ecust.edu.cn/admetsar1/predict/) (accessed on 23 May 2021) because toxicity is an essential factor to develop new drugs. Hence, Ames toxicity, carcinogenic properties, acute oral toxicity, and rat acute toxicity were predicted by admetSAR (East China University of Science and Technology, Shanghai, China).

Chemical Compounds from ZPFs
A total of 85 chemical compounds and its seven key chemical compounds in ZPFs were detected by the GC-MS analysis (Figure 2), and the name of compounds, PubChem ID, retention time (mins), and peak area (%) were enlisted in Table 1. Lipinski's rules accepted all 85 compounds (molecular weight ≤500g/mol; Moriguchi octanol-water partition coefficient ≤4.15; number of nitrogen or oxygen ≤10; number of NH or OH ≤5), and all chemical compounds satisfied with the criteria of "Abbott Bioavailability Score (>0.1)" through SwissADME. The TPSA (topological polar surface area) value of chemical compounds was also accepted ( Table 2).

Overlapping Targets between SEA and STP Related to Chemical Compounds
A total of 470 targets from Similarity Ensemble Approach (SEA) and 629 targets from SwissTargetPrediction (STP) were associated with 85 chemical compounds (Supplementary Table S1). The Venn diagram showed that 215 targets overlapped between the two public databases (Supplementary Table S1) ( Figure 3A).

Overlapping Targets between RA-Related Genes and the Final 101 Overlapping Targets
A total of 3377 targets associated with RA were selected by retrieving from DisGeNET, Online Mendelian Inheritance in Man (OMIM) databases and literature (Supplementary  Table S2). Venn diagram's result displayed 101 overlapping targets that were selected between 3377 targets related to RA and the 215 overlapping targets ( Figure 3B) (Supplementary Table S3).

Acquisition of a Key Target from PPI Networks
From STRING analysis, 99 out of 101 overlapping targets were directly associated with RA occurrence and development, indicating 99 nodes and 469 edges ( Figure 4). The two removed targets (CA1 and CA3) had no connectivity to the overlapping 101 targets. In PPI networks, Vascular Endothelial Growth Factor A (VEGFA) was the highest degree (42) and is considered a key target (Table 3).

Identification of Two Key Signaling Pathways from a Bubble Chart
The results of Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis unveiled that 19 signaling pathways were connected with 40 out of 101 targets (false discovery rate <0.05). The 19 signaling pathways were directly related to RA, indicating that these 19 signaling pathways might be significant pathways of ZPFs against RA. The description of 19 signaling pathways were shown in Table 4. Additionally, a bubble chart indicated that inactivation of MAPK signaling pathway and activation of PPAR signaling pathway might be key signaling pathways of ZPFs to alleviate RA ( Figure 5). Specifically, MAPK signaling pathway is associated with VEGFA (a hub target) in holistic PPI networks while PPAR signaling pathway is not related to VEGFA.

Construction of Signaling Pathway-Target-Compound Network
The signaling pathway-target-compound (STC) network is displayed in Figure 6. There were 19 pathways, 40 targets, and 63 compounds (122 nodes, 488 edges). The nodes stood for the total number of relationships between signaling pathways, targets, and compounds. The edges represented the relationship of three elements (19 pathways, 40 targets, and 63 compounds). The STC network suggested that the uppermost target is Protein Kinase C Alpha (PRKCA) with 14 degrees among 19 signaling pathways (Table 5).

KEGG Pathway Enrichment Analysis
KEGG pathway enrichment analysis [25] reveals that 19 signaling pathways are associated with the occurrence and development of RA. Out of 19 signaling pathways, MAPK signaling pathway and PPAR signaling pathway were most significantly related to RA (Figure 7). The red rectangles indicated core targets on the two signaling pathways. One of the MAPK signaling pathways, PRKCA (a species of PKC), is an excellent immunomodulatory target, which is linked to the onset of inflammatory arthritis, including RA [26]. It was subsequently shown that inactivation of PRKCA results in the inhibition of c-fos; eventually, the cascade leads to the amelioration of RA in an animal model [27]. One of the PPAR signaling pathways, FABP3, is located in the upstream region to regulate lipid metabolism. It was reported that defective lipid metabolism was observed in RA patients who are persistent to proinflammatory responses [28].

MDT of 6 Targets and 23 Chemical Compounds Associated with MAPK Signaling Pathway
From both SEA and STP databases, it was uncovered that FGF1 is associated with two chemical compounds, FGF2 with four chemical compounds, VEGFA with five chemical compounds, TNFRSF1A with two chemical compounds, PRKCA with 16 chemical compounds, and PLA2G4A with four chemical compounds. Out of 33 chemical compounds, 10 chemical compounds overlapped, and finally, 23 chemical compounds were identified on the MAPK signaling pathway.

Discussion
PPI network showed that the therapeutic efficacy of ZPFs on RA was associated with 99 targets. The KEGG pathway analysis enrichment of 99 targets revealed that 19 signaling pathways were related closely to the occurrence and progression of RA, suggesting that these signaling pathways might be remedial mechanisms of ZPFs to alleviate RA. The relationships of the 19 signaling pathways with RA are briefly discussed as follows. PPAR signaling pathway: the activation of PPAR is a good strategy to alleviate RA, which can suppress the inflammatory activity of NF-κB in fibroblast-like synoviocytes (FLSs) [30]. Relaxin signaling pathway: the combined relaxin with estrogen exerts anti-inflammatory effects by dampening neutrophil function [31]. Vascular Endothelial Growth Factor (VEGF) signaling pathway: VEGF expression level in RA patients increased significantly, compared with healthy groups; moreover, patients under RA for an extended period exerted higher VEGF expression level in serum [32]. Estrogen signaling pathway: the estrogen treatment might have an inhibitory effect on RA symptoms or delay in the onset of disease, and estrogen has anti-inflammatory activity in an animal test of RA [33]. Fc epsilon RI signaling pathway: Fc epsilon RI signaling pathway is associated with inflammatory etiology, antigen-induced autoimmune reaction, and control of lipid metabolic pathway, and damage in human joint [34]. Prolactin signaling pathway: Prolactin collaborates with other proinflammatory factors to stimulate macrophages through prolactin receptors which might be a potential therapeutic target in RA [35]. Sphingolipid signaling pathway: Sphingolipid expression level in RA elevated in the serum sample, consistent with known roles of sphingolipids related to inflammation [36]. Cyclic AMP (cAMP) signaling pathway: stimuli-induced cAMP production can exert proinflammatory effects in RA [37]. It implies that the activation of the cAMP level might drive inflammatory responses. AGE-RAGE (the receptor for advanced glycation end products) signaling pathway in diabetic complications: the expression of RAGE increased in RA patients, and IL-17 and IL-1β triggered RAGE production [37]. It suggests that the inhibition of RAGE might be a good target for RA treatment. Hypoxia-Inducible Factor-1(HIF-1) signaling pathway: synovial hypoxia is characterized by RA patients, leading to inflammation, cartilage destruction, and oxidative impairment [38]. Ras-associated protein-1 (RAP1) signaling pathway: The downregulation of RAP1 in RA induces ongoing inflammation due to the overproduction of free radicals [39]. Thyroid hormone signaling pathway: joint damages in thyroid gland abnormalities are due to hypothyroidism, and thyroid hormone plays a vital role in antioxidant modulation, as indicated by in vivo and in vitro tests [40,41]. Phospholipase D signaling pathway: phospholipase D1 (PLD1) is a mediator to induce proinflammatory cytokines with a reduction of the regulatory T (Treg) cell and recruitment of Th17 cell in collagen-induced arthritis mice [42]. Gonadotropin-releasing hormone (GnRH) signaling pathway: GnRH is a potent substance to alleviate inflammation in RA patients with a high level of GnRH [43]. Ras signaling pathway: T cells in RA patients show the Ras signalling pathway's overactivation-related deeply to inflammation [44]. T cell receptor signaling pathway: T cell recruitment to the inflammatory sites can induce chronic inflammation and autoimmunity [45]. Phosphoinositide 3-kinase-Akt (PI3K-Akt) signaling pathway: a report demonstrated that PI3K-Akt signaling pathway was excessively activated, aggravating in overexpression Bcl-2, Mcl-1, and FLIP to result in unbalanced apoptosis of synovial cells, which is associated with occurrence and progression of RA [46,47]. Calcium signaling pathway: it was reported that the cytoplasmic calcium concentration in RA naïve CD4 + T cells is increased significantly [48]. It implies that the calcium signaling pathway is involved intensely in inflammatory responses.
MAPK signaling pathway: a study demonstrated that andrographolide with antiinflammatory activities has potent anti-RA efficacy, inhibiting MAPK pathway [49]. This report is consistent with our result via network pharmacology analysis.
To sum things up, the relationship between compounds and targets through the network pharmacology concept was clarified, and ligands with the most excellent affinity out of selective ligands in ZPFs were considered as potential therapeutic candidates against RA, compared with existing ligands via MDT. The PPI network suggested that VEGFA is the highest degree of value (42) on the MAPK signaling pathway. However, the affinity of ZPFs compounds on VEGFA was invalid (> −6.0kcal/mol) [29]; furthermore, existing ligand (BAW2881: −7.6 kcal/mol) had better affinity than a compound (3,4-O-Isopropylidene-dgalactose: −5.3 kcal/mol) with the greatest affinity from ZPFs. The STC network indicated that PRKCA had the greatest degree of value (14) on the MAPK signaling pathway, and monoolein (PubChem ID: 5283468) had the greatest affinity (−6.7 kcal/mol), which was superior to an existing ligand (sphingosine: −5.5 kcal/mol). On PPAR signaling pathway, β-Caryophyllene (−8.6 kcal/mol) on PPARA, Stigmasta-5,22-dien-3-ol (−8.6 kcal/mol) on PPARD, and NSC402953 (−8.2 kcal/mol) on PPARG had better affinity than existing ligands on each target. Specifically, monoolein on FABP3 had the greatest affinity (−8.9 kcal/mol). Moreover, agonists of FABP3 were not reported yet. Furthermore, there were no valid ligands on MMP1 on the PPAR signaling pathway. Taken together, the most noticeable ligand of ZPFs against RA is monoolein with dual effects on RA. Monoolein might be an antagonist on the MAPK signaling pathway or an agonist on the PPAR signaling pathway ( Figure 10). A study demonstrated that monoolein treatment in LPS-stimulated primary murine bone marrow-derived dendritic cells (BMDCs) inhibited the activation of MAPK and NF-κB; moreover, it suppressed the production of NO and iNOS in RAW264.7 cells [50]. This suggests that the inactivation of the MAPK signaling pathway is a good strategy to ameliorate inflammation associated with bone damage. Another report indicated that the expression level of PPARG in FLSs of RA patients was significantly reduced compared with healthy FLSs [30]; this implies that the activation of the PPAR signaling pathway might be a therapeutic mechanism against RA. These reports are in line with our results. Additionally, we suggested that monoolein might be an antagonist in activating PRKCA on MAPK signaling pathway and might also be an agonist to activate FABP3 on the PPAR signaling pathway.

Conclusions
This study provides eight targets, seven compounds, and two key signaling pathways of ZPFs which show promise against RA, which might be useful for multiple target combination therapies on RA. Out of seven compounds from ZPFs, monoolein shows dual effects: inactivation of MAPK signaling pathway and activation of the PPAR signaling pathway. The key pharmacological mechanisms of ZPFs against RA might be to inhibit cytokine production in synovial cells by binding on PRKCA or FABP3.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/biology10080703/s1, Table S1: The list of 470, 629, and 215 targets from SEA, STP, and overlapping target proteins between SEA and STP, respectively, Table S2: The list of 3377 target proteins associated with RA, Table S3 Funding: This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.