Study on Hepatotoxicity of Rhubarb Based on Metabolomics and Network Pharmacology

Abstract

Background

Rhubarb, as a traditional Chinese medicine, is the preferred drug for the treatment of stagnation and constipation in clinical practice. It has been reported that rhubarb possesses hepatotoxicity, but its mechanism in vivo is still unclear.

Methods

In this study, the chemical components in rhubarb were identified based on UPLC-Q-TOF/MS combined with data postprocessing technology. The metabolic biomarkers obtained through metabolomics technology were related to rhubarb-induced hepatotoxicity. Furthermore, the potential targets of rhubarb-induced hepatotoxicity were obtained by network pharmacology involving the above components and metabolites. Meanwhile, GO gene enrichment analysis and KEGG pathway analysis were performed on the common targets.

Results

Twenty-eight components in rhubarb were identified based on UPLC-Q-TOF/MS, and 242 targets related to rhubarb ingredients were predicted. Nine metabolic biomarkers obtained through metabolomics technology were closely related to rhubarb-induced hepatotoxicity, and 282 targets of metabolites were predicted. Among them, the levels of 4 metabolites, namely dynorphin B (10–13), cervonoyl ethanolamide, lysoPE (18:2), and 3-hydroxyphenyl 2-hydroxybenzoate, significantly increased, while the levels of 5 metabolites, namely dopamine, biopterin, choline, coenzyme Q9 and P1, P4-bis (5ʹ-uridyl) tetraphosphate significantly decreased. In addition, 166 potential targets of rhubarb-induced hepatotoxicity were obtained by network pharmacology. The KEGG pathway analysis was performed on the common targets to obtain 46 associated signaling pathways.

Conclusion

These data suggested that rhubarb may cause liver toxicity due to its action on dopamine D1 receptor (DRD1), dopamine D2 receptor (DRD2), phosphodiesterase 4B (PDE4B), vanilloid receptor (TRPV1); transient receptor potential cation channel subfamily M member 8 (TRPM8), prostanoid EP2 receptor (PTGER2), acetylcholinesterase (ACHE), muscarinic acetylcholine receptor M3 (CHRM3) through the cAMP signaling pathway, cholinergic synapses, and inflammatory mediators to regulate TRP channels. Metabolomics technology and network pharmacology were integrated to explore rhubarb hepatotoxicity to promote the reasonable clinical application of rhubarb.

Keywords:

• rhubarb
• hepatotoxicity
• metabolomics
• network pharmacology

Introduction

Rhubarb is the dried root and rhizome of Rheum palmatum (Rheum palmatum L.), Tanghete rhubarb (Rheum tanguticum Maxim. ex Balf.) or medicinal rhubarb (Rheum officinale Bail.).¹ Rhubarb has the following effects: such as purging and attacking accumulation, clearing away fire and heat, removing dampness and relieving yellow, removing blood stasis and clearing menstruation. At the same time, rhubarb is the first choice for the treatment of stagnation and constipation in clinical practice.² Free anthraquinones and bound anthraquinones are the main effective ingredients of rhubarb. Free anthraquinones, including emodin, rhein, and aloe-emodin, are the antibacterial ingredients of rhubarb, and conjugated anthraquinones are the main diarrhea components of rhubarb, among which the bisanthrone glycoside has the strongest effect.^3,4 In clinical practice, rhubarb is widely used to treat constipation, acute pancreatitis, renal failure, and hyperlipidemia.^5–7 However, the misuse and abuse of rhubarb have endlessly emerged on account of the widespread clinical application of rhubarb in recent years. It has been found that rhubarb exhibits strong adverse reactions when used in large doses, among which liver toxicity of rhubarb is more prominent to become an unavoidable symptom.^8–10 However, the mechanism of rhubarb-induced hepatotoxicity and its correlative components are still unclear. Therefore, it is necessary to study rhubarb hepatotoxicity with accurate and rapid methods.

Metabolomics is a discipline that has grown up after proteomics and genomics.^11,12 The metabolic network of biological systems can be derived by observing the dynamic changes in metabolites before and after stimulation or perturbation of organism systems on the basis of various analytical tools.^13,14 Metabolomics is widely used in various scientific fields due to its high pathway integrity, high sensitivity, and high selectivity.¹⁵ Because the research concept of metabolomics has a certain similarity with the overall view of traditional Chinese medicine (TCM), it was rapidly introduced into the field of TCM research and then quickly infiltrated into the investigation of TCM toxicology, becoming one of the main means of screening toxicity markers and studying the mechanism of toxic effects.^16,17 Metabolomics is an important analysis method for studying the therapeutic effect and toxicity of drug.^18–20 Therefore, nontargeted metabolomics technology in the study is used to analyze the changes in endogenous metabolic markers after rhubarb-induced liver injury. Network pharmacology is a discipline that explains the process of disease development based on the principles of systems biology, and it uses the overall view of network databases to study the interaction between drugs and the body.²¹ Network pharmacology creates a novel model to rapidly identify the complex network relationship between drug targets and disease targets, and it has an important application value in researching the pathogenesis and therapeutic targets of complex diseases.²² Network pharmacology provides a new perspective for the analysis of drug effects. The connection and relationship of nodes in the biological network can be used to analyze the characteristics of interconnections, and the mechanism between drug effects and body changes can be further clarified.²³ Network pharmacology has made some progress in the study of combination drugs, diseases, and toxicity of TCM. Therefore, network pharmacology can be applied to the study of rhubarb hepatotoxicity, and the relevant targets of rhubarb-induced hepatotoxicity can be quickly obtained through virtual screening and analysis.

Metabolomics and network pharmacology are both important parts of systems biology.²⁴ Metabolomics explores the toxic effects of TCM through changes in the body’s endogenous metabolic molecules and related metabolic pathways. Network pharmacology then explores the toxic effects of TCM by screening ingredients and their targets, constructing a network of “component-target-pathways”, and it systematically and comprehensively reveals the toxic mechanism of TCM.²⁵ Drugs work by interacting with various molecular targets. Because most metabolites and drugs interact with multiple target proteins, the methods to predict new targets with existing endogenous metabolites and small molecule compounds have high biological and pharmacological value.^26,27 Network pharmacology links the toxicology of Chinese medicine with metabolomics, genes, and proteomics in modern research.²⁸ The combination of network pharmacology and metabolomics will better discover and describe the multicomponent, multitarget, and multipathway action mechanisms of Chinese medicine, which is conducive to the exploration of the action and toxicity mechanisms of Chinese medicine.²⁹ This is also a new idea for the modernization of Chinese medicine research, resulting in better development prospects.

In light of the above research status, we explored a new perspective combining metabolomics and network pharmacology to study rhubarb-induced hepatotoxicity. Based on the UPLC-Q-TOF/MS technology platform, metabolomics technology was used to preliminarily study changes in the metabolic markers of rhubarb-induced liver toxicity. Furthermore, the common targets of liver injury in the Comparative Toxicogenomics Database (CTD) and combined targets, which included rhubarb component targets and hepatotoxic metabolite targets, were recognized as potential targets for rhubarb-induced hepatotoxicity. The pathways involved in rhubarb-induced hepatotoxicity in vivo were predicted through potential targets, and the mechanisms of which were researched and analyzed. This research will provide a safety warning for rhubarb in future clinical applications.

Materials and Methods

Experimental Instruments and Reagents

An electric heating jacket (Beijing Zhongxing Weiye Instrument Co., Ltd.), rotary evaporator (Swiss Buchi Co., Ltd.), a circulating water multipurpose vacuum pump (Yuhua Instrument Co., Ltd.), a low-temperature high-speed centrifuge (Changsha Xiangyi Centrifuge Instrument Co., Ltd.), 1/10,000 balance (Shanghai Balance Instrument Factory), small vortex meter (Haimen Qilinbeier Instrument Manufacturing Co., Ltd.), Waters ACQUITY UPLC (Waters, USA), Waters Xevo TQD high-resolution mass spectrometer (Waters, USA Company), and ACQUITY UPLC BEH C18 column (2.1 × 100 mm, 1.7 μm, Waters Corporation, USA) were used in the present study.

Rhubarb (Anguo Pingchuan Pharmaceutical Co., Ltd., production batch number: 18080035), chrysophanol (Chengdu Herbpurify Co., Ltd.), aloe-emodin (Chengdu Herbpurify Co., Ltd.), rhein (Chengdu Herbpurify Co., Ltd.), emodin (Chengdu Herbpurify Co., Ltd.), 95% ethanol (Tianjin Concord Technology Co., Ltd.), physiological saline (Shijiazhuang Fourth Pharmaceutical Co., Ltd.), distilled water (Guangzhou Watsons Food and Beverage Co., Ltd.), chloraldehyde hydrate (Shanghai McLean Biochemical Technology Co., Ltd.), tissue fixative (Beijing Yili Fine Chemicals Co., Ltd.), chromatographically pure acetonitrile (Sigma-Aldrich), chromatographically pure methanol (US ROE), chromatographically pure formic acid (Sigma-Aldrich), and analytical pure ethanol (US ROE company) were used in the present study.

Preparation of Rhubarb Extract

The rhubarb medicinal material was crushed into a coarse powder through a 10-mesh sieve. Ten times the amount of 95% ethanol was added to the powder, which was then heated to reflux for 1.5 h and filtered with gauze. Eight times the amount of 95% ethanol was added to the medicine residue, which was then heated to reflux for 45 min and filtered with gauze. After the ethanol in the filter residue was evaporated, ten times the amount of distilled water was added, and the mixture was heated to reflux for 1 h. After filtering the residue with gauze, the filtrate and the alcohol extract were combined and concentrated under reduced pressure with a rotary evaporator to obtain the rhubarb extract.³⁰ Distilled water was added to the rhubarb extract to prepare a 1 g/mL drug solution (each 1 mL of the drug solution contained 1 g of crude drug). This drug solution was used for intragastric administration in rats. The prepared rhubarb extract was stored in 4 °C refrigerator.

Rhubarb Ingredient Detection

The abovementioned rhubarb extract was added to distilled water to prepare a 0.1 g/mL rhubarb extract (that is, 0.1 g crude drug per 1 mL of drug solution). The rhubarb extract was centrifuged at 13,000 rpm for 15 min at 4 °C, and the supernatant was collected through a 0.22 μm microporous membrane. All standard samples (0.5 mg/mL) were prepared in methanol aqueous solution through a 0.22 μm microporous membrane. Then, 200 μL of supernatant was transferred into the injection vial, and UPLC-Q-TOF/MS liquid quality sampling analysis was performed. The chromatographic conditions and mass spectrometry conditions are provided in the Supplementary Information.

Animal Experiment

Clean-grade male Wistar rats (Beijing Witton Lihua Experimental Animal Technology Co., Ltd.), weighing 200 ± 20 g, were kept in the Animal Laboratory of Tianjin Institute of Physical Education. After the rats were purchased, they were placed in a changing environment for 12 h day and night and then reared. The temperature in the cultivation room was 23 ± 2 °C, and the humidity was 35 ± 5%. The study was approved by the Animal Ethics Committee of Tianjin University of Traditional Chinese Medicine under permit number TCM-2012-078F01). All procedures were conducted in accordance with the Chinese national legislation and local guidelines. Rats were randomLy divided into the following four groups (n = 8): blank control group administered normal saline; low-dose group administered total rhubarb extract (1.5 g/kg/d); medium-dose group administered total rhubarb extract (6 g/kg/d); and high-dose group administered total rhubarb extract (12 g/kg/d). In animal experiments, all administration methods were gavage, and the rats were administered continuously for 20 days. The blank control group was given the same volume of saline as the administration group. The low dosage was calculated according to the maximum dosage of rhubarb in the “Chinese Pharmacopoeia”. In addition, dosing was performed according to the dosing schedule shown in Table 1.

Table 1 Experimental Groups and Dosage of Each Group

Groups	N	Dosage	Mode of Administration	Administration Time
Control group	8	10 mL/kg/d	Continuous saline/Gavage	20 days
Low-dose group	8	1.5 g/kg/d	Continuous administration/Gavage	20 days
Medium-dose group	8	6 g/kg/d	Continuous administration/Gavage	20 days
High-dose group	8	12 g/kg/d	Continuous administration/Gavage	20 days

Collection and Pretreatment of Biological Samples

Rats in each group fasted for 12 h after the last administration time, and whole blood was collected. Rats were first anesthetized with a 10% chloral hydrate solution, and blood was taken from the abdominal aorta after anesthesia.³¹ Approximately 5 mL of whole blood was collected, and 1.5 mL of the whole blood sample was placed in a heparinized test tube, which was slowly rotated to ensure even mixing with heparin. The supernatant was centrifuged at 3000 rpm and 4 °C for 12 min, and the supernatant was centrifuged again at 3500 rpm and 4 °C for 8 min. The plasma sample was obtained and stored in −80 °C refrigerator for metabolomics research. The remaining 3 mL of whole blood was placed in a centrifuge tube, and the supernatant was centrifuged at 3000 rpm and 4 °C for 15 min. After collecting blood from each group, rats were perfused with physiological saline.³² The liver tissues of the rats were removed, and the blood stains on the surface were washed with physiological saline and then soaked in tissue fixation solution for pathological section analysis.

Plasma Metabolomics

Plasma samples were frozen at −80 °C and thawed at room temperature. Plasma (100 μL) was placed in an Ostro sample preparation plate, and 300 μL of acetonitrile was added to each well. The samples were mixed vigorously three times and allowed to stand at low temperature for 5 min. Then, 2–8 in-Hg pressure suction filtration was used, and suction filtration was completed. The solution was then collected from the plate and placed in a centrifuge tube. The solution was dried with nitrogen, and 100 μL of ultrapure water was added to reconstitute the sample, which was then vortexed for 1 min, sonicated for 15 min in an ice water bath, and centrifuged at 13,000 rpm for 15 min at 4 °C. Supernatant (200 μL) was transferred into the injection vial, and UPLC-Q-TOF/MS liquid quality sampling analysis was performed.³³ The details for UPLC-Q-TOF/MS analysis conditions and mass spectrometry conditions are provided in the Supplementary Information.

Plasma samples (approximately 50 μL) were obtained from each group to prepare quality control (QC) samples, and the QC sample contained the biological information of all samples.³⁴ The QC sample was prepared following the same procedure as the tissue sample to verify the experimental accuracy, repeatability, and stability. The detailed QC methodology surveys are provided in the Supplementary Information.

Biochemical Detection and Histopathology

The contents of alanine aminotransferase (ALT) and aspartate aminotransferase (AST) in the serum of experimental rats were determined using a fully automatic biochemical analyzer. The collected liver tissues were stained with HE, and the pathological manifestations of the rat liver were observed under an optical microscope. Optical microscope examination and histological evaluation were performed to determine the damage of rhubarb to liver tissue. For the HE staining procedure, paraffin sections were dewaxed with xylene, hydrated, stained with hematoxylin for 10 mins, differentiated, stained with eosin, dehydrated, rendered transparent with xylene, and glued with neutral resin.

Target Screening and Online Pharmacology

Target Prediction of Rhubarb Composition

In this study, rhubarb chemical components were used for target prediction. The above chemical composition structure was obtained from the PubChem (https://pubchem.ncbi.nlm.nih.gov/) database, and the SDF file was downloaded and saved. The downloaded SDF file was used to obtain the target point corresponding to the rhubarb component in the Swiss Target Prediction database (http://www.swisstargetprediction.ch/).³⁵ The targets were imported into Cytoscape software to construct the “rhubarb component-target” interaction network diagram.

Target Prediction of Hepatotoxic Metabolites Caused by Rhubarb

Based on the hepatotoxic metabolites obtained in the metabolomics analysis, the chemical structure was obtained using the PubChem database and saved as a file in SPF format, and the corresponding target was obtained in the STP database. We summarized the obtained targets, merged targets, and deleted duplicates. We then obtained the “metabolite-target” and used Cytoscape software to construct an interaction network diagram.

Network Construction Analysis of Liver Injury Targets

The obtained rhubarb component targets and metabolite targets were crossed to obtain the rhubarb component metabolite targets. The intersection targets and the CTD database (http://ctdbase.org/) “chemical and drug-induced liver damage” targets were intersected to obtain the targets of hepatotoxicity caused by rhubarb. The STRING (https://string-db.org/) database was used to analyze the protein interactions of the above targets, and Cytoscape software was used to obtain protein interaction network maps.

Analysis of Related Channels

For the intersecting target proteins, the David database (https://david.ncifcrf.gov/) was used for GO functional enrichment analysis and KEGG pathway analysis. GO function analysis was mainly used to describe the function of gene products, including biological process, cell composition, and molecular function. The KEGG pathway was enriched to analyze the possible action pathways of the target and used to explore the mechanism of hepatotoxicity caused by rhubarb.

Data Processing

Masslynx (US Waters) software was used to collect the original data of each group and perform ion pair extraction, peak alignment, peak matching, and peak intensity correction on the spectrum. Targetlynx (US Waters) software was used to update the retention time and peak intensity of the category, and SIMCA-P14.0 statistical software was utilized for multivariate statistical analysis.³⁶ Unsupervised principal component analysis (PCA) was performed. Within the 95% confidence interval, outliers of each group were eliminated. Supervised partial least squares-discriminant analysis (PLS-DA) was performed in UV mode, and cross-validation was performed to confirm that the established model was reliable. In the results of PLS-DA, substances with VIP values > 1 were selected as candidate markers.³⁷ After the t-test, the candidate markers were selected to have significant differences as preliminary toxic biomarkers.

Results

Identification of Rhubarb

UPLC-Q-TOF/MS was used to qualitatively analyze the chemical constituents extracted from rhubarb. The base peak intensity (BPI) of ESI-MS is shown in Figure 1. Using information, such as mass spectrometry fragments and compound spectra, a total of 28 compounds were identified from rhubarb extracts by comparing standards and literature as shown in Table 2. The specific comparison information and atlas of the standard and sample are included in the Supplementary Information. The substances included anthraquinones, anthracnose, tannins, flavonoids, and simple phenolic acids. The ion fragmentation pathway for the identification of rhubarb components is included in the supporting information.

Table 2 Identification Information for Rhubarb Total Extract Composition

	Ingredient Name	Molecular Formula	Rt. (Min)	m/z (Theoretical Value)	m/z (Actual Value)	Ppm	Fragmentation Information
1	Epicatechin-(4β-8)-epicatechin-(4β-8)-catechin	C45H38O18	1.82	865.1980	865.1974	−0.69	577.1228, 483.0695, 169.0170, 407.0745, 313.0255, 125.0257
2	(+)-catechin	C15H14O6	3.02	289.0712	289.0698	−4.84	245.0748, 579.0983, 203.0747, 151.0403, 136.0273
3	Procyanidin B1 3-O-gallate	C37H30O16	17.17	729.1456	729.1481	3.43	651.1469, 577.1382, 407.0789, 364.0630, 289.0731, 163.0426
4	Procyanidin 3,3′-di-O-gallate	C44H34O20	7.68	881.1565	881.1566	0.11	729.1423, 559.1328, 577.1371, 541.0665, 407.0714, 289.0776, 169.0195, 125.0238
5	Sennoside A	C42H38O20	16.93	861.1878	861.1892	1.63	699.1359, 431.0953, 269.0416, 283.0689
6	Chrysophanol	C15H10O4	31.89	253.0501	253.0502	0.40	225.0563, 209.0601, 197.0291
7	Aloe-emodin	C15H10O5	28.90	269.0450	269.0452	0.74	240.0415, 211.0441, 183.0452
8	Rhein	C15H8O6	29.64	283.0243	283.0258	5.30	239.0348, 211.0394, 183.0448
9	Emodin	C15H10O5	30.91	269.0450	269.0452	0.74	241.0411, 225.0558, 197.0604, 183.0451
10	Gallic acid	C7H6O5	1.04	169.0137	169.0143	3.55	125.0246
11	Sennoside B	C42H38O20	12.48	861.1878	861.1873	−0.58	699.1273, 386.1023
12	Emodin-O-malonyl-glucoside	C24H22O13	26.71	517.0982	517.0979	−0.58	473.1086, 253.0508
13	Physcion-O-acetyl-glucoside	C24H24O11	9.41	487.1240	487.1239	−0.21	283.0240, 240.0522
14	Physcion	C16H12O5	29.10	283.0606	283.0614	2.83	269.0434, 240.0432, 225.0568
15	P-Hydroxybenzoic acid-O-galloyl-glucoside	C20H20O12	2.32	451.0877	451.0873	−0.89	451.0930, 289.0698, 281.0860, 151.0402, 143.0316
16	Proanthocyanidin B1	C30H26O12	1.49	577.1346	577.1353	1.21	483.0776, 289.0713, 407.0789, 169.0146
17	Proanthocyanidin B2	C30H26O12	1.74	577.1346	577.1339	−1.21	577.1334, 559.1176, 493.0526, 451.0984, 425.0836, 407.0782, 289.0691, 245.0802, 125.0254
18	Procyanidin B2 3ʹ-O-gallate	C37H30O16	3.22	729.1456	729.1443	−1.78	577.1407, 289.0725, 169.0140, 125.0233
19	Procyanidin B-5,3ʹ-O-gallate	C37H30O16	4.45	729.1456	729.1449	−0.96	577.1301, 407.0771, 289.0724, 169.0121, 125.0264
20	Aloe-emodin-8-O-β-D-glucopyranoside	C21H20O10	9.20	431.0978	431.0975	−0.70	269.0456, 240.0417, 225.0561, 195.0813, 169.0121
21	Emodin-1-O-β-D-glucoside	C21H20O10	18.14	431.0978	431.0977	−0.23	269.0465, 240.0434, 169.0159, 125.0226
22	Emodin-8-O-β-D-glucoside	C21H20O10	24.81	431.0978	431.0988	2.32	269.0451, 225.0549, 169.0119, 253.0509
23	Rhein-8-O-β-D-glucopyranose	C21H18O11	9.95	445.0771	445.0773	0.45	283.0237, 239.0356, 211.0396, 891.1628
24	Chrysophan-1-O-β-D-glucopyranose	C21H20O9	25.21	415.1029	415.1027	−0.48	253.0507, 225.0549
25	Gallic acid-3-O-(6ʹ-O-galloyl)-glucoside	C20H20O14	1.51	483.0775	483.0782	1.45	483.0780, 465.0664, 439.0826, 331.0663, 321.0301, 313.0497, 295.0297, 287.0656, 169.0137, 151.0013, 125.0237
26	Catechin-7-glucoside	C21H24O11	1.44	451.1240	451.1245	1.11	451.1241, 289.0717, 245.0817, 227.0677, 203.0705, 161.0662
27	Epicatechin-3-glucoside	C21H24O11	1.60	451.1240	451.1251	2.44	451.1229, 289.0698, 179.0323, 137.0249
28	(-)-Epicatechin-3-O-Gallate	C22H18O10	8.26	441.0822	441.0836	3.17	441.0860, 331.0469, 289.0740, 245.0804, 205.0542, 203.0715, 193.0124, 179.0388, 169.0152, 125.0251

Figure 1 Rhubarb total extract with negative ion mode BPI.

Rat Body Weight, Biochemical, and Pathological Results

The weight of each group of rats was weighed every day during the administration period. The change in body weight is shown in Figure 2A. The four curves in the figure represent the changes in the body weight of the normal control group (C), low-dose group (L), medium-dose group (M), and high-dose group (H) during the administration period. Throughout the entire administration period, the weight of the rats in each group showed an increasing trend. However, the weight of the rats in each dose group changed after the rhubarb total extract was administered. With the increase of the dose of rhubarb, the weight of rats increased more slowly. After 20 days of administration, the rats in the administration group weighed less than the control group. These findings indicated that the total extract of rhubarb has a certain degree of toxicity in rats.

Figure 2 Trend chart of body weight change of each group (A). Changes in the contents of ALT (B) and AST (C) in the serum samples of different groups. Liver pathological results (D) for the C group (a1-a2), L group (b1-b2), M group (c1-c2) and H group (d1-d2) with 200x magnification. n= 6, 7, 7, and 6, respectively.

Notes: (A) *P value of the treatment group compared with the control group is less than 0.05, and ^#P value of the treatment group compared with the control group is less than 0.01; (B and C) *P value <0.05 and **P value <0.01.
Abbreviations: C group, normal control group; L group, low-dose group; M group, medium-dose group; H group, high-dose group.

To evaluate the degree of damage of rhubarb to the liver, serum ALT and AST levels were used as biochemical indicators of liver damage. The results of the biochemical indicators detected by the experimental groups are shown in Figure 2B and C. The AST and ALT levels increased with increasing administered dose. During the 20-day administration process, the AST and ALT levels increased with increasing doses. The level of ALT increased significantly in the low and medium dose groups. The level of AST in the high-dose group increased significantly. These findings demonstrated that there was a certain degree of hepatotoxicity in rats after administration of the total extract of rhubarb.

The results of liver pathological sections are shown in Figure 2D. In the normal control group (a1-a2), the liver cells were regularly arranged with clear boundaries. The cytoplasm was loosely stained, and there was no obvious abnormality in the interstitium and no obvious inflammation. The cells in the low-dose group (b1-b2) were regularly arranged and clearly demarcated, but there was a small number of fibrotic cells as shown by black arrows. The fibrotic cells were accompanied by a small amount of inflammatory cell exudation as shown by red arrows. In the medium-dose group (c1-c2), a small amount of tissue damage was observed in the field of vision. A small number of cells disappeared, and a small amount of connective tissue hyperplasia and fibroblasts was observed as shown by black arrows. Additionally, a small amount of inflammatory cell exudation was observed with occasional granulocyte presence and lymph cell infiltration as shown by the red arrow. Liver tissue in the high-dose group (d1-d2) was extensively damaged. The connective tissue was relatively loose with a small amount of cell fibrosis as shown by the black arrow, and it was accompanied by inflammatory cell infiltration. In a large number of hepatocytes, the cytoplasm was loosely stained as shown by the red arrow, and a small amount of vacuoles, liquefied necrosis and steatosis was observed as shown by the yellow arrow. A small amount of nucleocytoplasm was observed in and around the venous cavity, and the nucleolus was obvious.

Metabolomics Results

PCA and PLS-DA Analysis

Using UPLC-Q-TOF/MS technology analysis, the four groups of rat plasma samples were subjected to metabolomics analysis. The plasma QC samples were analyzed by mass spectrometry in positive ion mode, and the resulting BPI diagram is shown in Figure 3.

Figure 3 Plasma QC sample with positive ion mode BPI.

PCA determined that the samples were within the 95% confidence interval (Figure 4A). The PCA results (PCA: R2X = 0.299, Q2 = 0.135) in the two-dimensional scatter plot were unable to distinguish three different doses. The samples in the dose groups and normal group had poor aggregation. There was a certain difference between the normal group and the different dose groups, which showed that the Chinese medicine rhubarb caused changes in endogenous metabolites in the body. The four groups of samples had differences in metabolic levels.

Figure 4 Two-dimensional PCA graphs of four groups (A). PLS-DA graphs of normal control group and toxic group (B). PLS-DA graphs of normal control group and low-dose group (C). PLS-DA graphs of normal control group and low-dose group (D). PLS-DA graphs of normal control group and low-dose group (E).

On this basis, PLS-DA analysis was used to mine potential data information and further screen for differential metabolites with significant changes. Four groups of PLS-DA models were established (PLS-DA: R2X = 0.206, R2Y = 0.488, and Q2 = 0.0741). The results showed that each group was separated into different areas (Figure 4B). PLS-DA analysis was then used to build models for the following groups: normal and low-dose groups; normal and medium-dose groups; and normal and high-dose groups (PLS-DA-L: R2X = 0.318, R2Y = 0.963, and Q2 = 0.55; PLS-DA-M: R2X = 0.226, R2Y = 0.983, and Q2 = −0.21; PLS-DA-H: R2X = 0.325, R2Y = 0.957, and Q2 = 0.414). Thus, the VIP > 1 substance was selected as a candidate liver toxicity marker. The specific results are shown in Figure 4C–E.

Determination and Quantification of Markers

On the basis of screening candidate markers with VIP > 1, a t-test group difference test was performed to identify substances with VIP > 1 and P value < 0.05 as liver toxicity biomarkers. The biomarkers were tested by standard curve quantitative analysis. Through a database search of differential metabolites and further identification of related literature, nine hepatotoxicity markers were finally determined. The specific information of the selected biomarkers is shown in Table 3.

Table 3 Ionic Information for Potential Biomarkers for Liver Toxicity

No.	Metabolites	Molecular Formula	Rt. (Min)	m/z (Measured Value)	m/z (Theoretical Value)	Error (ppm)	Adduct Form	Trend	VIP Value
1	Dopamine	C8H11NO2	2.49	176.0704	176.0687	9.66	M+Na	↓*	2.26
2	Biopterin	C9H11N5O3	5.77	276.0524	276.0499	9.06	M+K	↓*	1.89
3	Choline	C5H14NO	0.59	104.1077	104.1075	1.92	M+H	↓**	2.37
4	Dynorphin B (10–13)	C20H39N5O6	2.65	446.2906	446.2979	−16.36	M+H	↑*	1.94
5	Coenzyme Q9	C54H82O4	10.92	833.5842	833.5850	−0.96	M+K	↓*	1.97
6	Cervonoyl ethanolamide	C24H36O3	2.88	373.2741	373.2743	−0.54	M+H	↑*	2.19
7	P1,P4-Bis(5ʹ-uridyl) tetraphosphate	C18H26N4O23P4	10.93	812.9782	812.9836	−6.64	M+Na	↓*	1.91
8	LysoPE(18:2)	C23H44NO7P	4.91	478.2941	478.2934	1.46	M+H	↑**	2.36
9	3-hydroxyphenyl 2-hydroxybenzoate	C13H10O4	3.13	253.0497	253.0477	7.90	M+Na	↑**	3.02

Notes: *P value <0.05; and **P value <0.01. The upward arrow indicates that the level of metabolites in the toxic group is up-regulated compared to the control group. The downward arrow indicates that the level of metabolites in the toxicity group was lower than that in the control group.

Marker Analysis

The heatmap and ROC curve were used to analyze and optimize the markers. Table 3 shows that a total of 9 hepatotoxicity biomarkers were screened out. The dopamine, biopterin, choline, coenzyme Q9, and P1, P4-Bis (5ʹ-uridyl) tetraphosphate metabolic markers decreased significantly compared to the normal control group, and the dynorphin B (10–13), cervonoyl ethanolamide, lysoPE (18: 2), and 3-hydroxyphenyl 2-hydroxybenzoate metabolic markers increased significantly compared to the normal control group.

The hierarchical cluster analysis method (heatmap) was used to analyze the markers. Figure 5 shows the content differences of nine the metabolic markers in the hepatotoxicity group. Taking the variable information as the ordinate and the sample information as the abscissa, the depth of the color reflects the size of the variable value. The closer the bifurcation of the variable information in the vertical axis is, the higher the similarity of the substances; that is, they may be derived from metabolites of the same substance. In the present study, the heatmap colors indicated a significant difference between the normal control group and toxic group.

Figure 5 Heatmap of metabonomic data depicting the data structure of nine biomarkers.

Notes: The columns represent differentiating metabolites, and the rows indicate the samples. C1-C6, LI-L7, M1-M7 and H1-H7 belong to the normal control group, low-dose group, medium-dose group and high-dose group, respectively. The depth of the color represents the size of the variable.

The ROC curve was used to optimize the nine markers in Figure 6. The AUCs of three biomarkers were greater than 0.9, indicating that they had a good diagnostic ability for rhubarb-induced liver toxicity. The AUCs of six biomarkers were greater than 0.7, indicating good diagnostic capabilities as well. Information on the area under the curve of the nine metabolites is provided in the Supplementary Information.

Figure 6 ROC curve for evaluating the predictive ability of nine rhubarb hepatotoxicity biomarkers.

Network Pharmacology Analysis

Target Prediction of Rhubarb Composition

In this study, 28 chemical components of rhubarb were obtained by measurement. The structural formulas of these components were obtained through PubChem, and 242 targets of components were obtained by importing them into the STP database. The component-target interaction network of rhubarb constructed by Cytoscape software is shown in Figure 7.

Figure 7 Rhubarb “ingredient-target” network diagram.

Note: Orange represents the rhubarb component, and blue represents the target.

Target Prediction of Hepatotoxic Metabolites Caused by Rhubarb

Based on the UPLC-Q-TOF/MS metabolomics technology in the rat liver toxicity studies, nine hepatotoxicity metabolic biomarkers were screened. In total, 282 targets corresponding to nine metabolites were obtained from the STP database, and the “metabolite-target” interaction network was constructed as shown in Figure 8.

Figure 8 “Metabolite-target” network diagram.

Notes: Diamonds represent metabolic markers, and circles represent targets. The color from pink to blue represents the correlation from low to high.

Rhubarb Potential Hepatotoxicity Target Protein Network Analysis

Intersecting 242 rhubarb component targets with 282 metabolite targets generated 184 copredicted targets. After intersecting with the “chemical and drug-induced liver damage” targets in the CTD database, 166 common targets were obtained. Using these 166 targets as potential targets for hepatotoxicity of rhubarb, the target-protein interaction network map of 166 potential toxic targets of rhubarb was drawn using the String database (protein–protein interaction, PPI) as shown in Figure 9. The required interaction score was set to 0.400. The PPI network contained 166 nodes, 4247 edges, and an average node degree of 50.9. Each node in the figure represents a protein, and the connections between the nodes represent the interaction between the proteins. The common target may be related to the mechanism of hepatotoxicity caused by rhubarb. The network diagram of protein interaction relationships helps to more intuitively understand the functional relationships between proteins.

Figure 9 Rhubarb potential hepatotoxicity target PPI network.

Note: The PPI network contained 166 nodes, 4147 edges, an average node degree of 50.9 and enrichment p value <0.05.

GO Gene Enrichment Analysis and KEGG Pathway Analysis

To further discover the mechanism of hepatotoxicity caused by rhubarb, the 166 targets obtained above were analyzed for GO gene enrichment and KEGG pathway analysis using the DAVID database. The GO gene enrichment analysis of the targets was mainly divided into the following three parts: biological process (BP), cellular component (CC), and molecular function (MF). The GO analysis results showed that the biological process (BP) obtained a total of 149 enrichment results, mainly including the adenylate cyclase-inhibiting G-protein coupled receptor signaling pathway, adenylate cyclase-activating adrenergic receptor signaling pathway, phospholipase C-activating G-protein coupled receptor signaling pathway, sensory perception of pain, and locomotory behavior. Cell composition (CC) obtained 24 enrichment results, mainly including integral components of the plasma membrane, postsynaptic membrane, neuron projections, dendrites, and cell junctions. Molecular function (MF) obtained 34 enrichment results, mainly including serotonin binding, drug binding, G-protein coupled acetylcholine receptor activity, neuropeptide binding, and peptide binding. The enriched results were ranked by P value, and the R programming language was used to plot the analysis results of the top 10 rankings of the three parts as shown in Figure 10.

Figure 10 GO gene enrichment function analysis diagram.

Notes: The Y axis represents the name of the gene enrichment result, and the X axis represents the number of genes. The graphic color represents the significance of the enrichment, that is, the magnitude of the P value.

Through the enrichment analysis of target pathways related to KEGG gene annotation function analysis, 46 pathways were obtained. Among them, there were 46 pathways with P value <0.05. These pathways mainly included neuroactive ligand-receptor interactions, calcium signaling pathways, serotonergic synapses, cAMP signaling pathways, cocaine addiction, and cGMP-PKG signaling pathways. The results of the visual enrichment analysis, which were the first 20 channels of the P value in the KEGG analysis results, are shown in Figure 11.

Figure 11 KEGG pathway enrichment analysis diagram.

Notes: The Y axis represents the name of the pathway, and the X axis represents RichFactor, that is, the number of genes belonging to the pathway in the target gene/the number of all genes in the pathway in the background gene set. The size of the bubble area represents the pathway belonging to the target gene set and the number of genes. The bubble color represents the significance of enrichment, that is, the magnitude of the P value.

Discussion

Clinical Liver Injury of Rhubarb

In addition to being used in clinical practice for constipation, rhubarb is also widely used in the clinical treatment of chronic hepatitis, cholecystitis, diabetic nephropathy, renal failure, and other diseases.^38,39 However, in recent years, it has been reported that certain components in rhubarb can cause liver damage to the body, and long-term or unreasonable clinical application may cause drug-induced liver damage.⁴⁰ The literature has conducted statistics on cases of TCM causing DILI in hospitals in recent years.^41–44 In these studies, liver injury was determined based on the causal relationship between drugs and liver injury assessed by the RUCAM scale, which is recommended in the Guidelines for the Diagnosis and Treatment of Drug-induced Liver Injury. The proportion of DILI caused by taking Chinese herbal medicine and proprietary Chinese medicine was the highest in the cases, and the Chinese medicinal preparation containing Radix Polygoni Multiflori, Radix Notoginseng, and Radix et Rhizoma Rhei was the most likely to cause drug-induced liver injury.⁴¹ The medication time and dose should be strictly controlled in clinical use, and liver function should be closely monitored during medication.

Potential Toxic Ingredients in Rhubarb

In the present study, data mining was performed by summarizing the mass spectrum information of rhubarb components mentioned in the literature and combining with UPLC-Q-TOF/MS technology and a data postprocessing strategy, which identified 28 rhubarb components by comparing the ion fragment information of the components.^45–49 Aloe-emodin, emodin, rhein, chrysophanol and physcion were found to be the main anthraquinones in rhubarb. It has been reported that anthraquinone is an effective component of rhubarb and may also be its toxic component.⁵⁰ Gallic acid-3-O-(6ʹ-O-galloyl)-glucoside and procyanidin B-5,3ʹ-O-gallate are tannin components in rhubarb and are also components that cause hepatotoxicity in Polygonum multiflorum.^51,52 It has been demonstrated that gallic acid is a substance that causes hepatocyte damage and may also be a potentially toxic ingredient triggering hepatotoxicity in rhubarb.⁵³ However, the hepatotoxicity of other components in rhubarb is not clear. Therefore, the intersection targets of component targets, metabolite targets and liver injury targets in the database were used as potential targets of liver toxicity to reverse explore the components.⁵⁴ The results showed that some targets of the 28 components were in the liver toxicity intersection targets.

Target Prediction of Metabolite Markers Related to Rhubarb-Induced Hepatotoxicity

The rhubarb-induced hepatotoxicity in rats can be markedly observed through biochemical and pathological results. In addition, nine biomarkers related to hepatotoxicity caused by rhubarb were obtained through metabolomics. Among them, dopamine, biopterin, choline, coenzyme Q9, cervonoyl ethanolamide and lysoPE (18:2) have been demonstrated to be related to liver injury in the literature.^55–60 A previous study has shown that the incidence of liver tissue lesions and plasma aminotransferase in rats after dopamine treatment are reduced and that the survival rate is improved, suggesting that the increase in dopamine may reduce fulminant liver injury in mice.⁶¹ At the same time, experiments have shown that in liver injury caused by liver ischemia, the excitatory neurotransmission caused by norepinephrine and dopamine is suppressed.⁶² These results are consistent with the results of the decrease in dopamine content after liver injury in this experiment. In the present study, the choline content decreased significantly. A previous study has used methionine-choline deficiency and a high-fat diet to induce fatty liver injury, and some animal models of chronic liver injury have been generated using methionine-choline deficiency in the experiments.^63,64 In a previous study, the damage to liver cells mediated by oxidative stress was prevented through the addition of coenzyme Q9, which has antioxidant activity.²⁷ These findings also indirectly suggested that the decrease in coenzyme Q9 content in the results is correct. Biomarkers are endogenous changes in the body when affected by external factors, and they can be used as abnormal signal indicators at different biological levels.⁶⁵ Therefore, the rhubarb-induced liver toxicity of the body was reflected by changes in metabolic markers, and the targets related to metabolic markers were predicted via network pharmacology.

Potential Targets of Rhubarb Hepatotoxicity and Protein–Protein Interactions

The rhubarb component targets, hepatotoxicity metabolic marker targets, and liver injury-related targets in the CTD database were crossed to obtain the potential targets of rhubarb hepatotoxicity that were used for protein interaction analysis. Protein–protein interactions (PPIs) are the basis of intracellular biological processes.⁶⁶ The protein–protein interaction (PPI) network was constructed using the STRING database for potential hepatotoxic targets of rhubarb in the results. The interaction network among proteins is useful for deciphering the molecular mechanism of biological functions. In addition, the protein interaction diagram also contains systematic information and related signaling pathways.⁶⁷ Through node data analysis, we identified that DRD1, DRD2, PDE4B, TRPV1, TRPM8, TRPA1, PTGER2, and CHRM3 are important nodes. GO function analysis and KEGG pathway enrichment analysis were performed on the potential targets of rhubarb-induced hepatotoxicity after protein interaction, and pathways enriched by potential targets of hepatotoxicity were analyzed.

Rhubarb-Induced Hepatotoxicity Analysis of Common Potential Targets in the KEGG Prediction Pathway

cAMP Signaling Pathway and Liver Injury

cAMP is a key second messenger molecule that regulates various cell functions, including inflammation and cell damage, by affecting gene/protein expression and function. Dopamine receptors (DRs) belong to the G-protein coupled receptor family (GPCR). Both dopamine receptor D1 (DRD1) and dopamine receptor D2 (DRD2) are expressed in hepatocytes, and their receptor agonists can increase cyclic adenosine phosphate (cAMP).⁶⁸ Increased cAMP inhibits apoptosis induced by the JNK or PI3K/Akt pathway through the key conversion factor of guanine nucleotides to protect liver cells.⁶⁹ Stimulating DRD1 negatively regulates the activation of NLRP3 inflammatory bodies through cAMP, thereby regulating the occurrence of inflammation.⁷⁰ In addition, DR exists in a variety of immune cells, such as regulatory T cells (Tregs), effector T cells (Teffs), and NK cells. DRD1 involves the inhibition of CD4⁺ and CD8⁺ T cells. Treg cells and Teff cells produce dopamine, possibly through DRD2, to inhibit the function of their cells so that Teff cells are activated and proliferate.⁷¹ Cyclic nucleotide phosphodiesterase (PDE) regulates cAMP and cGMP levels as well as signals through its hydrolysis, thereby affecting cAMP and cGMP-dependent processes.⁷² The PDE4 family (PDE4B) is widely expressed in cAMP-specific PDEs, accounting for the majority of intracellular cAMP hydrolytic activity, and it is used as a target for various inflammatory diseases. The ingredients of rhubarb may inhibit DRs and PDE to promote the occurrence of inflammation and aggravate liver damage, resulting in liver toxicity.

Cholinergic Synaptic and Liver Injury

Acetylcholine (ACh) is a neurotransmitter widely distributed in the central nervous system (CNS).⁷³ When ACh is released in a synaptic cleft, it can bind to two different types of receptors as follows: the ionic nicotinic acetylcholine receptor (nAChR) and the metabolic muscarinic acetylcholine receptor (mAChR).⁷⁴ The activation of mAChRs by ACh is relatively slow, but through the subtypes of mAChRs (M1-M5), it can directly change the cellular homeostasis of phospholipase C, inositol triphosphate, cAMP, and free calcium. The key cholinergic receptor, muscarinic subtype 3 (Chrm3), is a specific gene that induces liver injury. The Chrm3 muscarinic receptor inhibits hepatocyte proliferation and increases hepatic nodules, apoptosis, and fibrosis. At the same time, Chrm3 regulates hepatocyte proliferation and apoptosis in the liver injury response.⁷⁵

TRP Channel Inflammatory Mediators and Liver Injury

TRPV1, TRPM8, and TRPA1 in the TRP channel can be expressed in a subset of nociceptive ganglion neuronal cell bodies (DRGs). These channels can be indirectly regulated by inflammatory mediators (such as PGE2, NGF, and proinflammatory cytokines) produced during liver tissue damage. Prostaglandin E receptor 2 (EP2) and prostaglandin E receptor 4 (EP4) are immunomodulatory receptors with prostaglandin E2 (PGE2), which is a ligand. PGE2 inhibits the expression of cirrhotic macrophages through EP2. After liver injury, the level of PGE2 increases, and EP2 exerts an immunosuppressive effect on cirrhotic macrophages.⁷⁶ Under cell death conditions, the induction of the cyclooxygenase 2 gene (COX2) increases the secretion of PGE2, which is more closely related to liver failure and the severity of liver disease.⁷⁷ EP2 also plays an important role in the early stage of liver failure. As EP2 further blocks PBMCs, liver failure is negatively correlated with EP2 expression. In addition, the expression of EP2 is negatively correlated with MIP-1β. The blockade of EP2 increases the secretion of MIP-1β, indicating that the downregulation of EP2 leads to an increase in the liver failure index, MIP-1β.⁷⁸ It has also been reported that during acute liver injury or liver failure, hepatocytes die and release DAMPs in the pathway to trigger further injury.^79,80

Conclusion

In the present study, we found that rats treated with rhubarb extract showed a certain degree of liver damage. The metabolomics results indicated that nine metabolic biomarkers were related to hepatotoxicity caused by rhubarb. On the basis of the network pharmacology results, it can be inferred that rhubarb may cause liver toxicity due to its action on DRD1, DRD2, PDE4B, TRPV1, TRPM8, TRPA1, PTGER2, ACHE, and CHRM3 through the cAMP signaling pathway, cholinergic synapses, and inflammatory mediators that regulate TRP channels. This study provides a reference for research on rhubarb-induced hepatotoxicity and provides a basis for the safe application of rhubarb. However, the mechanism between the target and the metabolite pathway needs to be further studied.

Abbreviations

UPLC-Q-TOF/MS, ultrahigh-performance liquid chromatography-tandem quadrupole/time-of-flight mass spectrometry; DRD1, dopamine D1 receptor; DRD2, dopamine D2 receptor; PDE4B, phosphodiesterase 4B; TRPV1, vanilloid receptor; TRPM8, transient receptor potential cation channel subfamily M member 8; TRPA1, transient receptor potential cation channel subfamily A member 1; PTGER2, prostanoid EP2 receptor; ACHE, acetylcholinesterase; CHRM3, muscarinic acetylcholine receptor M3; TCM, traditional Chinese medicine; ALT, alanine aminotransferase; AST, aspartate aminotransferase; GO, gene ontology; PCA, principal component analysis; KEGG, Kyoto Encyclopedia of Genes and Genomes; PLS-DA, partial least squares discriminant analysis; PCA, principal component analysis; RSD, relative standard deviation; VIP, importance in projection; QC, quality control; ROC, receiver operating characteristic; AUC, area under curve.

Disclosure

The authors declare no conflicts of interest in this work.

Acknowledgment

This work was supported by the National Natural Science Foundation of China No. 81873194, Tianjin Talent Development Special Support Project for High Level Innovation and Entrepreneurship, and the National Natural Science Foundation of China NO. 81903933.