New bitongling, a traditional Chinese medicine (TCM) compound, is formulated using herbs such as Cinnamomi ramulus (Gui Zhi), Saposhnikoviae radix (Fang Feng), Ephedrae herba (Ma Huang), Sinomenii caulis (Qing Feng Teng), Aconiti radix (Chuan Wu), and Vespae nidus (Feng Fang), sourced from the pharmacy at Jiangsu Provincial Hospital of TCM. Prior to preparation, all raw herbs underwent stringent quality control measures, including morphological identification and visual inspection to verify the authenticity and purity of each herb. Additionally, each herb was standardized based on established pharmacopeial standards (China Pharmaceutical Association, 2020), ensuring that they met predefined criteria for key active components to guarantee consistency across different batches.

The preparation process involves several key steps: Aconiti radix (Chuan Wu) is first decocted in double-distilled water at 100°C for 30 min. Subsequently, the remaining herbs are added and the mixture is decocted for an additional 45–60 min. The resulting liquid is then filtered to obtain the herbal extract. This decoction process is repeated to ensure thorough extraction of active constituents. The combined extracts are concentrated and centrifuged at 5,000 rpm for 30 min to eliminate residues, yielding a final concentration of 1.22 g/mL. The decoction is then aliquoted and stored at −20°C for future use. All preparation batches underwent quality control checks, including visual inspection and adherence to standardized preparation protocols, to confirm the consistency and potency of the extracts.

For the positive control, Tripterygium wilfordii tablets containing 10 mg of the active ingredient were used. These tablets, produced by Zhejiang Deyongde Pharmaceutical Co., Ltd. (Batch number: 2204111B), were ground into a fine powder and dissolved in distilled water to create a suspension. The suspension was standardized to a final concentration of 8 mg/mL through precise measurement and preparation techniques to ensure the stability and uniformity of the active components for experimental use.

In this study, we selected 6-week-old female Sprague–Dawley (SD) rats, weighing between 180–200 g, sourced from Zhejiang Vital River Laboratory Animal Technology Co., Ltd. (License No. SCXK(Zhe)2021-0006). The rats were acclimatized for 1 week before the experiment in a specific pathogen-free (SPF) facility at the Affiliated Hospital of Nanjing University of Chinese Medicine (License No. SYXK(Su)2017-0069), maintained under controlled conditions (22 ± 2°C, 45–60% humidity). They had free access to SPF-grade feed (Synergy Bio, 1010086) and water. The study adhered to National Institutes of Health (NIH) guidelines for the care and use of laboratory animals and received approval from the Ethics Committee of the Affiliated Hospital of Nanjing University of Chinese Medicine (Approval No. 2022 DW-37-01).

For arthritis induction, bovine type II collagen (CII) (Beijing Bolide Biotech, 20021) was dissolved in 0.1 mol/L glacial acetic acid to a concentration of 4 g/L, stirred on ice, and stored overnight at 4°C. This solution was then emulsified with an equal volume of complete Freund’s adjuvant (CFA) (Sigma, F5881) to create a 2 g/L collagen II emulsion. The rats were randomly assigned to either a normal group or a collagen-induced arthritis (CIA) model group. The CIA induction protocol involved intradermal injections of 0.3 mL of the CII emulsion on day 1, followed by a booster dose on day 7. The primary immunization consisted of a 1:1 mix of bovine type II collagen and CFA, administered at a 2 g/L concentration. The booster immunization used a similar 1:1 mix with incomplete Freund’s adjuvant (IFA) (Sigma, F5506). Arthritis severity was evaluated on day 14 using a visual scoring system for each paw (0–4, with a maximum total score of 16): 0 for normal joints, 1 for slight swelling or erythema, 2 for moderate erythema and swelling, 3 for pronounced erythema and swelling, and 4 for severe swelling, erythema, or ankylosis (Trentham et al., 1977).

The rats with successfully induced arthritis were randomly allocated into different treatment groups, along with a blank control group selected from the normal group. The groups received the following treatments via oral gavage for 28 days (Li et al., 2021): blank control group received distilled water; model group received distilled water; low dose new bitongling group received 2.75 g/kg/day of new bitongling; medium dose new bitongling group received 5.5 g/kg/day of new bitongling; high dose new bitongling group received 11 g/kg/day of new bitongling; and the positive control group received 0.008 g/kg/day of Tripterygium wilfordii solution. The selected doses for the low, medium, and high new bitongling (NBTL) groups were based on dose extrapolation from equivalent human doses using the standard body surface area (BSA) conversion formula (Nair and Jacob, 2016; Reagan-Shaw et al., 2008). These doses were further refined through preliminary experiments reported in our prior study, which demonstrated the efficacy and safety of these specific concentrations in a similar model (Li et al., 2021). The rats were monitored daily for changes in appearance, behavior, food intake, and excretion, with body weight measurements taken weekly. After 28 days, the rats were anesthetized for blood serum collection and then euthanized for the collection of ankle and knee synovial tissue.

We selected 6–8 week-old SPF-grade SD rats and randomly divided them into the traditional Chinese medicine (TCM) group and the blank control group. The TCM group received an oral gavage of new bitongling (NBTL) at the optimal dose for in vivo rat experiments (high dose: 11 g/kg/day, with a drug concentration of 1.22 g/mL). The blank control group was administered an equal volume of physiological saline. Both treatments were administered daily for seven consecutive days. One hour after the final administration, the rats were anesthetized, and blood samples were collected. The blood was allowed to clot at 4°C for 2 h to facilitate serum and clot separation. The samples were then centrifuged at 3,000 rpm for 10 min to separate the serum from the clot. The supernatant, which is the serum, was carefully collected. The collected serum was then inactivated at 56°C to eliminate any potential biological activity. This step is performed after serum separation. Following inactivation, the serum was filtered for sterilization, aliquoted, and stored at −80°C for subsequent experiments (Pan et al., 2023; Pei et al., 2024).

We used ultra-performance high-performance liquid chromatography coupled with time-of-flight mass spectrometry (UPHLC-MS-TOF) to identify the prototypical components and metabolites of new bitongling (NBTL) in rat serum. Serum samples were prepared by protein precipitation and filtration. Chromatographic separation was performed on an acquity UPLC BEH C18 column, followed by mass spectrometric analysis with an electrospray ionization source in both positive and negative ion modes. Data processing included peak detection, alignment, and compound identification through database matching.

Ion peak data were filtered based on various conditions, with missing values removed. Positive and negative ion data were processed separately and converted to a long format for plotting, retaining only NBTL medicated serum group data. Vertical offsets were applied to avoid overlapping annotations. The “ggplot2” package (Villanueva and Chen, 2019) was used to visualize the data, creating plots for positive and negative ions to illustrate the prototypical components and metabolites.

Synovial tissue samples were fixed in 4% paraformaldehyde (PFA) (Aladdin, C104188; Sinopharm, 30525-89-4), dehydrated through graded ethanol solutions (Sinopharm, 100092683), and embedded in paraffin. Sections were cut at 4 microns, dewaxed, and rehydrated. Hematoxylin staining (Zhuhai Beisuo, BA4097) was used to highlight cell nuclei, followed by eosin staining (Zhuhai Beisuo, BA4099) to highlight cytoplasm. The sections were then dehydrated with anhydrous ethanol (Sinopharm, 100092683), cleared with xylene (Sinopharm, 1330-20-7), and mounted with neutral balsam (Sinopharm, 10004160). Hematoxylin and eosin (HE) staining was evaluated for synovial cell arrangement, hyperplasia, and inflammatory cell infiltration, with scores ranging from 0 (no lesion) to 4 (very severe) based on the extent of disorganization, hyperplasia, and infiltration. The average score was calculated to assess tissue damage severity across samples (Guo et al., 2017).

Serum samples stored at −80°C were thawed to room temperature for enzyme-linked immunosorbent assay (ELISA) analysis. The following reagents and consumables were used: rat interleukin-6 (IL-6) ELISA kit (CSB-E04640r), rat interleukin-17 (IL-17) ELISA kit (CSB-E07451r), and rat tumor necrosis factor alpha (TNF-α) ELISA kit (CSB-E11987r) from Cusabio. Biotin-labeled antibodies and HRP-labeled streptavidin were diluted with their respective diluents at a ratio of 10 μL to 990 μL and kept on ice. The wash buffer was prepared by mixing 1 mL of concentrate with 24 mL of deionized water. Standard proteins were prepared through serial dilutions. The ELISA procedure involved adding standard proteins and serum samples to a coated 96-well plate, followed by incubation at 37°C. After 2 h, 100 μL of biotin-labeled antibody was added and incubated for another hour. The plate was then washed five times with 200 μL of wash buffer per wash, followed by the addition of 100 μL of HRP-labeled streptavidin and a 1-h incubation. The wash step was repeated, and 90 μL of 3,3′,5,5′-tetramethylbenzidine (TMB) substrate was added and incubated in the dark. Finally, 50 μL of stop solution was added, and the absorbance was measured at 450 nm using a microplate reader.

In this study, western blot analysis was performed on synovial tissue samples. Key reagents included 10X Transfer Buffer, 10X Electrophoresis Buffer (pH 8.3), and 10X Tris Buffered Saline with Tween-20 (TBST) Buffer. We prepared working solutions such as 10% sodium dodecyl sulfate (SDS), 10% ammonium persulfate (APS), and Tris-HCl buffers at pH 8.8 and 6.8, along with stacking and separating gels of varying concentrations for sodium dodecyl sulfate-polyacrylamide gel electrophoresis (SDS-PAGE). The synovial tissues were minced, ground in liquid nitrogen, and lysed using radioimmunoprecipitation assay (RIPA) buffer with a protease inhibitor cocktail (Sigma, P8340). Protein concentrations were determined using the bicinchoninic acid (BCA) protein assay kit (Beyotime, P0010). The samples were resolved on SDS-PAGE gels and transferred to polyvinylidene difluoride (PVDF) membranes (Millipore, ISEQ00010). Membranes were blocked with 5% skim milk (BD, 232100) in TBST, then incubated overnight with primary antibodies in 1% bovine serum albumin (BSA)/TBST (Sigma-Aldrich, V900933). The primary antibodies used were glyceraldehyde 3-phosphate dehydrogenase (GAPDH) (Zen-Bio, 200306-7E4), vascular endothelial growth factor (VEGF) (Santa Cruz, sc-53463), vascular endothelial growth factor receptor 1 (VEGFR1) (Proteintech, 19003-1-AP), vascular endothelial growth factor receptor 2 (VEGFR2) (Abcam, ab32152), hypoxia-inducible factor 1-alpha (HIF-1α) (Santa Cruz, SC-13515), and apoptosis-associated speck-like protein containing a CARD (ASC) (Santa Cruz, sc-514414). Secondary antibodies included HRP-conjugated anti-mouse (Beyotime, A0216) and anti-rabbit (Beyotime, A0208), as well as fluorescent antibodies Alexa Fluor 488 Donkey Anti-Rabbit IgG (H + L) (Abcam, ab150061) and Alexa Fluor 488 Donkey Anti-Mouse IgG (H + L) (Abcam, ab150105). Detection was performed using an enhanced chemiluminescence (ECL) substrate (Tanon, 180-5001) and visualized with a Tanon5200 Imaging System.

For quantitative polymerase chain reaction (qPCR) detection of mRNA and miRNA in rats and cells, we use specific kits and reagents. RNA is extracted using RNAiso Plus from Takara, combined with Alladin’s chloroform, isopropanol, ethanol, and diethylpyrocarbonate (DEPC) (C128130, I112011, E111989, D105557). This involves homogenizing in TRIzol (Thermo Fisher Scientific, 15596026), adding chloroform, and undergoing centrifugation to isolate RNA, which is then precipitated and washed before dissolving in DEPC water for concentration measurement.

For reverse transcription, the PrimeScript^™ RT kit with gDNA Eraser from Takara is used for mRNA, while miRNA in rats employs the Ribobio’s Bugle-Loop^™ RT kit (Guangzhou RiboBio Co., Ltd.) (Zhou et al., 2015), following the specific protocols of each. Quantitative PCR is performed using SYBR^® Premix Ex Taq^™ from Takara, starting with a 95°C denaturation, followed by 40 cycles and melting curve analysis.

Expression levels are calculated using the 2 − ΔΔCT method, comparing Ct values between target and reference genes (Livak and Schmittgen, 2001). Rat mRNA primers are detailed in Table 1. Due to intellectual property constraints, rat miRNA primer sequences from the Ribobio kit are not disclosed. Cell miRNA primers are provided in Table 2.

Rats were anesthetized with 10% Sodium Pentobarbital (Wengjiang Reagent, CAS: 57-33-0) via intraperitoneal injection. Fecal samples were collected under sterile conditions and either used immediately or stored at −80°C. Microbial DNA was extracted using the E.Z.N.A.^® Soil DNA Kit (Omega Bio-tek, United States), with quality checked by agarose gel electrophoresis. Concentration and purity were measured with a NanoDrop2000 spectrophotometer (Thermo Scientific, United States). The V3–V4 region of the 16S rRNA gene was amplified using primers 338F and 806R. The PCR mix included TransStart FastPfu Buffer, dNTPs, primers, and DNA polymerase, in a final volume of 20 μL. The PCR process began with denaturation at 95°C for 3 min, followed by 27 cycles at 95°C, 55°C, and 72°C, each for 30 s. The final step was a 10-min extension at 72°C, followed by cooling at 4°C. PCR products were separated on a 2% agarose gel, then purified and quantified. Sequencing libraries were created using the NEXTFLEX Rapid DNA-Seq Kit, which included adapter ligation, magnetic bead purification, and PCR amplification. Sequencing was performed using the Illumina PE300/PE250 platform. The raw paired-end data underwent quality control with Fastp (Chen et al., 2018) (v0.19.6) and were merged using FLASH (Magoč and Salzberg, 2011) (v1.2.11). Clustering into operational taxonomic units (OTUs) was done with UPARSE (Edgar, 2013; Stackebrandt and Goebel, 1994) (v7.1) at 97% similarity, with chimeric and mitochondrial sequences removed. To balance sequencing depth, each sample was normalized to 20,000 reads, achieving a Good’s coverage of 99.09%. OTUs were classified with the RDP Classifier (Wang et al., 2007) (v2.11) at a 70% confidence threshold, referencing the Silva 16S rRNA Gene Database (v138), followed by community composition analysis.

The “data table” package (Barrett et al., 2021) was used to read the OTU table, set OTU ID as row names, and extract sample and taxonomy tables. The “amplicon” package (Callahan et al., 2016) calculated microbial relative abundance with a threshold of 1 × 10⁻⁴.

The “ropls” package (Thévenot et al., 2015) constructed the partial least squares discriminant analysis (PLS-DA) model, and the “ggplot2” and “ggsci” (Xiao, 2023) packages visualized the results, including sample score and variable importance plots.

To observe gut microbial abundance changes, the total abundance of each species was calculated and sorted. The top 10 species by abundance and the total abundance of remaining species were displayed using stacked bar plots with “ggplot2” package. The “aplot” package (Yu, 2023) added group information.

For LEfSe (linear discriminant analysis effect size) analysis, the OTU table was loaded with identifiers as row names. A sample table and formatted taxonomy table were created. The “microeco” package (Liu C. et al., 2021) integrated these tables into a Microtable object. LEfSe analysis was conducted using the trans_diff$new function, with bar plots and a cladogram generated for taxonomic differences. Correlation analysis with the “corrplot” package (Wei et al., 2021) produced a correlation heatmap among gut microbiota.

We utilized the “tidyverse” (Wickham et al., 2019) packages for data preprocessing. The support vector machine-recursive feature elimination (SVM-RFE) algorithm, implemented via the “e1071” package (Meyer et al., 2023), was employed for feature selection, identifying influential features for classification. K-fold cross-validation evaluated the feature set effectiveness, with error rates and model accuracy visualized to assess the impact of feature subsets. The best-performing subset was selected based on error rates. A random forest model, created with the “randomForest” package (Liaw and Wiener, 2002), evaluated feature importance, and its generalizability was tested with five repetitions of 10-fold cross-validation. The optimal feature subset was selected based on importance scores and error rates.

Next, key microbiota identified from the SVM-RFE and random forest analyses were extracted and compared using the “VennDiagram” package (Chen and Boutros, 2011), with overlapping microbiota exported for receiver operating characteristic (ROC) analysis. The “pROC” package (Robin et al., 2011) was used for ROC analysis, calculating the area under the curve (AUC) values. An AUC >0.7 was considered indicative of strong discriminative power for key microbiota. ROC curves were plotted and saved, and data for high AUC microbiota were exported.

Finally, the “data table” and “tidyverse” packages processed expression data of target microbiota. Boxplots were generated with the “ggsignif (Ahlmann-Eltze and Preisenergie, 2021), “ggplot2,” and “ggpubr” (Kassambara, 2023) packages to illustrate abundance differences among groups. The stat_compare_means function added p-values, identifying significant microbiota associated with NBTL treatment for rheumatoid arthritis.

RNA was extracted from tissue samples using TRIzol^® Reagent (Invitrogen) and RNAiso Plus (Takara, 9108). Concentration, purity, and integrity were assessed with a NanoDrop 2000 spectrophotometer and agarose gel electrophoresis. Quality was further verified using an Agilent 2100 Bioanalyzer, ensuring ≥1 μg of total RNA, a concentration of at least 50 ng/μL, and an OD 260/280 ratio within 1.8 to 2.2. diethylpyrocarbonate (DEPC) water (Alladin, D105557) was used in all preparations to ensure RNase-free conditions. mRNA was extracted using Oligo(dT) magnetic beads, then fragmented into ~300 bp pieces with a buffer for Illumina Sequencing. Using the PrimeScript^™ RT reagent Kit with gDNA Eraser (Takara, RR047A), double-stranded cDNA was synthesized from fragmented mRNA with random primers. The cDNA then underwent end repair, A-tailing, and Y-shaped adapter ligation using the Truseq^™ RNA Sample Prep Kit (Illumina). Adapter-ligated cDNA was purified, size-selected, and amplified by PCR to generate the final sequencing library. The libraries were quantified with the QuantiFluor^® dsDNA System (Promega), and sequencing was conducted on the Illumina HiSeq X Ten platform using its reagent kit. Sequencing reads were analyzed for base distribution and quality, followed by quality control, reference genome alignment, and expression level analysis, with results reported in transcripts per million (TPM).

We used the “ropls” package to create a PLS-DA model, extracting results like the model summary and VIP (variable importance in projection) values. Sample distributions and confidence ellipses were then plotted using the “ggplot2” package. Kyoto Encyclopedia of Genes and Genomes (KEGG) gene sets were acquired using the “msigdbr” package (Dolgalev, 2022) to prepare gene lists for gene set variation analysis (GSVA). The analysis was carried out with the “GSVA” package (Hänzelmann et al., 2013), followed by differential analysis using GraphPad Prism.

For weighted gene co-expression network analysis (WGCNA), we processed and cleaned the gene expression data. The “WGCNA” package (Langfelder and Horvath, 2008, 2012) was then used for clustering analysis and the removal of outlier samples. The optimal soft threshold was determined, and a co-expression network was constructed to identify modules. These modules were analyzed for their correlation with phenotypic traits, and relevant plots were generated.

Next, module genes were analyzed for protein–protein interactions (PPI) using the STRING database (Szklarczyk et al., 2019).¹ Cytoscape (Shannon et al., 2003) was used to visualize the interaction data, and the cytoHubba plugin (Bader and Hogue, 2003) identified hub genes through degree analysis. Following this, we read and cleaned the gene expression data, imported the top 10 genes based on degree ranking, and standardized the GSVA results for the vascular endothelial growth factor (VEGF) pathway. The “vegan” package (Oksanen et al., 2022) was used for correlation analysis, and the results were visualized with the “ggcor” package (Huang et al., 2020).

We imported and cleaned the synovial PCR results and microbiota data. Spearman correlation tests were conducted to examine relationships between these datasets, and scatter plots were generated. Lastly, we imported synovial PCR data, selected samples, conducted correlation tests, and plotted scatter diagrams to display the relationship between miRNA and VEGF expression.

Load the “data table,” “mlr3” (Lang et al., 2019), “mlr3learners” (Lang et al., 2023), “mlr3pipelines” (Binder et al., 2021), and “keras” (Chollet and Allaire, 2023) packages, and use the “reticulate” package (Kalinowski et al., 2024) to access the Keras library (Chollet, 2015) from a Python environment. Specify the working directory, load the dataset, and split it into 80% for training and 20% for testing. Standardize the training data and apply principal component analysis (PCA) for dimensionality reduction. Construct a sequential deep neural network (DNN) model with two hidden layers, using rectified linear unit (ReLU) activation and the Adam optimizer. The model was compiled using mean squared error (MSE) as the loss function and mean absolute error (MAE) as the performance metric. It was trained with a 20% validation split and evaluated on the test set. Plot comparisons of predicted versus actual values, calculate the R² value, monitor the loss curve during training, and display a graphical representation of the model structure.

Human umbilical vein endothelial cells (HUVECs) were sourced from the American Type Culture Collection (ATCC) and stored in liquid nitrogen. Before culture, cells were quickly thawed in a 37°C water bath and transferred to a 15 mL centrifuge tube containing 10 mL of pre-warmed Roswell Park Memorial Institute (RPMI)-1640 medium (Gibco, C11875500BT) supplemented with 10% fetal bovine serum (FBS) (Gibco, 10099-141C) and 2 mM L-glutamine. After centrifugation at 1,000 rpm for 5 min, the supernatant was discarded, and the cell pellet was resuspended in fresh pre-warmed medium and transferred to a pre-coated T25 culture flask. HUVECs were maintained at 37°C in a 5% CO₂ atmosphere, with medium changes every 2–3 days. Cells were passaged at 80–90% confluence using 0.25% trypsin-EDTA (Gibco, Wisent 325-043-EL) for 1–2 min at 37°C, followed by centrifugation at 1,000 rpm for 5 min. Digestion was halted with an equal volume of RPMI-1640 containing 10% FBS, and cells were resuspended in fresh complete medium and redistributed into new culture flasks. For cryopreservation, cells were harvested, resuspended in CELLSAVING serum-free cryopreservation solution (C40100), aliquoted into cryovials, and stored at −80°C.

To screen for the effective concentration of NBTL-containing serum, HUVECs were cultured to an appropriate density and seeded into 96-well plates at 5 × 10³ cells per well in RPMI-1640 medium with 10% FBS. Cells were treated with various concentrations of NBTL-containing serum (0, 2, 4, 8, 12, 16, 20%) or an equal volume of blank rat serum for 24 h at 37°C and 5% CO₂. It is important to note that the intervention was carried out solely with rat serum and RPMI-1640 basal medium, without the addition of FBS. Each concentration was tested in triplicate, followed by the addition of 10 μL of Cell Counting Kit-8 (CCK-8) reagent (APExBIO, K1018) and 90 μL of RPMI-1640 medium, with further incubation for 1–4 h. Absorbance was measured at 450 nm using a microplate reader.

HUVECs (1 × 10⁵ cells per well) were cultured in 6-well plates until they reached confluence. A scratch was then created using a sterile 200 μL pipette tip, followed by treatment with either NBTL-containing serum at its optimal concentration or blank rat serum. Cell migration and scratch closure were evaluated by capturing images at the beginning (0 h) and after 24 h using an inverted microscope.

HUVECs (1 × 10⁵ cells per well) were seeded in 6-well plates and labeled with 10 μM carboxyfluorescein diacetate succinimidyl ester (CFSE) (Beyotime, C1031). After 15 min of incubation at 37°C in the dark, the cells were rinsed with phosphate-buffered saline (PBS) (Nanjing Shengxing, SN331). The cells were then exposed to either control rat serum or NBTL-containing serum for 24 h, collected, and resuspended in PBS. Flow cytometry was used to assess cell proliferation, utilizing a 488 nm excitation and 518 nm emission.

HUVECs (1 × 10⁵ cells per well) were first cultured in 6-well plates and treated with either control rat serum or NBTL-containing serum for 24 h. The cells were then transferred to Matrigel-coated 24-well plates (1.5 × 10⁵ cells per well) and incubated for 4 h. Tube formation was documented and quantified using ImageJ software with the Angiogenesis Analyzer plugin.

Cells were seeded at 1 × 10⁵ density into 3.5 cm glass-bottom culture dishes and incubated overnight. The next day, the media was replaced with control or drug-containing serum based on CCK8 results, followed by a 24-h incubation. Cells were fixed with 4% paraformaldehyde (PFA) for 10 min, then washed with pre-cooled PBS and permeabilized using 0.5% Triton X-100 on ice. Blocking was performed with 3% bovine serum albumin (BSA) for 30 min at room temperature. Cells were incubated overnight at 4°C with primary antibodies (VEGF, VEGFR2, HIF-1α, VEGFR1) diluted in PBS with 1% BSA. After washing, secondary antibodies (1:400 dilution) were applied for 1 h at room temperature in the dark. Nuclei were stained with 100 ng/mL DAPI for 10 min, followed by washes. Finally, an anti-fade reagent was added, and the cells were imaged using a fluorescence microscope.

In this study, the analysis and processing of gut microbiota data were primarily conducted using R software (R Core Team, 2021) (version 4.3.0). For the correlation analysis between synovial PCR results and characteristic gut microbiota, R software (version 4.3.1) was employed. The construction and training of the deep neural network (DNN) model utilized R software (version 4.3.0), while other R analyses were mainly performed using version 4.3.2. Grayscale value analysis of western blot results, scratch assay migration area analysis, and angiogenesis assay analysis were all conducted using ImageJ software. CFSE flow cytometry analysis was carried out with FlowJo software (version 10.8.1). Further data processing and graphical plotting were completed using GraphPad Prism software (version 9.5.0).

For statistical comparisons, one-way analysis of variance (ANOVA) was used for comparisons among multiple groups with a single factor, and two-way ANOVA was employed for comparisons involving two factors. Paired t-tests were used for within-group comparisons, while unpaired t-tests were employed for comparisons between independent groups. Before performing statistical tests, normality was assessed using the Shapiro–Wilk test, and homogeneity of variance was evaluated using Levene’s test. If these assumptions were met, parametric tests (ANOVA and t-tests) were applied. For non-normally distributed data, non-parametric tests (e.g., Mann–Whitney U test) were used.

In all statistical analyses, significance levels were set as follows: p-values under 0.05 were considered statistically significant, those below 0.01 were highly significant, and p-values less than 0.001 or 0.0001 were deemed extremely significant.

Using UPHLC-MS-TOF, we conducted an in-depth analysis of the prototype components and metabolites of new bitongling (NBTL) after its entry into the bloodstream.

In the positive ion mode (Figure 1A and Table 3), the primary prototype components identified included zingerone (Rt. min. 5.029, m/z: 177.0907), wogonoside (Rt. min. 7.219, m/z: 461.1054), and uracil (Rt. min. 2.764, m/z: 113.0346), among others. In the negative ion mode (Figure 1B and Table 4), the main components detected were terpineol-4 (Rt. min. 8.989, m/z: 153.12506), succinic anhydride (Rt. min. 13.701, m/z: 99.00613), and salicyluric acid (Rt. min. 6.295, m/z: 194.04202), among others.

Furthermore, the analysis of NBTL metabolites in the bloodstream revealed several key compounds. In the positive ion mode (Figure 1C and Table 5), we identified pyrimidine (Rt. min. 13.328, m/z: 81.04522), pyridine (Rt. min. 0.539, m/z: 80.04992), and myrcene (Rt. min. 4.501, m/z: 137.13211), among others. In the negative ion mode (Figure 1D and Table 6), the primary metabolites detected included D-galactose (Rt. min. 13.649, m/z: 179.0527), arabitol(D) (Rt. min. 13.649, m/z: 151.0578), and (+)-triptophenolide (Rt. min. 11.431, m/z: 311.1644), among others.

This detailed UPHLC-MS-TOF analysis successfully identified the prototype components and metabolites of NBTL in the bloodstream, providing crucial data to support further research into the pharmacological mechanisms of NBTL.

The animals underwent a one-week acclimation period to adjust to the experimental environment. Following this, two immunizations were performed on day 0 and day 7 to establish an arthritis model. The experimental animals were randomly divided into six groups: control, model, low-dose new bitongling (NBTL-L), medium-dose new bitongling (NBTL-M), high-dose new bitongling (NBTL-H), and Tripterygium wilfordii tablets (TWT). Treatment for each group began on day 14 and continued until day 42 (Figure 2A).

Compared to the model group, all new bitongling groups and the TWT group showed reductions in arthritis index, especially the NBTL-H group (p < 0.01) (Figure 2B). Body weight changes during the experiment did not differ significantly across the groups (Figure 2C). Photographs taken at various time points revealed typical arthritis symptoms in the model group, which were noticeably reduced in the new bitongling and TWT groups (Figure 2D).

Histopathological examination of the ankle joints in the model group revealed significant structural damage along with marked infiltration of inflammatory cells. In contrast, both the NBTL-H and TWT groups showed a notable reduction in these pathological changes (p < 0.05 and p < 0.0001) compared to the model group (Figures 2E,F).

Following new bitongling treatment, TNF-α and IL-17 levels were notably decreased in the NBTL-L, NBTL-M, and NBTL-H groups (p < 0.01 and p < 0.0001), with comparable reductions seen in the TWT group (Figures 2G,H). Regarding IL-6, significant reductions were noted in the NBTL-M and NBTL-H groups (p < 0.01 and p < 0.001), but not in the NBTL-L group (Figure 2I).

Western blot analysis showed significantly increased ASC protein expression in the synovial tissues of the model group. After new bitongling treatment, ASC expression showed a dose-dependent decline, with the greatest reduction observed in the NBTL-H group (p < 0.001). The TWT group also showed a similar decreasing trend (Figures 2J,K). Quantitative polymerase chain reaction (qPCR) results further confirmed that ASC mRNA levels in the synovial tissues were significantly reduced in the NBTL-H, and TWT groups (p < 0.0001) (Figure 2L).

In this study, we employed PLS-DA to deeply analyze the gut microbiota structure across various taxonomic levels: order (Figure 3A), family (Figure 3B), genus (Figure 3C), species (Figure 3D), and operational taxonomic units (OTU) (Figure 3E). Experimental data points are marked in purple, red, and green to represent the control group (control), model group (model), and high-dose new bitongling treatment group (NBTL-H), respectively. Each subplot reveals the main directions of data variation through the first and second principal components, with the corresponding variance percentages labeled on the axes. The principal component analysis results clearly demonstrate the distribution characteristics of the three groups within the principal component space, with the NBTL-H treatment group showing significant separation at multiple taxonomic levels.

We also comprehensively analyzed the microbial abundance distribution at the order (Figure 3F), family (Figure 3G), genus (Figure 3H), and species (Figure 3I) levels. The data is displayed using distinctly colored bar charts, highlighting the abundance of each microbial group. The control, model, and NBTL-H groups are marked by different color bands at the chart base, aiding in the comparative analysis of gut microbiota across various taxonomic levels. Notable differences were observed in the proportions of certain microbial groups. For instance, at the order level (Figure 3F), the abundance of Bacteroidales decreased in the NBTL-H treatment group compared to the control and model groups. At the family level (Figure 3G), Prevotellaceae abundance decreased in the treatment group relative to the model group. At the genus level (Figure 3H), UCG-005 abundance significantly declined in the high-dose treatment group. At the species level (Figure 3I), Lactobacillus johnsonii showed increased abundance in the high-dose treatment group relative to the other two groups.

Additionally, we applied LEfSe to detect significant microbial taxa differences among the control, model, and NBTL-H groups. The linear discriminant analysis (LDA) scores indicated the degree of intergroup differences, with higher scores reflecting greater differences. Analysis revealed that f_Mycoplasmataceae, g_Prevotellaceae_Ga6A1_group, s_Metamycoplasma_sualvi, s_uncultured_bacterium_g_norank_f_Eggerthellaceae, and s_uncultured_bacterium_g_Prevotellaceae_Ga6A1_group exhibited high LDA scores, indicating significant abundance differences among groups (Figure 3J). Furthermore, the phylogenetic tree generated based on LDA scores (Figure 3K) illustrated the evolutionary relationships of these significantly different microbial taxa. Each node’s size represents the significance of the corresponding taxon, and color differences reflect different experimental groups.

Correlation heatmap analysis (Figure 3L) revealed significant correlations of f_Mycoplasmataceae with various other microbial taxa, especially with g_Prevotellaceae_Ga6A1_group in the NBTL-H treatment group. Additionally, changes in the correlation patterns of s_Metamycoplasma_sualvi with s_uncultured_bacterium_g_norank_f_Eggerthellaceae and s_uncultured_bacterium_g_Prevotellaceae_Ga6A1_group might reflect ecological interactions or similar physiological metabolic pathways. Notably, the altered correlations of s_Metamycoplasma_sualvi post-treatment suggest that NBTL may impact specific gut microbial interactions. Furthermore, significant correlation changes for s_uncultured_bacterium_g_Prevotellaceae_Ga6A1_group indicate that this uncultured taxon may undergo niche alterations under drug influence.

The SVM-RFE analysis results demonstrate a gradual improvement in the model’s 5-fold cross-validation accuracy (Figure 4A), which reaches its peak at 0.867 with 20 features. Simultaneously, the 5-fold cross-validation error (Figure 4B) decreases to its lowest point of 0.133 at 20 features, indicating optimal model performance at this point. The random forest analysis results show a steady decline in cross-validation error (Figure 4C) as the number of variables increases, reaching its lowest at 24 variables, demonstrating the model’s excellent generalization ability. The variable importance assessment (Figure 4D) highlights significant contributions to the model’s predictive power from species such as s_Metamycoplasma_sualvi and s_uncultured_bacterium_g_norank_f_Eggerthellaceae. The Venn diagram (Figure 4E) reveals a consensus between SVM-RFE and random forest in identifying 19 key microbial taxa, including f_Mycoplasmataceae and g_Prevotellaceae_Ga6A1_group, as potential targets for new bitongling’s regulation of the gut microbiome.

For the 19 key microbial taxa identified by machine learning, ROC analysis assessed their sensitivity and specificity across various classification thresholds. Six taxa had area under the curve (AUC) values exceeding 70%, indicating strong discriminative power. These include g_Prevotellaceae_Ga6A1_group and s_uncultured_bacterium_g_Prevotellaceae_Ga6A1_group with AUCs of 81.25%, s_Metamycoplasma_sualvi with an AUC of 93.75%, f_Mycoplasmataceae with an AUC of 90.62%, g_norank_f_Erysipelotrichaceae with an AUC of 77.34%, and s_uncultured_bacterium_g_norank_f_Eggerthellaceae with an AUC of 87.89% (Figure 4F).

Box plots illustrate the relative abundance differences of these taxa among the control, model, and NBTL-H. In the model condition, the relative abundance of g_Prevotellaceae_Ga6A1_group (Figure 4I) and s_uncultured_bacterium_g_Prevotellaceae_Ga6A1_group (Figure 4L) was markedly higher than in the control group but decreased significantly following NBTL-H treatment. In contrast, the abundance of f_Mycoplasmataceae (Figure 4G), s_Metamycoplasma_sualvi (Figure 4J), and s_uncultured_bacterium_g_norank_f_Eggerthellaceae (Figure 4K) was notably lower in the model group, yet it significantly increased after NBTL-H treatment. However, g_norank_f_Erysipelotrichaceae (Figure 4H) showed significant differences between the NBTL-H treatment and model groups, but not between the control and model groups, resulting in its exclusion.

Considering the high discriminative power of the ROC curves and the significant differences among the three groups, we identified five taxa as characteristic microbial taxa influenced by NBTL: f_Mycoplasmataceae, g_Prevotellaceae_Ga6A1_group, s_uncultured_bacterium_g_Prevotellaceae_Ga6A1_group, s_Metamycoplasma_sualvi, and s_uncultured_bacterium_g_norank_f_Eggerthellaceae.

Transcriptome analysis results (Figure 5A) through PLS-DA show that PC explain 24.6 and 7.59% of the total variation, respectively, clearly distinguishing between the control, model, and NBTL-H groups. gene set variation analysis (GSVA) (Figure 5B) indicates a significant increase in the VEGF signaling pathway in the model group compared to the control group (p < 0.001), while a significant downregulation is observed in the NBTL-H group compared to the model group (p < 0.05).

Further investigation into the impact of NBTL-H on gene expression patterns was conducted using WGCNA. The soft threshold power analysis (Figure 5C) identified the optimal threshold as 20, optimizing the scale-free topology of the network. In the gene clustering dendrogram (Figure 5D), different modules were distinguished by color. The module-trait relationship heatmap (Figure 5E) revealed that the MEyellow module had a significant positive correlation with the NBTL-H group (r = 0.53, p = 0.02), suggesting that NBTL affects the VEGF pathway by regulating genes within this module, which is crucial for understanding the drug’s mechanism.

We also constructed a PPI network for the NBTL-H module genes (Figure 5F), revealing complex interactions among the drug module genes. Using Cytoscape, we visualized the top 10 genes with the highest degrees within the PPI network (Figure 5G), indicating their potential key roles in the mechanism of NBTL. Spearman analysis (Figure 5H) revealed a strong positive association between Rps9 and VEGF pathway activity, indicating a potential involvement of Rps9 in NBTL’s mechanism via the VEGF pathway. Although other genes (such as Rps16, Rpl11, etc.) did not show significant correlations with the VEGF pathway (p > 0.05), they might still be involved in the pathway through more complex mechanisms.

The study also assessed the effects of different doses of NBTL on VEGF mRNA expression (Figure 5I). VEGF mRNA levels were markedly elevated in the model group, whereas all NBTL doses, particularly the high dose, led to a notable reduction. The positive control group treated with Tripterygium wilfordii tablets (TWT) also showed a similar reduction. Additionally, Spearman correlation analysis revealed significant associations between NBTL-H treatment and characteristic gut microbiota. Specifically, the abundance of f_Mycoplasmataceae (Figure 5J, r = −0.76, p = 0.017), s_Metamycoplasma_sualvi (Figure 5L, r = −0.82, p = 0.0066), and s_uncultured_bacterium_g_norank_f_Eggerthellaceae (Figure 5M, r = −0.83, p = 0.0056) showed a negative correlation with VEGF mRNA expression, while g_Prevotellaceae_Ga6A1_group (Figure 5K, r = 0.76, p = 0.017) and s_uncultured_bacterium_g_Prevotellaceae_Ga6A1_group (Figure 5N, r = 0.84, p = 0.0044) showed a positive correlation.

Our findings (Figures6A,H) revealed that in the model group, the levels of miR-20a-5p and miR-223-3p were significantly reduced, likely contributing to increased angiogenesis in the disease state. Treatment with NBTL, especially at high doses (NBTL-H), significantly restored the expression levels of these miRNAs. Likewise, expression recovery was observed in the positive control group receiving tripterygium glycosides, showing a similar trend. Spearman correlation analysis further demonstrated that, in the high-dose NBTL group, miR-20a-5p and miR-223-3p were significantly negatively correlated with VEGF mRNA expression (Figure 6B, r = −1, p = 5.6 × 10⁻⁶; Figure 6I, r = −0.92, p = 0.0013), suggesting these miRNAs may inhibit VEGF expression and, consequently, angiogenesis in the disease context.

Additionally, Spearman correlation analysis indicated that in the NBTL-H group, the abundance of g_Prevotellaceae_Ga6A1_group (Figure 6D, r = −0.82, p = 0.0066) and s_uncultured_bacterium_g_Prevotellaceae_Ga6A1_group (Figure 6G, r = −0.82, p = 0.0066) were negatively correlated with miR-20a-5p, suggesting that the reduction of these bacterial groups may promote the expression of miR-20a-5p. Conversely, the abundance of s_Metamycoplasma_sualvi (Figure 6E, r = 0.83, p = 0.0053) showed a positive correlation with miR-20a-5p, indicating that the increase of this beneficial bacterium might enhance the production of miR-20a-5p. Although the abundance of f_Mycoplasmataceae (Figure 6C, r = 0.27, p = 0.48) and s_uncultured_bacterium_g_norank_f_Eggerthellaceae (Figure 6F, r = 0.44, p = 0.24) was positively correlated with miR-20a-5p, the correlations were not significant, but they still hold potential biological relevance.

Similarly, in the NBTL-H group, g_Prevotellaceae_Ga6A1_group (Figure 6K, r = −0.78, p = 0.013) and s_uncultured_bacterium_g_Prevotellaceae_Ga6A1_group (Figure 6N, r = −0.76, p = 0.017) both exhibited a strong negative correlation with miR-223-3p. Conversely, s_Metamycoplasma_sualvi (Figure 6L, r = 0.71, p = 0.032) was positively correlated with miR-223-3p. Although the correlations of f_Mycoplasmataceae (Figure 6J, r = 0.45, p = 0.22) and s_uncultured_bacterium_g_norank_f_Eggerthellaceae (Figure 6M, r = 0.66, p = 0.053) with miR-223-3p were not significant, they still possess potential biological importance.

We initially employed a deep neural network (DNN) to predict how NBTL inhibits VEGF-driven angiogenesis by analyzing five characteristic gut microbiota features. The training process (Figure 7A) showed a rapid decline and subsequent stabilization in the model’s loss and error rates, with validation data performing slightly worse than the training data, indicating good generalization. Further optimization (Figure 7B) revealed an R² of 0.671 for the training set and a high R² of 0.958 for the validation set, with a loss of 0.01807 and a MAE of 0.12241, underscoring the model’s robust predictive capability regarding NBTL’s regulation of the VEGF pathway. The model architecture, from input to output layers, effectively captured complex relationships within the data through multiple layers of feature extraction (Figure 7C).

Building on the DNN model’s predictions, we further validated NBTL’s inhibitory effects on VEGF signaling using a rat model. Western blot analysis of VEGF, VEGFR1, VEGFR2, and HIF-1α protein expression demonstrated a dose-dependent inhibitory effect of NBTL across low, medium, and high dosage groups (Figures 7D–G). The high-dose group showed a highly significant suppression of VEGF, VEGFR1, VEGFR2, and HIF-1α expression (p < 0.05 and p < 0.0001), with the medium dose also significantly reducing VEGF, VEGFR1, and HIF-1α levels (p < 0.05 and p < 0.0001) (Figure 7H). In conclusion, NBTL not only exhibited strong potential in inhibiting the VEGF signaling pathway in predictive models but also demonstrated substantial efficacy in rat models, particularly at higher doses, supporting its potential use in treating rheumatoid arthritis.

In this study, we utilized the CCK-8 assay to evaluate the effects of NBTL medicated serum on the viability of HUVECs. The results indicated a dose-dependent decrease in cell viability as the concentration of NBTL medicated serum increased from 0 to 20%, demonstrating an inhibitory effect on cell growth (Figure 8A). At the highest concentration (20%), cell viability significantly decreased relative to the 0% group (p < 0.05). A comparison between the control group (20% blank serum) and the NBTL medicated serum group (20%) further confirmed that NBTL significantly reduced cell viability (p < 0.05) (Figure 8B), establishing 20% as the optimal inhibitory concentration for subsequent experiments.

Next, we assessed the impact of NBTL medicated serum on HUVEC migration using a scratch assay. After 24 h of incubation, the control group displayed significant wound closure, while the NBTL-treated group showed limited closure, indicating a substantial inhibition of cell migration (Figures 8C,D, p < 0.001). This suggests that NBTL may exert its anti-angiogenic effects, in part, by inhibiting cell migration.

We also conducted a tube formation assay to investigate the effect of NBTL medicated serum on angiogenesis. Microscopic observations revealed that HUVECs in the control formed more tube-like structures than in the NBTL group, where these structures were notably fewer (Figures 8E–H, p < 0.05 and p < 0.01), confirming its anti-angiogenic properties.

Additionally, we analyzed the effect of NBTL medicated serum on HUVEC proliferation using CFSE staining and flow cytometry. The results demonstrated a marked decrease in cell proliferation in the NBTL group relative to the control (Figures 8I–N, p < 0.01), supporting its role in inhibiting cell growth.

Quantitative PCR (qPCR) analysis assessed the expression of miR-20a-5p and miR-223-3p, two miRNAs linked to angiogenesis inhibition. The data showed that NBTL-treated serum significantly upregulated these miRNAs (Figures 8O,P, p < 0.05 and p < 0.01), indicating a possible mechanism through which NBTL inhibits angiogenesis.

Finally, immunofluorescence staining demonstrated that NBTL medicated serum significantly downregulated the expression of VEGF, VEGFR1, VEGFR2, and HIF-1α proteins, which are key components of the VEGF signaling pathway in HUVECs (Figures 8Q–T).

In conclusion, this study systematically demonstrates that NBTL medicated serum significantly inhibits HUVEC proliferation, migration, and angiogenesis through multiple mechanisms, offering a potential therapeutic application in anti-angiogenesis treatment.