Integrated network pharmacology analysis and in vitro cell experiments to reveal the mechanisms of Ligusticum chuanxiong Hort. in the treatment of non-alcoholic fatty liver disease

Background: Non-alcoholic fatty liver disease (NAFLD) is a liver disease of unknown etiology. A traditional Chinese medicine Ligusticum chuanxiong Hort. (CX), it has been used about 2,000 years. Until now, the mechanism of action of CX on NAFLD remains unclear. Method: We first tested the toxicity of CX to AML12 cells with CCK-8. In vitro cell models of NAFLD were made using free fatty acid, and used Oil Red O staining tested lipid droplets. Then the active compounds of CX were collected from TCMSP and literatures, and SwissTargetPrediction, Search Tool for Interacting Chemicals, Encyclopedia of Traditional Chinese Medicin, Bioinformatics Analysis Tool for Molecular mechANism of Traditional Chinese Medicine database were used to predict the targets of the compounds. DRUGBANK, Online Mendelian Inheritance in Man, Therapeutic Target Database, DisGeNET and GeneCards database were used to predict the targets of NAFLD. Use Venn diagram to obtained the intersection targets by, and analyzed the protein-protein intersection network. Use Kyoto Encyclopedia of Genes and Genomes and Gene Ontology to forecast the function of intersection genes. Molecular docking was used to evaluate the interaction between hub gene and active ingredients. Finally, use western blotting to determine the effects of CX on PPARA, PPARG, IL1B and TNF proteins. Result: CX can reduce the production of AML12 cell lipid droplets. A total of 15 chemical components were identified from CX. Folic acid, chrysophanol and sitosterol were the main components of CX against NAFLD. ALB, TNF, PPARG and PPARA proteins were the main targets of CX in the treatment of NAFLD. PPAR signaling pathway and fatty acid degradation were closely related to anti-NAFLD. Molecular docking results shows that folic acid was the main active ingredient of CX for NAFLD treatment, and TNF is the main potential target. The cellular NAFLD model showed that CX up-regulated the expression of PPARA and PPARG protein and down-regulated inflammatory factor IL-1B and TNF expression. Conclusion: Our study suggests that CX has a therapeutic effect on NAFLDA, which may be related to the PPAR pathway and the reduction of inflammatory cytokines.


Introduction
Fatty liver (FAL) refers to abnormal lipid accumulation in parenchymal cells of the liver. Liver is the representative organ of ectopic fat. Therefore, obesity is the main cause of FAL. Normally, 5 percent of cells contain triglycerides, but in a FAL, an imbalance in fat accumulation and consumption leads to excess production of triglycerides in liver cells. Cells with less than 50% triglycerides are called fatty infiltrates, and cells with more than 50% triglycerides are called FAL [1,2]. One of the types of FAL is non-alcoholic fatty liver disease (NAFLD) [3][4][5]. NAFLD has become a global public health problem and is considered the most common cause of chronic liver disease worldwide [6,7]. NAFLD is a common disease associated with metabolic syndrome, there is hepatic steatosis not associated with high alcohol intake, it is the common cause of asymptomatic liver enzyme elevation [8][9][10]. NAFLD covers a range of diseases, from simple fat deposits in the liver to the more severe form of non-alcoholic steatohepatitis (NASH), and finally to fibrosis, cirrhosis and hepatocellular carcinoma [11][12][13]. NAFLD is also associated with extrahepatic manifestations such as chronic kidney disease, cardiovascular disease, and sleep apnea [14]. The exact cause of NAFLD is unknown, but it may be caused by insulin resistance combined with increased oxidative stress in liver cells. There is no specific treatment for NAFLD, the main treatment is weight loss [15]. Therefore, there is an urgent need to find a drug that can effectively prevent or treat NAFLD.
Chuanxiong Rhizoma, the dried rhizome of Ligusticum chuanxiong Hort. (CX), is derived from umbellifer plants [16]. CX is a traditional Chinese medicine, it contains a large of active substances, including polysaccharide, alkaloids, phenolic acids and phthalide lactones [17]. The compositions of CX are similar but have different contents in rhizome, fibrous roots and leaves [18]. CX can regulate NRf2, NF-ĸB, caspase-3, and protein kinase pathways to exert anti-apoptotic, antioxidant and anti-inflammatory effects [19]. CX polysaccharide has antioxidant, immune-promoting, antitumor and antibacterial activities [20]. In addition, CX can also be combined with other traditional Chinese medicine to form a formula and exert multiple therapeutic effects. For example, CX can combined with Hedysarum multijugum Maxim [21] and Borneolum [22] to combat stroke.
Based on systems biology theory, a more comprehensive approach to studying molecules and cells in biology [23], network pharmacology is a subject and a branch of systems and network medicine and its treatment. It has changed the way we define, diagnose and treat disease. It analyzes of biological systems networks and design multi-target drug molecules by selects specific signal nodes. Network pharmacology is characterized by predictability and emphasis on the regulation of multiple signaling pathways [24,25]. CX can treat osteoarthritis through the core targets of ALB, AKT3, IGF1, caspase-5, MAPK14, ANXA3 and MAPK1 based on the prediction of network pharmacology [26]. And CX may treat diabetic nephropathy through steroid hormone, estrogen, and thyroid hormone pathways [27]. However, it is still unclear whether CX has therapeutic effect on NAFLD.
In this study, the method of network pharmacology was used to predict the active ingredients and targets of CX in the treatment of NAFLD, and related target enrichment analysis was conducted. Then the experiment was used to verify, providing scientific basis for studying the mechanism of action of CX treated NAFLD.

Preparation of CX extract
Weigh 60 g of CX decoction piece coarse powder, add 600 mL 80% ethanol for reflux extraction for 1.5 h, and collect filtrate. 400 mL 80% ethanol was added to the CX residue for reflux extraction for 1.5 h, and filtrate was collected. Combine filtrate twice and drain. Concentrate the filtrate under vacuum pressure and remove the ethanol. The concentrate was then transferred to a 200 mL volumetric bottle and the purified water was filled to 200 mL. Finally, the extract of CX was 0.3 g/mL.

Cell culture
The study used the AML12 cells thoroughly. Cultured AML12 cells in DMEM/F-12 medium, which contains 10% FBS, 1% penicillin and streptomycin. The cell will incubate in a humidified atmosphere of 5% CO 2 at 37°C. The medium was changed every other day, and the cells were passaged when they reached the exponential growth phase.

CCK-8 determination of CX on cell viability
Cell viability treated with different CX concentrations was measured using the CCK-8. In 96-well plates, 1×10 4 cells/well cells were be cultured. When the cells reached 70-80% confluence, treated cells with different concentrations (5,10,20,40,80, and 160 μg/mL) of CX. After incubation of 24 h, the media was removed, and 100 μL CCK-8 buffer was added per well and incubated for an additional 30 minutes. Measured the absorbance of each well at 450 nm.

Construction of cell models and Oil Red O staining
When AML12 cells grew to 80%, they were given 1 mmol/L free fatty acid (FFA) mixture (oleic acid: palmitic acid = 2:1) and cultured with Dulbecco's modified eagle's medium medium containing 10% FBS for 24 h to construct the cell model of NAFLD degeneration [28]. Divided the cells to normal group, model group and CX group. After model establishment, the CX group was treated with 80 μg/mL CX extract for 24 hours. Oil Red O Stain: removed cell culture medium, washed twice with PBS, and fixed with Oil Red O solution for 30 min. Soak in 60% isopropyl alcohol for 30 seconds. Discard isopropyl alcohol, add the newly prepared Oil Red O dyeing solution, and soak for 20 min. Discard the dye and rinse with 60% isopropyl alcohol for 30 seconds. Wash away excess dye with distilled water. Redyeing the core for 2 min by Mayer hematoxylin staining solution. Add oil red O buffer for 1 min, discard, add distilled water to cover cells and observe under microscope.

The active compounds of CX
The active compounds of CX were screened from Traditional Chinese Medicines Systems Pharmacology Database and Analysis Platform (TCMSP, https://old.tcmsp-e.com/tcmsp.php) based on oral bioavailability ≥ 30% and drug-likeness ≥ 0.18. In addition, we searched the relevant literature of CX to supplement the active compounds [29][30][31][32][33][34][35]. Get the PubChem ID, formula, molecular weight (g/mol) and structure of the compound by PubChem database (https://pubchem.ncbi.nlm.nih.gov/). Targets prediction of CX Submit a manuscript: https://www.tmrjournals.com/im The active compounds obtained from CX were searched in TCMSP to obtain related targets gene name, and the gene name was converted to Gene Symbol by UniProt database (https://www.uniprot.org/). The SDF format of these compounds was input into Swiss Target Prediction database (http://www.swisstargetprediction.ch/) for target prediction, and the species was selected as "Homo sapiens", with the probability ≥ 0.6 strictly required. Enter the compounds names into Search Tool for Interacting Chemicals database (http://stitch.embl.de/) for target prediction, select "Homo sapiens" and confidence score ≥ 0.6. Search the names of compounds in Encyclopedia of Traditional Chinese Medicine (ETCM, http://www.tcmip.cn/ETCM/), and select the confidence score > 0.8, collect related targets of compounds. PubChem ID of the compounds were input into Bioinformatics Analysis Tool for Molecular mechANism of Traditional Chinese Medicine (http://bionet.ncpsb.org.cn/batman-tcm/), selected score cutoff ≥ 20 and adjusted P-value ≥ 0.05 to forecast the target of these compounds.

The intersection targets of CX and NAFLD Intersection targets of CX and NAFLD were obtained using Bioinformatics & Systems
Biology online mapping site (http://bioinformatics.psb.ugent.be/). Meanwhile, Venn diagram was use to demonstrate the relationship between CX and NAFLD targets.

PPI network construction and screening hub genes
The intersection targets of CX and NAFLD was analyzed in the Search Tool for the Retrieval of Interacting Genes/Proteins database (https://cn.string-db.org/) [36]. Input intersection targets in a multi-protein manner, set the organisms as "Homo sapiens" and selected a medium confidence (minimum required interaction score of 0.4) to construct the protein-protein interaction (PPI) network diagram [37]. Then, the PPI results was visualized by the Cytoscape software (version 3.7.2). The hub genes of CX anti-NAFLD were identified using Cytohubba plug-ins. The top 10 genes generated by maximum neighbor component method were considered as the hub genes and construct the core targets PPI by Cytoscape software [26].

Drug-compounds-target-disease (DCTD) network construction
Enter the "NETWORK" and "TYPE" files of CX-compounds, compounds-intersection targets, and intersection targets-NAFLD to Cytoscape software to construct the DCTD network diagram. Calculating the node degree, and evaluated the importance of each node in the network.

Enrichment analysis of gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG)
Through enrichment analysis of GO and KEGG to explore the function and pathway of CX anti-NAFLD in Metascape database (https://metascape.org/). We input targets as "H.sapiens" and setting the P-value cutoff as 0.01, min enrichment as 1.5, min network size and min overlap as 3 for enrichment analysis, including cell composition (CC), biological process (BP), molecular function (MF) and KEGG pathways. The results of enrichment, gene count and P-value were input into bioinformatics software (http://www.bioinformatics.com.cn/). In this study, we presented the top 10 items of CC, BP, MF and top 20 items of KEGG. Cytoscape software was used to visualize the relationship between KEGG pathway and targets.

Molecular docking
Download the 3D structure of compounds at PubChem database. Chem3D Pro 14.0 software was used to adjust the structure of these compounds to achieve the minimum binding energy. The Entry ID of hub genes is searched in UniProt database. In the Research Collaboratory for Structural Bioinformatics Protein Data Bank database (https://www.rcsb.org/) download the 3D structure Protein Data Bank format of target protein. PyMOLWin software was used to remove water and ligands from proteins. Use Discovery Studio CDOCKER software to hydrogenate the proteins and evaluate the interactions between the proteins and compounds.
Western blotting verification 5×10 6 AML12 cells were cultured in culture dishes. After modeling and CX processing, the upper medium was removed and washed 3 times with pre-cooled PBS. Total protein was obtained by adding RIPA lysis buffer in culture dish and cracking at 4°C for 30 min. After detecting the protein concentration with the BCA protein detection kit, each protein was diluted to the same concentration. Then, use SDS-PAGE protein loading buffer to denature the protein. Separated the prepared proteins by SDS-PAGE adhesive. Then transferred the protein to activated PVDF membranes. Using 5% fat-free milk incubate PVDF membranes at room temperature for 1 h. Then, using corresponding primary antibody PPARA, PPARG, IL1B and TNF (1:1000) incubate PVDF membranes at 4°C overnight. The PVDF membrane was washed 3 times with tris buffered saline with tween-20 for 10 minutes each time after incubation. Finally, use horseradish peroxidase-conjugated antibody (1:2000) incubate PVDF membrane for 1h at room temperature. Using enhanced chemiluminescence kit visualize protein bands. Set GAPDH (1:50000) as internal reference protein. Use Image J software to analysis the gray value of each blot.

Statistical analysis
Presented data as mean ± standard deviations. Student's t-test or one-way analysis of variance was used for statistical comparisons in GraphPad Prism 8 software (GraphPad Software Inc., La Jolla, CA). Highly significant was presented as ** P < 0.01 and ## P < 0.01. Statistically significant was presented as * P < 0.05.

Effects of CX on cell viability and Oil Red O staining
Use CCK-8 to detected the influence of CX on cell viability. As Figure  1A shows, it was found that CX extract had no effect on cell viability at 5, 10, 20, 40 and 80 μg/mL (P > 0.05), but CX extract at 160 μg/mL significantly reduced cell viability (P < 0.01). Therefore, CX 80 μg/mL was selected as the dosage for the follow-up experiment.
In vitro cell models of NAFLD were constructed using FFA. Detected the inhibitory effect of CX on FFA modeling by Oil Red O staining. As shown in Figure 1B, compared with the normal group, lipid droplets in the cytoplasm of the model group increased in volume, color and quantity. However, these effects can be inhibited by CX.
Next, 134 compound targets were obtained from TCMSP database, 13 from Swiss Target Prediction database, 30 from Search Tool for Interacting Chemicals database, 123 from ETCM database, and 229 from Bioinformatics Analysis Tool for Molecular mechANism of Traditional Chinese Medicine. A total of 529 targets corresponding to 15 compounds were obtained. After the deletion of duplicate targets, 370 targets remained.

Disease targets of NAFLD
A total of 208 disease targets related to NAFLD were obtained from the Online Mendelian Inheritance in Man database, 5 from DRUGBANK and 16 from TTD. A total of 903 targets with relevance score ≥ 10 were obtained from GeneCards and 157 targets with score ≥ 0.1 were obtained from DisGeNET. A total of 1289 disease targets related to NAFLD were obtained, and 1107 were left after the deletion of 182 duplicate targets.

PPI network analysis and screen hub gene
The cross targets of CX and NAFLD were represented by Venn diagram ( Figure 3A). PPI networks were used to demonstrate protein-to-protein interactions. Using the degree to demonstrate the number of connections a protein has with other proteins [38]. In the PPI network ( Figure 3B), each node represents a protein. The size and color of the node are related to the degree. The node changes from small to large, and the color changes from blue to red, indicating that the degree of the node increases. Interactions between proteins are represented by edges. The edges thickness and color represent the combined score between nodes. The edges go from thin to thick and from blue to red, representing the combined score between nodes from low to high. Ten hub genes were screened using Cytohubba plug-ins, including ALB, TNF, PPARG, PPARA, IL1B, VEGFA, ESR1, CYP3A4, and CYP2E1 (Table 2). Protein interaction analysis of 10 hub genes showed that ALB was the most core target ( Figure 3C).

DCTD network analysis
Visualization analysis of the relationships among CX, NAFLD, intersection targets, and compounds by Cytoscape software (Figure 4). The green triangle represents NAFLD, the blue hexagon represents common targets, the red oval represents active compounds, and the yellow diamond represents CX. The degree value of the node is related to the color depth and node size in hexagon and oval. The larger the node and the darker the color, the higher the degree value. Among them, the degree values of folic acid, chrysophanol and sitosterol were the highest, which may be the main components of the anti-NAFLD effect of CX. PTGS2, ESR1, RXRA and AR may be the core targets of CX compounds.

Enrichment analyses of Go and KEGG
GO and KEGG analyses were performed on the intersection targets, and the results are shown in Figures 5 and Figure 6. In the bar chart, the GO and KEGG results were represented by the Y axis and the gene counts were represented by the X axis. The color changes from purple to red, indicating that P-value changes from high to low. In bubble plot, the bubble size represents gene count, the X-axis shows the enrichment fraction GO and KEGG results. And the Y axis represents the enriched GO and KEGG results. The change of bubble color from green to red, demonstrating that P-value changes from high to low. In GO analysis, the first 10 items of molecular functions, cellular components, and biological processes were ranked according to the gene counts. In Figure 5, the results were presented. Oxidoreductase activity, nuclear receptor activity, carboxylic acid binding, tetrapyrrole binding and steroid binding were mainly response to molecular functions. The main cellular components are membrane raft, ficolin-1-rich granule lumen, secretory granule lumen, receptor complex and endoplasmic reticulum lumen. The biological processes were related to monocarboxylic acid metabolic process, regulation of hormone levels, steroid metabolic process, gland development and response to lipopolysaccharide.
In addition, multiple signaling pathways were enriched from 96 cross genes according to KEGG analysis ( Figure 6). In this study, the top 20 enrichment pathways were selected. The enrichment pathways included PPAR signaling pathway, fatty acid degradation and Non-alcoholic fatty liver disease.    In the KEGG pathway-targets interaction diagram (Figure 7), protein targets represented as purple nodes, and KEGG pathways represented as red nodes. The degree values of NFKB1, PPARG and MAPK1 were the highest, indicating that these three targets were involved in multiple KEGG pathways.

Molecular docking of active compounds with hub gene
The 15 active compounds were docked with 10 hub genes. The LibDockScore is used to evaluate the affinity of a compound with a target. With a higher LibDockScore, the compound binds to the protein were more stable and the greater the likelihood of interaction. Table 3 shows the LibDockScore of compound and protein. It was found that 5 compounds had the highest LibDockScore with TNF protein. From the perspective of compounds, except IL1B protein, folic acid had the highest LibDockScore with other 9 proteins among all compounds. We show the complexes with the highest docking fraction with 15 compounds (Figure 8-11). It indicates these compounds and proteins can be identified through carbon hydrogen bond, conventional hydrogen bond, pi-alkyl, unfavorable donor-donor, amide-pi stacked and unfavorable acceptor-acceptor .

Effect of CX on protein expression
According to KEGG results, the treatment of NAFLD by CX may be related to the PPAR signaling pathway, so we used western bolt to detect the expressions of PPARA and PPARG. The inflammation-related factors protein expressions of IL1B and TNF were detected. The results showed that the expression levels of PPARA and PPARG proteins were lower in FFA treated cells, but the expression levels of TNF and IL1B proteins were increased. However, these changes could be reversed after CX treatment ( Figure 12).

Discussion
In the world, about a quarter of people suffer from NAFLD. With the rising levels of type 2 diabetes, obesity, and metabolic syndrome, the incidence of NAFLD is increasing [39]. Hepatic steatosis, balloon-like and lobular inflammation with or without fibrosis can be used to diagnose NAFLD. The gold standard for diagnosis of NAFLD is histological estimation of liver biopsy [40]. At present, the prevalence rate of NAFLD is high, the diagnosis is difficult, the pathogenesis is complex, and the standard treatment methods are lacking [41]. A traditional Chinese medicine, CX, has been proven to prevent liver and blood vessel damage in rats following ovariectomies through antioxidant and vasorelaxant compounds [42]. Moreover, Supercritical fluid CO 2 extract of CX inhibited d-galactose induced nuclear factor protein expression up-regulation, and alleviated liver and kidney injury in mice [43]. It indicated that CX has a hepatoprotective effect, but the therapeutic effect of CX on NAFLD remains unclear. We first evaluated the activity of CX. In vitro NAFLD cell model was established by FFA, and it was found that CX inhibited lipid droplet formation in AML12 cells. Then, we further studied the potential mechanism of CX on NAFLD by network pharmacological focused on CX. Different from previous reports related to CX and network pharmacology [26,27,36], in addition to using TCMSP platform to screen the active components of CX, we also collected the bioactive CX components reported in the literature to supplement the active components of CX. Protein-protein intersection results shows that ALB, TNF, PPARG, PPARA, IL1B, VEGFA, CYP3A4, CYP2E1 and ESR1 was the hub genes of CX treated NAFLD. ALB is the most highly expressed protein in liver. It is specifically expressed in mature liver cells and plays a role in maintaining body nutrition and osmotic pressure. It was found that high-fat diets decreased the expression of ALB protein in the liver [44,45]. APOA1 expression and molecular transport promoted by PPARA, PPARD and PPARG could form the PPAR-APOA1 signaling pathway, in the pathogenesis of NAFLD it playing a beneficial role [46]. NAFLD covers a range of diseases, from steatosis to NASH [47]. In NASH, inflammatory cytokines TNF and IL-1B plays an important role [48]. VEGFA is an endothelial growth factor that promotes endothelial cell proliferation, migration and differentiation. The degree of angiogenesis in NASH was proportional to the degree of liver fibrosis compared with FAL alone or healthy liver [49,50]. The cytochrome P450 protein family mainly exists in the endoplasmic reticulum of hepatocytes, and plays a role in promoting drug metabolism. Increased activity of CYP2E1 and decreased activity of CYP3A4 were observed in rats with fatty liver induced by high fat diet [51]. In another experiment, however, the activity of CYP3A and CYP2E1 decreased by 30% and 20%, respectively, in rats fed a high-gluc-idic fat-free diet that increased liver fat by 235% [52]. It is may be due to the different induction and progression of the disease. Liver is the target organ of sex hormones and can express ESR1 and androgen receptor. The major differences in liver disease incidence were related to differences in sex hormone levels and expression of sex hormone-specific genes [53]. Therefore, these 10 hub genes play an important role in the treatment of NAFLD.
Enrichment analysis of GO showed that the main biological process of CX target was monocarboxylic acid metabolic process. Extramitochondrial fatty acids regulate extramitochondrial fatty acid oxidation and fatty acid binding protein pathways in liver [54]. The main molecular function of CX target was oxidoreductase activity. NAFLD is closely related to oxidative stress. In liver, lipid metabolism disorders lead to lipid accumulation, affecting different reactive oxygen species generators, including NADPH oxidase, mitochondria and endoplasmic reticulum [55]. The KEGG pathway enrichment including metabolism of xenobiotics by cytochrome P450, insulin resistance and fatty acid degradation. Cytochrome P450-derived eicosanoids have been shown to play an important protective role in a variety of inflammatory liver diseases [56]. The accumulation of FFA in liver cells may be caused by insulin resistance in the context of NAFLD, resulting in lipotoxicity, inflammation, and fibrosis [57].
Novel compounds of therapeutic significance can be identified using molecular docking. And molecular docking can predict ligand-target interactions at the molecular level [58]. We found that CX active compounds have high affinity for target. 1/3 of the compounds had the highest binding capacity with TNF. It is may related to the arginine, proline, alanine, isoleucine, lysine, and tyrosine residues on TNF can bind to the compound. This is of great significance for preventing NAFLD from developing into NASH. 90% of hub genes had the highest binding capacity with folic acid. This may be due to its chemical structure including a pteridine ring, p-aminobenzoic acid, and one or more gamma-linked glutamate residues [59]. The carboxyl group of glutamate residues can be dehydrated and condensed with amino acid residues of receptor proteins to form peptide bonds. The atoms in the peptide bond are in resonance and are highly stable [60]. This may be one of the reasons why folic acid has higher binding energy. In the following experiments, we can consider folic acid content as an indicator to establish a method for extracting CX active components for verification.
KEGG results showed that CX could treat NAFLD through PPAR signaling pathway. PPARA, PPARG, TNF and IL-1B are also the hub genes we predicted for CX treatment of NAFLD. Therefore, use western blotting to detect the four proteins expression after CX treatment. The results showed that CX could increase PPARA and PPARG protein expression while decrease the TNF and IL-1B protein expression.

Conclusions
In this study, we first found that CX can inhibit the increase of AML12 cell lipid droplets caused by FFA. Then, using network pharmacology with molecular docking methods to determine the molecular mechanism of CX therapy for NAFLD. Finally, using western blotting to verify the effect of CX on AML12 protein expression. Our study suggests that CX has a therapeutic effect on NAFLD.