Discovery and Validation of Traditional Chinese and Western Medicine Combination Antirheumatoid Arthritis Drugs Based on Machine Learning (Random Forest Model)

The combination of traditional Chinese medicine (TCM) and Western medicine is a promising method for treating rheumatoid arthritis (RA). Combining the two fully exploits the advantages of Western and TCM to treat RA and has the potential to greatly improve the therapeutic effect on RA. In this study, we developed a combination drug training set by using 16 characteristic variables based on the characteristics of small molecules of TCM ingredients and Food and Drug Administration-certified combination drug data downloaded from the DrugCombDB database. Furthermore, we compared the prediction and classification abilities of five models: the k-nearest neighbors, naive Bayes, support vector machine, random forest, and AdaBoost algorithms. The random forest model was selected as the classification and prediction model for Western and TCM and Western combination drugs. We collected data for 41 small molecules of TCM ingredients from the Traditional Chinese Medicine Systems Pharmacology database and 10 small molecule drugs commonly used in anti-RA treatment from the DrugBank database. Combinations of Western and TCM for anti-RA treatment were screened. Finally, the CellTiter-Glo method was used to determine the synergy of these combinations, and the 15 most predicted drug combinations were carried out experimental verification. Myricetin, rhein, nobiletin, and fisetin had high synergy with celecoxib, and rhein had high synergy with hydroxychloroquine. The preliminary findings of this study can be further applied for practical clinical anti-RA combined treatment strategies and serve as a reference for clinical treatment of RA with integrated Western and TCM.


Introduction
Modern medical research has demonstrated that many diseases often present complex pathogenic mechanisms, and consequently, the past methods for treating diseases with drugs are no longer adequate. In recent years, drug combination has been reported to offer clear advantages in treating tumors, AIDS, cardiovascular diseases, osteoarthritis, and other complex diseases; it has been widely implemented in clinical practice [1][2][3][4], with treatment methods based on combined drugs demonstrating advantageous synergistic effects, or a "1 + 1 ≥ 2" effect. The therapeutic effects of two or more drug combinations are greater than those achieved when each drug is used alone. Therefore, developing and evaluating potentially synergetic combination drugs, whether for development of new drugs or clinical treatment optimization, are essential. However, most existing combination drugs have been developed on the basis of limited clinical experience or experimental strategies. Verified combination drugs with synergistic therapeutic effects currently comprise a small portion of listed drugs. Therefore, making full use of knowledge of existing drug combinations, discovering new drug combinations based on well-accepted and understood drugs, and establishing a model for evaluating the synergistic effects of drug combinations are key concerns that must be addressed.
Rheumatoid arthritis (RA) is a chronic, systemic, inflammatory disease of the synovial and other joints. RA is a chronic progressive autoimmune disease [5], and its underlying pathogenesis is currently unclear. However, studies have reported that the incidence of RA increases with age due to infection and inflammatory mediators, and approximately 0.3%-1.0% of the global population is affected by the disease every year [6,7]. At present, the drugs used in clinical treatment of RA are mainly nonsteroidal antiinflammatory drugs (NSAIDs), disease modifying antirheumatic drugs (DMARDs), and biological macromolecular therapy. Because RA is complex, several studies have explored different combination therapy strategies for its treatment [8][9][10]. However, some therapeutic drugs, such as methotrexate, may produce unfavorable side effects [11,12]. Therefore, if sole use of anti-RA combination therapy has the potential to produce unfavorable side effects, longterm use in patients with RA may do greater harm. Traditional Chinese medicine (TCM) has long been used to treat RA, and several studies have found that combined anti-RA TCM and a compound prescription have not only antiinflammatory, analgesic, immune regulation, and multilevel and multilink therapeutic effects but also the advantages of high safety, few adverse reactions, and low cost. Such a combination has gained attention in studies of RA treatment and has increasingly attracted international attention [13,14].
Therefore, we propose that the combination of Western medicine and TCM may be a promising method for improving the efficacy of anti-RA treatments and reducing side effects. Through their combination, the advantages of both Western medicine and TCM can be fully exploited to improve the therapeutic effects of RA treatment. However, because numerous active ingredients in TCM and hundreds of listed Western medicines have been reported to be effective in anti-RA treatments, the potential drug combinations are in the thousands or even tens of thousands. This multitude of combinations cannot feasibly be evaluated or verified through clinical practice and trials; therefore, the most practical method for identifying potential combinations would be to use drug data for existing combination drugs to model, evaluate, and screen potential anti-RA combinations of Western medicine and TCM and then experimentally validate these models. Therefore, in this study, we developed a classification and prediction model for combining Western and TCM drugs based on characteristics obtained through small molecule research on TCM and existing data on combined drugs. This model was developed using machine learning modeling followed by manual screening of anti-RA combination drugs. Through experimental verification, we identified effective combinations of Western and TCM drugs for treating RA; these findings can serve as a reference for clinical treatment of RA by using such combined treatment. The flow chart illustrating the process of drug discov-ery and verification for anti-RA integrative medicine through machine learning modeling is displayed in Figure 1.

Data Collection and Collation.
The DrugCombDB database (http://drugcombdb.denglab.org/) is a comprehensive combinatorial drug database [15] containing integrated drug combinations from various data sources. To ensure the reliability of the research data, we first downloaded the Food and Drug Administration-(FDA-) provided known drug combination dataset from the database (up to May 31, 2019). The molecular formula, Chemical Abstracts Service (CAS) number, canonical simplified molecular-input lineentry system (SMILES), and other chemical informatics data corresponding to each combination of drugs and their action targets were collected from the DrugBank database (https:// go.drugbank.com/). Furthermore, we searched "Rheumatoid Arthritis" on the Traditional Chinese Medicine Systems Pharmacology (TCMSP) database (https://tcmspw.com/) to retrieve target Chinese medicine ingredient data reported to be relevant to treating RA. We obtained the main small molecules chemical information from the PubChem database (https:// pubchem.ncbi.nlm.nih.gov/). In addition, several common anti-RA drugs were selected from the DrugBank database along with those discussed in the literature and those recommended by clinicians. All small molecule data for the TCM ingredients (collected from the TCMSP) and for known anti-RA Western drugs were combined in pairs to form the final prediction sample set. Canonical SMILES of all drugs were collected from the DrugBank and PubChem databases, and the 1024-dimensional molecular fingerprints for each of the small molecule drugs were calculated using DRAGON software (version 7.0). The fingerprint similarity (S FP ) between two small molecule drugs was then calculated using Tanimoto coefficient-based similarity:

Construction of Combined
where A is the 1024-dimensional molecular fingerprint feature vector corresponding to the small molecule drug d A , k Ak is the vector length, B is the 1024-dimensional molecular fingerprint feature vector corresponding to the small molecule drug d B , kBk is the vector length, and A•B is the vector inner product.

Drug Combination Characteristics Based on Small
Molecule ADME Similarity. To predict the adsorption, distribution, metabolism, and excretion (ADME) characteristics of each drug, Canonical SMILES for all small molecule drugs were imported into the SwissADME online system (http:// www.swissadme.ch/) [16]. We selected 12 indicators as characteristic indexes for our ADME similarity calculations: 2  BioMed Research International   molecular weight, heavy atoms, aromatic heavy atoms, fraction Csp3, rotatable bonds, H-bond acceptors, H-bond  donors, molecular refractivity, topological polar surface area, Silicos-IT Log P, and Silicos-IT LogSw. All ADME indexes were first normalized, and the Euclidean distance was then used to calculate the similarity value S ADME of ADME: where x A i is the ADME index corresponding to the drug small molecule d A and x B i is the ADME index corresponding to the drug small molecule d B (n = 12).

Drug Combination Characteristics Based on Sequence
Similarity of Small Molecular Targets. The action targets of all drugs were obtained from DrugBank and the TCMSP, and the protr package (version 1.6.2) based on R language was used to obtain the sequence of action protein targets of all small molecule drugs. The sequence similarity between the two target proteins, S Seq hd A , d B i, was calculated using the Smith-Waterman algorithm. The maximum, minimum, median, and mean values of the sequence similarity S Seq h d A , d B i between all target proteins of the d A and d B were included in the characteristic variables as S Seqmax , S Seqmin , S Seqmed , and S Seqmean , respectively.

Drug Combination Characteristics Based on GO
Functional Similarity of Small Molecular Targets. To calculate the gene ontology (GO) functional similarity S GO hd A , d B i of the two drugs d A and d B , we used the GOSemSim package (version 2.14.2). Similarly, the minimum, maximum, median, and mean values of the GO functional similarity S GO hd A , d B i of the two drugs d A and d B were included in the characteristic variables as S GOmin , S GOmax , S GOmed , and S GOmean , respectively.

Drug Combination Characteristics Based on Pathway
Similarity of Small Molecule Targets. We used the cluster-Profiler package (version 3.14.3) to perform a Kyoto Encyclopedia of Genes and Genomes (KEGG) Pathway enrichment analysis on all drug small molecule targets. We then screened the enrichment pathways of all targets based on p value < 0.05 and q value < 0.05. A vector, P = ðp 1 , p 2 , ⋯,p n Þ, consisting of all pathways was constructed according to the enriched pathways, fp k jk = 1, 2,⋯,ng (where n is the number of pathways). If a target was enriched in a pathway p k in the target set T A of the drug small molecule d A , then the value was set as p A k = 1 at the corresponding position of the drug small molecule A pathway vector P A , and if not, it was set as 0: The target pathway vector P B of the drug small molecule d B can be similarly obtained using a string similarity formula to calculate the pathway similarity value (S Pathway ) of two drugs: where kP A k and kP B k are path vector lengths and P A · P B is the vector inner product.
2.2.6. Drug Combination Characteristics Based on Distance of Small Molecule Targets in Human PPI Network. The distance and network proximity of small molecule targets of different combinations of drugs on the human protein-protein interaction (PPI) network were calculated to measure the synergy characteristics of combination drugs within the PPI network: where d ða, bÞ is the shortest path distance between two protein targets, a and b, within the PPI network (calculated using the igraph package [version 1.2.6]), hd AB i = 1/kAk × kBk∑ a∈A, b∈B dða, bÞ, kAk is the number of targets for the drug small molecule d A , and kBk is the number of targets for the drug small molecule d B . The human PPI background network used in the calculation was derived from the literature [17]; the PPI network contained a total of 16677 proteins and 243603 interactions.
The Tanimoto coefficient-based similarity, Euclidean distance, and cosine similarity calculations used in this paper were calculated using the philentropy package (version 0.4.0) based on R language (version 4.0.2).

BioMed Research International
Furthermore, we used established sample training sets to compare the classification prediction effects of the five machine learning models. First, the sample set was divided into prediction training sets and internal validation sets at a 7 : 3 ratio. Then, the prediction training set was input into the five machine learning models, and the machine learning models were trained using three sets of 10-fold crossvalidation. The trained machine learning models were evaluated using the internal validation sets, and the model with the best classification prediction effect was selected as the final drug prediction model for the Western and TCM drug combinations.

Evaluation of Machine Learning Models.
In this study, we defined precision rate as the proportion of positive samples correctly predicted by the model out of all positive pre-dictions, recall rate as the proportion of positive samples correctly predicted by the model out of all predictions, accuracy as the proportion of positive and negative samples to all samples correctly predicted by the model, and F1-score as the overall score of the model when precision rate and recall rate were given equal weight (calculated using Formulas (6)-(9), respectively). A higher score indicated a better classification effect.
where true positive (TP) indicates a sample is positive and was predicted to be positive, false negative (FN) indicates a sample is positive but was predicted to be negative, false positive (FP) indicates a sample is negative but was predicted to be positive, and true negative (TN) indicates a sample is negative and was predicted to be negative.
In the plot of the receiver operating characteristic (ROC) curve, the horizontal axis represents the false positive rate (FPR) and the vertical axis represents the true positive rate (TPR). The area under the ROC curve is the area surrounded by the ROC and its lower coordinate axis and is represented as AUROC in this paper. In addition, because the precision and recall rates of the model are not fixed, when the classification threshold of the model is adjusted, the precision and recall rates continue to change. The precision-recall curve represents the relationship between the precision rate of the ordinate and the recall rate of the abscissa with the change of the classification threshold p value; the area beneath this curve is represented as the area under the precision-recall curve (AUPRC).        BioMed Research International counted and inoculated into 96-well culture plates at 2000 cells/100 μL/well; 90 μL of culture medium (including serum) was added to each well for overnight culture, and 10 μL of drug solution was added. After a 48 h culture, CTG was used to detect the proliferation of RAW264.7 cells; CTG (100 μL) was added to each well, the content of which was incubated at room temperature in a dark environment for 10 min. The chemiluminescence value of each well was measured at 500 ms using enzyme calibration.

Combination Index Measurement.
In this study, according to the principle of the median effect-specifically by employing the Chou-Talalay model [1]-the synergistic effect between each set of drugs was analyzed, and the effects of integrated Western and TCM drugs on RA cell lines were evaluated using the combination index (CI). The CI indicates the degree of drug synergy, addition, or antagonism for any given drug combination: where when CI = 1, the two drugs have an additive effect, when CI < 1, the two drugs have a synergistic effect, and when CI > 1, the two drugs are antagonistic. D A and D B represent the concentrations of the two drugs, d A and d B , when administered in combination, and ðD A Þ x and ðD B Þ x represent the concentrations of the two drugs, d A and d B , when they have been administered alone and the combined drug inhibition rate has been reached.
In Formula (11), f a and f u are the inhibition and survival rates, respectively, when the two drugs are combined, f a = 1 − f u , D m is the median dose, and m is the curve coefficient (calculated according to half-maximal inhibitory concentration (IC 50 ) theory).
In this study, CompuSyn software (http://www .combosyn.com/) was used to calculate the CI values of the drug combinations to evaluate whether each screened combination of Western and TCM drugs had synergistic therapeutic effects.
2.6. Statistical Analysis. Data were expressed as the mean ± standard deviation, and a t-test and Mann-Whitney nonparametric tests were performed using IBM SPSS Statistics (version 26.0).

Collection and Collation of Combined Drug Data.
We obtained 946 pairs (comprising 816 drugs) of combinatorial drug data reported by the FDA from the DrugCombDB database. We further downloaded all existing small molecule drug data from the DrugBank database. After name matching and sorting, some small molecule drug data without key information, such as targets, were removed. In total,

Construction of Combination Drug
Training Set. First, using the method described in Section 2.2, we obtained the molecular fingerprint similarity value (S FP ), the ADME index similarity value (S ADME ), and the minimum, maximum, median, and mean of the similarity (S Seqmin , S Seqmax , S Seqmed , and S Seqmean , respectively) of the drug target sequence between the two small molecule drugs. The minimum, maximum, median, and mean values of the GO functional similarity of the drug targets were S GOmin , S GOmax , S GOmed , and S GOmean , respectively. The KEGG pathway similarity value (S Pathway ); the minimum, maximum, median, and mean of the shortest path of different drug targets on the human PPI network (S PPImin, S PPImax , S PPImed , and S PPImean , respectively); and the proximity (S np ) of the different drug targets were also obtained. Thus, we obtained a total of 16 characteristic variables to predict combinations of Western and TCM drugs.
In this study, we evaluated 16 characteristic variables between two drugs for 546 pairs of combined drugs and removed those for which characteristic variables could not be calculated or were missing. We obtained a total of 404 positive samples of combined drugs (comprising 488 drugs) reported by the FDA. Due to a lack of effective noncombination drug data, we combined the 488 drugs in pairs, removed the existing combination drug samples and those missing data, and obtained a total of 45603 noncombination drugs as negative samples (Supplementary 2; positive and negative sample training dataset for all drug combinations).
Finally, to construct the training set for the TCMcombined drug samples, because the number of noncombined drugs should be much larger than that of combined drugs and to avoid large deviation in positive and negative samples, we set the ratio of positive to negative samples as 1 : 5; the number of positive samples was 404, and the number of negative samples randomly extracted from the 45603 noncombined drug samples was 404 × 5 = 2020. Through

Analysis and Screening of Characteristic Variables in
Combination Drug Training Set. A t-test and a Mann-Whit-ney nonparametric test were used to determine whether significant differences existed between the 16 characteristic variables in the combined drug group and those in the noncombined drug group (p < 0:05). The results revealed that S FP , S ADME , S Seqmin , S GOmax , S GOmed , S GOmean , and S PPImax did not significantly differ between the combined drug and noncombined drug groups (p > 0:05). Significant differences 9 BioMed Research International were found in the S Seqmax , S Seqmed , S Seqmean , S Pathway , S PPImin , S PPImed , S PPImean , and S np between the two statistical tests (p < 0:05). Statistical differences were found in S GOmin in the nonparametric tests (p < 0:05). All statistical tests were performed in IBM SPSS Statistics (version 26.0), and their data are presented in Table 1. The violin diagrams of the distribution of characteristic variables in different groups were completed using ggplot2 (version 3.3.5) and are displayed in Figure 3.
An rfe function using the caret package was used to screen the 16 characteristic variables. A random forest model was used as the selection function, which was f uncti ons = rf Funcs, and the method was cross-validation, which was method = " cv: " The others parameters were set to default. The training set was repeated 1000 times to extract feature variables from 1000 randomly generated training set samples, and each extracted feature variable set was recorded. The frequency of each feature variable in the 1000 feature selections was counted; the results are displayed in Figure 4.
As illustrated in Figure 4, eight characteristic variables, namely, S Pathway , S GOmin , S GOmax , S Seqmin , S Seqmax , S Seqmean , S PPImean , and S np , had appearance frequencies of over 70% after 1000 characteristic screenings. S Pathway and S GOmin were identified as key characteristic variables for predicting combined drugs. The other eight indicators had appearance frequencies below 70%; for example, S FP had a rate of less than 10%. Thus, we removed the group of indicators with frequencies lower than 70% and included the top eight indicators as feature variables for subsequent machine learning.

Comparison and Selection Results of Machine Learning
Models. The five included machine learning models were subjected to three replicates of 10-fold cross-validation under 1000 random samples to evaluate the predictive classification performance of different models. The results are displayed in Table 2, and the ROC and PR curves of the five models are presented in Figure 5.
Through evaluation and comparison of the predictive ability of the different models, we found that the classification index scores of the random forest model were generally better than those of the other prediction classifiers. Therefore, the random forest model was adopted as the combined drug prediction classifier for this study. The parameter was represented as ntree = 300.

Prediction
Results of Integrative Medicine. We entered "Rheumatoid Arthritis" into the TCMSP database as a search term and identified 41 small molecule types in TCM ingredients that have been used to treat RA. We collected the target treatment information for each small molecule and verified it using the Uniprot database (https://www .uniprot.org/). We collected the molecular formula, Pub-Chem CID, CAS number, and canonical SMILES data of each small molecule from PubChem ( Table 3).
The information for drugs commonly used in RA treatment was also obtained from DrugBank. In this study, acetaminophen, celecoxib, choline salicylate, diclofenac, hydroxychloroquine, indometacin, leflunomide, methotrexate, naproxen, sulfasalazine, and 10 other common drugs were selected, as presented in Table 4.
The aforementioned 41 small molecule TCM (collected from TCMSP) and 10 common anti-RA Western medicine   The results revealed that out of 1000 random sample predictions, 34, 62, and 114 pairs of Western and TCM drugs were predicted to be combination drugs more than 500, 100-500, and 1-100 times, and 200 groups were predicted to be noncombination drugs. In this study, we considered predictions of more than 100 times to be potential anti-RA Western and TCM combination drugs, with higher numbers of positive sample predictions indicating a greater potential for drug combination. Therefore, we selected the 15 groups that were predicted more than 700 times for subsequent cell experiments. The combination drug information is displayed in Table 5 (Supplementary 4; prediction results of all Western and TCM combinations).
3.6. Experimental Verification Results. Because RAW264.7 cells were used for experimental verification, we first determined whether a single small molecule had an inhibitory effect on RAW264.7 cells. We then tested the Western and TCM drug combinations. Finally, the CI value was calculated using CompuSyn software to determine whether the small molecule combinations of the two drugs were synergetic.
We then applied the method described in Section 2.5. Based on observations of the proliferation effect of the 15 small molecules on the RAW264.7 cells, myricetin, nobiletin, fisetin, celecoxib, and hydroxychloroquine significantly reduced the survival rate of the RAW264.7 cells. The IC 50 values of the myricetin, nobiletin, fisetin, celecoxib, and hydroxychloroquine were> 100, 23.93, 31.19, >100, and 79.43 μmol/L, respectively, as presented in Figure 6.
We further tested combinations of Western and TCM containing these five small molecules (a total of seven pairs) based on the prediction results of the machine learning model and calculated the CI values, as displayed in Table 6. We found that of the seven pairs of Western and TCM combinations, five of the pairs had CI values less than 0.7, indicating that these drugs have good synergy.

Discussion
Machine learning methods have been widely used in many fields of research, such as text classification [23], sentiment analysis [24], breast cancer diagnosis [25], web page classification [26], and image generation [27]. Many researchers have used computational and modeling methods to identify potential combinatorial drugs [17,[28][29][30], and the use of machine learning to predict drug combinations has notably increased in recent years [31][32][33][34]. However, nearly all of this research has focused of combinations of Western drugs. Combinations of Western and TCM drugs have rarely been reported, especially within the context of RA treatment. TCM is a complex system composed of many small molecules with synergistic effects. In recent years, the research of TCM small molecule is still the focus of modern pharmacology. We mainly focused on the prediction of the combination of small molecule with clear pharmacological effects in TCM and Western medicines that have been approved drugs. Therefore, in this study, we examined the small molecule characteristics of TCM ingredients that differed from those in Western combination drugs. We developed a method for complete modeling and prediction to screen anti-RA combination drugs.
For feature variable construction, we selected drug structures and ADME features. However, the results of the statistical analyses and feature variable screening revealed that 12 BioMed Research International drug structure similarity and ADME characteristics had no significant effect on the classification of combination drugs. Notably, sequence similarity, GO functional similarity, and the distance of drug targets within the PPI network of the small molecule drug targets were found to be related to the combination of drugs, with the similarity value of the drug target pathway (S Pathway ) being identified as a key feature in combination drugs, both statistically and in feature screening. Therefore, the synergistic mechanism of two or more drugs within a pathway may serve as a key reference for future developments of combined drugs. In constructing the sample training set, we were unable to collect sufficient noncombined drug data (namely the number of negative samples) due to limited research data. Therefore, in the original data, the number of positive samples was only 1% that of negative samples, which created a serious imbalance in the number of positive and negative samples. If the model had been constructed using the original datasets, the disproportionate number of negative samples would have caused overfitting, and the predictions would have been biased toward the classifications of the larger number of samples. This would have greatly reduced the normalizing ability of the models. However, simply reducing the sample size may have excluded potential drug combinations in the 45603 samples.
To effectively and reasonably overcome this challenge, we proposed a method in which multiple training sets were constructed using random sampling performed 1000 times. For each instance of modeling, the positive samples remained unchanged, and a limited number of negative samples were randomly selected from the total negative samples; 1000 training sets were randomly composed. The machine learning model was trained by and predicted the 1000 training sets. The final potential combination drugs were then determined on the basis of the percentage of each sample predicted to be a combination drug (positive) in the 1000 learning sessions. However, because the number of noncombination drugs is generally much larger than that of combination drugs but large deviations in numbers of positive and negative samples would cause prediction deviation in the machine learning model, in this study, the ratio of positive and negative samples was set to 1 : 5. Through these changes, we were able to prevent imperfect data sampling of combined drugs.
In this study, from the prediction and evaluation results of five different machine learning classification models in Rosmarinic acid+celecoxib >100 + 100 CI > 1 F1 score, AUROC, and AUPR values of RF are the highest among the five machine learning models. AdaBoost performs best on recall, followed by SVM and RF. The prediction scores of k-NN and NB on these six machine learning evaluation indexes were lower, suggesting that these two models may not be suitable for combination drug prediction. The prediction performance of SVM and AdaBoost is second only to that of RF. Therefore, through the evaluation and comparison of the predictive ability of different models, we found that the RF classification index scores are generally better than other predictive classifiers. Therefore, the RF model is used as the final combined drug prediction classifier, where the parameter ntree = 300. We used the RF classifier to predict 410 groups of Western and TCM drug combinations to treat RA and thus preliminarily predicted 96 pairs of potential Western and TCM drug combinations to treat RA. We used RAW264.7 cells to experimentally verify the effects of a single drug small molecule on the proliferation of RAW264.7 cells. We then verified the effects of combined drugs on the cells. The results revealed that five of the seven selected pairs of combined Western and TCM drugs had high synergistic effects, indicating that our model had high accuracy and practicability. The five optimal combinations were Rhein and Celecoxib, Nobiletin and Celecoxib, Myricetin and Celecoxib, Fisetin and Celecoxib, and Rhein and Hydroxychloroquine.
Notably, celecoxib was found to have a higher synergistic therapeutic effect when combined with several small molecule TCM ingredients, namely, myricetin, rhein, nobiletin, and fisetin. In addition, rhein did not demonstrate an inhibitory effect on RAW264.7 cell proliferation when applied alone. However, when combined with celecoxib (CI < 0:7), the combined inhibition rate of 100 μmol/L rhein and 20 μmol/L celecoxib was 35.68%, which was higher than the rate of each of the two alone. Rhien had high synergy with hydroxychloroquine.
Notably, triptolide combined with Western drugs, such as methotrexate, has previously been reported to produce favorable results in past experiments [35]. However, the number of predicted combinations containing triptolide in our model was not high; it predicted triptolide combined with diclofenac (75 times), naproxen (75 times), and celecoxib (24 times) but did not predict triptolide combined with methotrexate (0 times). This may be due to triptolide's reported clinical toxicity [36]. Although triptolide has been clinically reported to improve treatment outcomes, it may not be the best option for drug combinations and may even lead to more systemic side effects [37]. Our prediction results may support this hypothesis. The screening of potential combination of TCM and Western medicine based on the machine learning combined with experimental validation methods can improve efficiency and reduce research costs.
In conclusion, we constructed combination drug characteristics based on existing combination drug data and characteristics from TCM drug research. Key characteristic indexes were screened, and a random forest classification model was used to predict potential combinations of Western and TCM drugs for anti-RA treatment. Through exper-imental verification, we preliminarily demonstrated that this method of combination therapy for anti-RA treatment can be applied to clinical practice. In the future, we will conduct vivo experiments and clinical trials to validate the findings of this study. In addition, our machine learning model for integrated Western and TCM drug combination can be applied for predicting not only anti-RA drug combinations but also TCM combinations for other diseases.
Abbreviations ADME: Adsorption, distribution, metabolism, and excretion AUPRC: Area under the precision-recall curve CAS: Chemical Abstracts Service CI: Combination index DMARDs: Disease modifying antirheumatic drugs FDA: Food and Drug Administration FN: False negative FP: False positive FPR: False positive rate GO: Gene ontology KEGG: Kyoto Encyclopedia of Genes and Genomes k-NN: k-nearest neighbors NB: Naive Bayes NSAIDs: Nonsteroidal anti-inflammatory drugs PPI: Protein-protein interaction RA: Rheumatoid arthritis RF: Random forest ROC: Receiver operating characteristic SMILES: Canonical simplified molecular-input line-entry system SVM: Support vector machine TCM: Traditional Chinese medicine TCMSP: Traditional Chinese Medicine Systems Pharmacology TN: True negative TP: True positive TPR: True positive rate.

Data Availability
The data used to support the findings of this study are included within the article.