A network pharmacology‐based approach to explore the effects of Chaihu Shugan powder on a non‐alcoholic fatty liver rat model through nuclear receptors

Abstract The pathogenesis of non‐alcoholic fatty liver disease (NAFLD) is still not fully understood, and currently, no effective pharmacotherapy is available. Nuclear receptors (NRs) are important biological participants in NAFLD that exhibit great therapeutic potential. Chaihu Shugan powder (CSP) is a traditional Chinese medicine (TCM) formula that has a wide therapeutic spectrum including NAFLD, but the effective components and functional mechanisms of CSP are unclear. We adopted a network pharmacology approach using multiple databases for Gene Ontology (GO) enrichment analysis and the molecular complex detection (MCODE) method for a protein‐protein interaction (PPI) analysis, and we used molecular docking method to screen the NR targets and determine the corresponding CSP components. The screening results were validated through a NAFLD rat model that was used to explain the possible relationship between CSP and NAFLD. Finally, we screened PPARγ, FXR, PPARα, RARα and PPARδ as target genes and quercetin, kaempferol, naringenin, isorhamnetin and nobiletin as target compounds. The five components were detected through high‐performance liquid chromatography‐mass spectrometry (HPLC‐MS), the results of which aligned with the docking experiments of PPARγ, PPARα and PPARδ. After CSP intervention, the NAFLD rat model showed ameliorated effects in terms of bodyweight, hepatic histopathology, and serum and liver lipids, and the mRNA levels of PPARγ, FXR, PPARα and RARα were significantly changed. The results from this study indicate that CSP exhibits healing effects in an NAFLD model and that the network pharmacology approach to screening NR targets and determining the corresponding CSP components is a practical strategy for explaining the mechanism by which CSP ameliorates NAFLD.


| INTRODUC TI ON
Non-alcoholic fatty liver disease (NAFLD) is a metabolic stress injury to the liver that is closely related to insulin resistance and is based on genetic susceptibility without a history of excessive alcohol consumption. The NAFLD disease spectrum includes non-alcoholic simple fatty liver (NAFL), non-alcoholic steatohepatitis (NASH), and cirrhosis and hepatocellular carcinoma (HCC). 1,2 NAFLD has an estimated global prevalence of 24%, indicating that it is the most common cause of liver disease worldwide. 3 Currently, there is no approved pharmacotherapy for NAFLD. Lifestyle modification, control of comorbid metabolic illness and some pharmacological options, such as vitamin E, pioglitazone and pentoxifylline, are the treatment options. 4 Nuclear receptors (NRs) are ligand-activated transcription factors, of which 48 members constitute the largest family of transcription regulatory mechanisms. 5 In the progression of NAFLD, NRs are important biological participants involved in insulin resistance and fatty acid (FA) synthesis,triglyceride (TG), total cholesterol (TC) and glucose (GLU) metabolism; bile acid homeostasis; drug detoxification; inflammation; regeneration; fibrosis; tumour formation; etc. 6,7 Some NRs are sensitive to lipid metabolic processes which may play important roles in NAFLD. 7  Current studies have reported that CSP has therapeutic potential for the treatment of NAFLD because it ameliorates insulin resistance, lipid peroxidation, steatohepatitis, etc, [8][9][10] but the chemical components and functional mechanisms of CSP have not been fully studied. Therefore, an exploration into the underlying mechanism of CSP in NAFLD is necessary.
The network pharmacology approach is a big data integration method based on a large number of database resources and statistical algorithms that are used to observe the synergy of multiple components, targets and mechanisms of disease and can be applied to broad-perspective analyses of TCM drugs. [11][12][13] As with other natural medicines, CSP contains abundant chemical components that indispensably correspond to many target genes and biological mechanisms. 14 Therefore, it is very difficult to fully clarify the effective components and functional target range and specifically to elucidate the therapeutic mechanisms of CSP. Therefore, we adopted a network pharmacology approach to screen NR targets and corresponding components of CSP and conducted verification experiments on a rat model of NAFLD to further elucidate the relationship between CSP and NAFLD.

| Components and target gene screening of CSP
The TCMSP, BATMAN-TCM, TCMID and TCMGeneDIT databases were used to identify the chemical components of CSP. DrugBank and PubChem were used to search for the information of chemical compounds and target genes of CSP. NCBI, GeneCards and DisGeNET were used to search for the therapeutic target genes of NAFLD. UniProt and Ensembl were used to convert protein names to gene names (ENSG identifiers). The target genes at the intersection of CSP and NAFLD (C&N) were regarded as potential target genes through which CSP influences NAFLD, and the corresponding chemical compounds of the C&N candidates were thought to be possible therapeutic components that affect NAFLD. Then, the NRs were screened out from the C&N candidates and regarded as potential therapeutic NR targets of CSP in NAFLD. Similarly, the CSP compounds corresponding to these NRs were regarded as possible therapeutic agents for NAFLD treatment.

| High-potential target screening of CSP
As there are differences in data sources, data storage and statistical algorithms, to identify high-potential NR targets, we integrated multiple databases and conducted Gene Ontology (GO) enrichment analysis using a comprehensive approach, which included the biological process (BP), molecular function (MF) and cellular component (CC) categories. From the results of the BP, MF and CC analysis of each database, we selected the top thirty clusters based on statistical significance and then ranked the total NR frequency in all these clusters in descending order.
Finally, the top five NRs were selected as the high-potential target genes of CSP. For GO enrichment analysis, the smaller P-value of a cluster was considered to indicate greater biological significance of CSP treatment of NAFLD; hence, the NRs included in these clusters were regarded as more effective potential therapeutic targets. In addition, we also integrated multiple databases to conduct a signalling pathway enrichment analysis of the C&N groups and then selected the top thirty pathways derived from each database based on statistical significance and determined the expression of the top 5 NRs in each pathway. Protein-protein interaction (PPI) analysis of the C&N data was conducted by STRING, and the highly interconnected regions in the PPI network were found through the molecular complex detection algorithm (MCODE) plugin in Cytoscape with the parameters set to the default values. Cytoscape was also used to visualize the topological structure of the interaction network. Each node in the network represented a target or compound, and each line represented the connection of the target to the compound or the target to the target.

| High-potential component screening of CSP
As the top 5 NRs may function through more than one compound of CSP, it is necessary to conduct an ADME pharmacodynamic screen based on oral bioavailability (OB) > 30% and druglikeness (DL) > 0.18. After the compounds that corresponded to top 5 NRs were ranked by frequency, the top 5 compounds were selected as the final target CSP components. Similar to the top 5 NRs, the top 5 compounds could be considered high-potential therapeutic components of CSP for treating NAFLD, presupposing that they function through the top 5 NRs. at 80°C and cooled to room temperature before gavage (the CSP calories were calculated based on 4 kcal/g carbohydrate). The rats were isolated and fed for one week while they adapted to the environment and then randomly divided again into a normal control (NC) group (n = 10), high-fat diet (HFD) control group (n = 10) and CSP intervention (CSP) group (n = 10) according to the random number table method. There were 3-4 rats fed in one cage for a total of 8 weeks. The NC group rats were fed freely with normal chow, and the HFD and CSP group rats were fed freely with high-fat chow. All the rats had free access to drinking water. In addition, the CSP group rats were given a weight-corresponding dosage of CSP (10 mL/kg), and the NC and HFD group rats were given the same dosage of deionized water. The room temperature was 22 (±2)°C, the relative humidity was 55 (±5)%, light and dark conditions were alternated every 12 hours, and the rat cages were cleaned every 3 days.

| Biochemical detection
After 8 weeks of feeding, the rats were fasted for 12 hours, but had free access to drinking water, and then were intraperitoneally anesthetized with pentobarbital sodium (40 mg/kg). After the rats were unconscious, abdominal aortic blood was collected and placed at room temperature for 30 minutes, and then, the serum was col-

| Quantitative real-time PCR
Total RNA was extracted from the liver tissue using TRIzol reagent (15596026, Thermo Fisher). Reverse transcription reactions were conducted according to the instructions for the PrimeScript reagent kit (RR047A, TaKaRa). qPCR amplification was conducted using a 25 μL reaction system with a PCR amplifier (iQ5, Bio-Rad) according to the instructions of the qPCR reagent kit (RR820A, TaKaRa).
The relative quantification of RNA, in number of folds, was calculated by the 2 −ΔΔCt method. The gene sequences were acquired from GenBank, and GAPDH was used as the internal reference. The primers were synthesized by Shanghai General Biotech Co., Ltd., (primer information is listed in Table S2).

| Automatic Western blotting
The sample protein was extracted from liver tissues after the cells were lysed by RIPA buffer (89901, Thermo Scientific) at 4°C for 20 minutes and centrifuged for 10 minutes (4°C, 11180.0g). The secondary antibodies were anti-PPARγ (ab209350, Abcam), anti-PPARα (ab24509, Abcam), anti-FXR (NBP2-16550, Novus), anti-RARα (ab254098, OriGene) and anti-PPARδ (ab23673, Abcam). The protein was quantified completely by automatic Western blotting through a Wes automated system (ProteinSimple). The dynamic linear range of the standard protein was diluted in a gradient from 0.125 to 4 mg/mL. The Western blotting process was conducted according to the instructions for the Wes system, which differs from the traditional method. The machine was run with default parameters, and imaging and analysis were performed using Compass software (ProteinSimple).

| HPLC-MS
To detect the top 5 compounds in CSP, the CSP granules, which were from the same batch produced for rat feed, were precisely weighed according to the ratio described in the Animals and treatments sec- formic acid water-acetonitrile (96:4), the column temperature was 40°C, the flow rate was 450 μL/min, and the injection volume was 5 μL. Gas-assisted electrospray ionization (ESI) was used for mass spectrometry (MS), and the detection method was based on multiion reaction measurement (MRM).

| Molecular docking
PubChem and PDB were used to find the chemical and conformational information of the relevant proteins and small-molecule compounds.
The keywords associated with the proteins were determined, and the candidates were selected according to the following criteria: a crystal resolution of <3 Å for protein structures obtained by the X-ray crystal diffraction method, a well-defined type of protein and a protein with small-molecule ligands. The AutoTools software was used to remove the redundant protein chains, ligands and water molecules with hydrogenation before running docking experiments. The AutoGrid software was used to calculate the energy lattice points with the grid box coordinates of 20 × 20 × 20. AutoDock Vina was used to simulate the docking condition between proteins and small molecules. The Schrodinger software was used to analyse the preferential conformation and map the simulation (databases are listed in Table S1).

| Statistical method
In addition to the built-in statistical algorithms in the databases and/or software used for network pharmacology, all the statistical calculations of the one-way analysis of variance (ANOVA) were performed using GraphPad Prism 8.0 software. The experimental results are presented as the means ± standard deviation, and P < .05 indicates significance.  Radix Bupleuri was the herb that contributed to the highest proportion of collected compounds and targets (Figures 2 and 3).   36 interactions and a score of 4.235; the cluster containing FXR had 4 targets, 5 interactions and a score of 3.333; and the cluster containing RARα had 3 targets, 3 interactions and a score of 3.000 ( the cluster score represents the core density of the node and the topologically adjacent nodes; a higher score represents a more concentrated cluster), suggesting that the top 5 NRs may function through a multitarget synergistic mechanism with C&N targets (Figure 4).

| Rat general status
The average daily calorie intake per rat indicated that the NC group rats ingested the most and the CSP group rats ingested the least, but the differences were not significant. The bodyweight and Lee index of the rats in the HFD group were significantly higher than those of the NC group and CSP group, and for the CSP group, the F I G U R E 2 The amount of targetcompound interactions for each CSP herb. A, The CSP compound amount in each herb. B, The CSP targets amount in each herb. C, The C&N targets amount in each herb. D, The C&N compounds amount in each herb. E, The target amount of NRs among the C&N network in each herb. F, The compound amount of NRs among the C&N network in each herb. Above all the UpSet plots, the horizontal bars represent the total quantity and proportion of each herb, and the vertical bars represent the independent quantity of each herb and the quantities that intersect with that of the other herbs bodyweight and Lee index were slightly lower than those of the NC group, but the differences were not significant. The variable tendencies in liver weight, perirenal fat weight and peritesticular fat weight in rats were relatively similar among that rats; in the HFD group, these weights were significantly higher than they were in the NC or CSP group, and in the CSP group, they were higher than those in the NC group, but only the difference in liver weight was significant Figure 6

| Pathological observation of the liver
In general, the liver colour, blood flow, morphology, and the texture and tactility condition of liver tissues in the NC group were better than those in the HFD group or CSP group, whereas they were better in the CSP group than they were in the HFD group. The results of HE and Oil red O staining showed that the histological structure, cell morphology and cytoplasmic structures of the NC group were significantly more normal than those of the HFD and CSP groups, and only few lipid droplets were observed. In the HFD group, many deposited lipid droplets and an irregular cell morphology were observed. Compared with these parameters in the HFD group, the lipid droplet quantity was lower and the cell structure morphology was more regular in the CSP group. Under an electron microscope, for most of the mitochondria in the NC group, contours of the internal cristae were visible and a small amount of lipid was distributed in the cytoplasm. The HFD group had lipid droplets with a diameter of 2-3 μm in the cytoplasm. The nuclei were concentrated; the mitochondria were diluted, swollen and distorted; the internal cristae were blurred; the endoplasmic reticulum quantity was reduced; necrotic hepatocytes and blurred structures were apparent; and the cell membrane had almost disappeared. In the CSP group, the lipid droplets were observed in hepatocytes, but the diameters F I G U R E 3 The NR target-compound network of each herb in CSP. In the NR networks, the red nodes represent the targeted NRs, the blue and green nodes represent corresponding compounds, and the edges represent the connection between targets and compounds F I G U R E 4 GO enrichment and PPI analysis. A, GO enrichment analysis of BP, MF and CC functions. All the target genes in the GO enrichment analysis were significant (P < .05). B, PPI analysis. The edge represents the connection between two proteins, and the size of the node represents the total number of connecting edges of the lipids were mostly <1 μm. The mitochondria were swollen, but the internal cristae were still apparent, and the changes in the endoplasmic reticulum were significantly less dramatic than those in the HFD group. The lipid area calculation also showed that the NC group had significantly less lipid area than the HFD or CSP group, with that of the CSP group lower than that of the HFD group ( Figure 6).

| Biochemical parameters of the serum and liver homogenate
In addition to LDL-c, the other serum indicators were relatively consistent. The results of the NC and CSP groups indicated significantly lower levels than those in the HFD group, with some differences significant. According to the results from the liver homogenate, the TC and TG levels in the NC group were significantly lower than those in the HFD and CSP groups, whereas those in the HFD and CSP groups were relatively similar, with no significant differences ( Figure 7).

| The mRNA expression of the top 5 NRs in the liver
Compared with those in the NC group and CSP group, PPARγ and RARα were expressed at significantly higher levels in the HFD group with significance. FXR, PPARα and PPARδ were expressed at lower levels in the HFD group, but only the differences in FXR and PPARα were significant among groups. In contrast to those

| The protein expression of the top 5 NRs in the liver
FXR, PPARα and PPARδ were expressed at lower levels in the HFD group than in the NC group or CSP group, whereas the expression in the CSP group was more similar to that in the NC group than to that in the HFD group. The RARα expression level was higher in the HFD group and lower in the NC group and the CSP group. Generally, although no significance was found among the protein expression levels, the proteins were expressed following the approximate trend as the mRNA expression ( Figure 7).

| Network pharmacology and CSP
The network pharmacology method is essentially a combination of clustering algorithms and network topology. As the complex data interaction relationship is well presented by the visualized node interactions, this method is often used to analyse TCM herbs. [11][12][13][14] Based on the 730 compounds and 917 corresponding target genes found, CSP has a large quantity of chemical components and corresponding targets. If the compounds and the targets only had one-to-one relationships, then they would form nearly a thousand pairs of interacting connections. It is widely known that a target generally functions through several synergistic mechanisms with the related targets and that there may be many interactions between the targets. Therefore, it was not surprising that 50 928 GO enrichment clusters, 6435 pairs of PPI connections and 2791 pathways were found. In addition, most of the CSP compounds have pharmacological functions through multiple NAFLD targets. Therefore, the whole therapeutic process from compounds to targets and finally to NAFLD was composed of a very large and complex hierarchical crossover compound-target network ( Figure 9). With the current technology level, verification of the whole potential target network is difficult, but pinpointing the factors with high research potential is possible. For example, significant differences could be used as a screening condition after enrichment to determine the clusters that were more significant. Because of the large number of screened targets, the biological mechanism covered almost all aspects of the development of NAFLD, and biological functions could also be used as screening conditions, which contribute to a more precise focus on the highly significant targets with specific biochemical functions.

| The effects of CSP on NAFLD
The development of NAFLD includes a wide pathological process of mechanisms involving insulin resistance, oxidative stress, endoplasmic reticulum stress, apoptosis, lipotoxicity, inflammation and the microbiome. 15 CSP is widely used in China and other East Asian countries, and the treatment spectrum involves NAFLD, 8 NASH,9 functional dyspepsia, 16 depression, 17 atherosclerosis, 18 etc, although the specific components and targets of CSP that are most important in NAFLD are currently unclear. All the components and targets of CSP revealed in this study seem to cover most of the biological mechanisms of NAFLD,therefore, we inferred that CSP has a broad impact on many factors closely related to NAFLD. Obesity and lipid accumulation are important cofactors in the development of NAFLD. 15,19,20 The alleviation of overweight and lipid accumulation of the liver, perirenal and peritesticular tissue were apparently found in the CSP when compared to those in the HFD group. According to the pathological observations, fewer lipid droplets were found in the CSP group than in the HFD group. The pathological results also F I G U R E 6 Rat general status and pathological observations. A, Rat general status. The bodyweight in the HFD group was higher than it was in the NC group (P < .0001) and the CSP group (P < .0001); the liver weight in the HFD group was higher than it was in the NC group (P < .0001) and CSP group (P = .0059), and it was lower for the NC group than it was for the CSP group (P < .0010); the perirenal fat weight in the HFD was higher than it was in the NC group (P = .0010) and CSP group (P = .0025); the peritesticular fat weight in the HFD group was higher than that it was in the NC group (P = .0116) and CSP group (P = .0309); and the Lee index of the HFD was higher than it was for the NC group (P < .0001) and CSP group (P < .0001) (n = 10). B, Morphology of the liver. C, HE-stained liver tissue. D, Oil red O-stained liver tissue. E, Transmission electron microscopy observation of the liver tissue. F, Oil red O-stained area. The NC group had a lower lipid area than did the HFD group (P < .0001) and CSP group (P < .0001), and the lipid area for the CSP group was lower than it was for the HFD group (P = .0032) (n = 6) showed that the lipid area of the HFD group was 22.86%, while that of the CSP group was only 17.96%, which was 4.9% lower. According to relevant reports, the criterion of at least 5% of the hepatocytes showing steatosis, as indicated by large cysts, large pure droplets or mixtures of small droplets and large droplets, can be used as the basis for an NAFLD diagnosis. 21 Therefore, in our study, the NAFLD rat model was successfully established, with obviously better liver adipose pathology amelioration in the CSP group. In addition, in the HFD group and CSP group, no obvious inflammations were observed in the liver tissue sections, which indicated that the model was simple fatty liver without a sign of NASH. The biochemical results indicated that the AST, ALT and ALP levels in the CSP group were significantly lower than those in the HFD group. Currently, no recognized biochemical diagnostic indicators are available for NAFLD, AST, ALT and ALP can only be used only as the auxiliary diagnostic indicators for NAFLD or for differential diagnosis with NASH and fibrosis. 3,22,23 Abnormal blood lipids are characterized by an increase in TG and LDL-C, and a decrease in HDL-C levels. 3,24,25 GLU is also an important biochemical factor in NAFLD and an important substrate in the synthesis of lipids. 26,27 In the present study, the serum levels of TG, TC, HDL-C and GLU in the CSP group were significantly lower than those in the HFD group. The biochemical and pathological results indicate that CSP has the prospect of lowering bodyweight and liver and blood lipids in NAFLD patients. can inhibit the differentiation of adipocytes that is induced by rosiglitazone, a PPARγ agonist, and reduce the development of obesity and the fatty degeneration of the liver caused by a high-fat diet and leptin deficiency. 43 Nobiletin can prevent the inflammation, insulin resistance, dyslipidaemia and fatty liver caused by a high-fat diet. 44 According to their functions, the top 5 compounds may pharmacologically act on NAFLD and may function through PPARs. Although these compounds can bind with PPARα, PPARδ and PPARγ well, we cannot overlook the fact that the screened compounds and target F I G U R E 7 Biochemical measurements and the mRNA and protein expression of the top 5 NRs. A, The results of biochemical parameter measurements. Among the serum indicators, ALT in the HFD group was higher than it was in the NC group (P = .0004) and CSP group (P = .0222); the AST in the HFD group was higher than it was in the NC group (P = .0006) and CSP group (P = .0149); ALP in the HFD group was higher than it was in the NC group (P < .0001) and CSP group (P = .0013), and it was lower in the NC group than in the CSP group (P = .0017); in the HFD group, the TC level was higher than it was in the NC group (P < .0001) and CSP group (P = .0003), and in the NC group, it was lower than it was in the CSP group (P < .0001); in the HFD group, the TG level was higher than it was in the NC group (P = .0011) and CSP group (P = .0476); in the HFD group, the GLU level was higher than it was in the NC group (P < .0001) and CSP group (P = .0002); and in the HFD group, the HDL-c level was higher than it was in the NC group (P < .0001) and CSP group (P = .0022), and in the NC group, it was lower than it was in the CSP group (P < .0001) (n = 10). Among the liver indicators, in the NC group, the TC level was lower than it was in the HFD group (P = .0276) and CSP group (P = .0427); in the NC group, the TG level was lower than it was in the HFD group (P = .0427) and CSP group (P = .0448) (n = 6). B, The mRNA expression of the top 5 NRs in the liver. PPARγ in the HFD group was expressed at a higher level than it was in the NC group (P = .0072) and CSP group (P = .0489); FXR in the HFD group was expressed at a lower level than it was in the NC group (P = .0222) and CSP group (P = .0331); PPARα in HFD group expressed lower than NC group (P = .0143) and CSP group (P = .0082); and RARα in the HFD group was expressed at a higher level than it was in the NC group (P = .0001) and CSP group (P = .0012) (n = 10). C, The expression of the top 5 NRs in the liver was measured by automatic Western blot analysis

RE S E ARCH
In this study, CSP showed positive effects on some factors closely related to NAFLD, such as bodyweight and lipid accumulation.
The network pharmacology approach with multiple databases, as adopted in this study, was useful in screening potential targets, and the screened compounds and genes found exhibited therapeutic potential for NAFLD. The final five compounds that were screened out are flavonoids, which can combine with the screened-out target genes PPARγ, PPARα and PPARδ, and thus show potential as drugs for NAFLD. This study also provided a methodological exploration based on network pharmacology. Although we tried to rescreen the C&N targets through enrichment analysis, the target frequency presented in clusters, and ADME screening, the biological function classification was still necessary to identify targets precisely. In vivo and in vitro experiments need to be conducted in the next step to investigate the specific pharmacological effect of the top 5 compounds on NAFLD.
In addition, various high-throughput screening methods, such as sequencing and genomics or proteomics analyses, can be combined with

CO N FLI C T O F I NTE R E S T
The authors declared that they have no conflicts of interest related to this study. Zhang provided suggestion for revising the manuscript. All authors read and approved the final manuscript.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data used to support the findings of this study are available from the corresponding author upon reasonable request.