Skip to main content

Plasma metabolomic profiling in patients with rheumatoid arthritis identifies biochemical features predictive of quantitative disease activity

Abstract

Background

Rheumatoid arthritis (RA) is a chronic, autoimmune disorder characterized by joint inflammation and pain. In patients with RA, metabolomic approaches, i.e., high-throughput profiling of small-molecule metabolites, on plasma or serum has thus far enabled the discovery of biomarkers for clinical subgroups, risk factors, and predictors of treatment response. Despite these recent advancements, the identification of blood metabolites that reflect quantitative disease activity remains an important challenge in precision medicine for RA. Herein, we use global plasma metabolomic profiling analyses to detect metabolites associated with, and predictive of, quantitative disease activity in patients with RA.

Methods

Ultra-high-performance liquid chromatography-tandem mass spectrometry (UPLC-MS/MS) was performed on a discovery cohort consisting of 128 plasma samples from 64 RA patients and on a validation cohort of 12 samples from 12 patients. The resulting metabolomic profiles were analyzed with two different strategies to find metabolites associated with RA disease activity defined by the Disease Activity Score-28 using C-reactive protein (DAS28-CRP). More specifically, mixed-effects regression models were used to identify metabolites differentially abundant between two disease activity groups (“lower”, DAS28-CRP ≤ 3.2; and “higher”, DAS28-CRP > 3.2) and to identify metabolites significantly associated with DAS28-CRP scores. A generalized linear model (GLM) was then constructed for estimating DAS28-CRP using plasma metabolite abundances. Finally, for associating metabolites with CRP (an indicator of inflammation), metabolites differentially abundant between two patient groups (“low-CRP”, CRP ≤ 3.0 mg/L; “high-CRP”, CRP > 3.0 mg/L) were investigated.

Results

We identified 33 metabolites differentially abundant between the lower and higher disease activity groups (P < 0.05). Additionally, we identified 51 metabolites associated with DAS28-CRP (P < 0.05). A GLM based upon these 51 metabolites resulted in higher prediction accuracy (mean absolute error [MAE] ± SD: 1.51 ± 1.77) compared to a GLM without feature selection (MAE ± SD: 2.02 ± 2.21). The predictive value of this feature set was further demonstrated on a validation cohort of twelve plasma samples, wherein we observed a stronger correlation between predicted and actual DAS28-CRP (with feature selection: Spearman’s ρ = 0.69, 95% CI: [0.18, 0.90]; without feature selection: Spearman’s ρ = 0.18, 95% CI: [−0.44, 0.68]). Lastly, among all identified metabolites, the abundances of eight were significantly associated with the CRP patient groups while controlling for potential confounders (P < 0.05).

Conclusions

We demonstrate for the first time the prediction of quantitative disease activity in RA using plasma metabolomes. The metabolites identified herein provide insight into circulating pro-/anti-inflammatory metabolic signatures that reflect disease activity and inflammatory status in RA patients.

Background

Rheumatoid arthritis (RA) is a chronic, autoimmune inflammatory disease primarily affecting the small diarthrodial joints and other organ systems [1,2,3,4] that can eventually lead to bone/cartilage erosion, joint deformity, loss in mobility, and organ damage [5]. Known to be associated with a variety of factors, such as genetic susceptibility [6], age [7], sex [8], smoking status [9], and dietary habits [10], RA is diagnosed in nearly 5 per 1000 adults worldwide, and women are 2 to 3 times more likely to develop RA than men [5]. In addition, as is the case with many complex and progressive disorders, patients with RA exhibit vast heterogeneity in clinical symptoms (e.g., joint inflammation, swelling, pain, stiffness) [11] and in responses to methotrexate and other disease-modifying anti-rheumatic drugs (DMARDs) [12]. Furthermore, immune cells (mainly, B cells, T cells, and macrophages) and cytokines are known to be implicated in RA pathogenesis [13]. For example, Haringman et al. observed that the abundance of macrophages in synovial tissue was positively correlated with disease activity [14], and Chung et al. identified significant differences in the levels of multiple cytokines (e.g., IL-6, IL-11, LIF) between RA and healthy controls [15]. In this regard, further understanding of the pathophysiological mechanisms that drive either progression or remission in RA disease activity would be important for identifying prognostic factors and developing more effective treatments [5, 16].

Having practical measures of disease activity is essential for determining the course of RA treatment and for monitoring patient response [3]. To this end, several studies have suggested strategies to quantify (or categorize) RA disease activity by using clinical and inflammatory core components, which include, but are not limited to, the number of tender and swollen joints, erythrocyte sedimentation rate (ESR), serum C-reactive protein (CRP) levels, and patients’ pain levels [17,18,19,20]. Among these various strategies, the modified Disease Activity Score that considers 28 joints (DAS28) with either ESR (DAS28-ESR) or CRP (DAS28-CRP) is currently one of the most well-recognized and recommended measures in RA [20].

An emerging area of RA research is in using high-throughput metabolomic profiling approaches, which comprehensively measure all small-molecule biochemicals in a biological specimen (e.g., plasma, serum, urine, synovial fluid) to enable biomarker discovery and novel insights into the biochemical processes governing disease pathophysiology [11, 21,22,23]. In particular, recent studies have demonstrated the promise of using such metabolomic technologies on patient-derived biospecimens for classifying patients with RA according to their disease activity categories [21, 24, 25], and for identifying metabolic signatures predictive of treatment response [26,27,28,29]. For instance, Teitsma et al. used metabolomic profiling in serum samples from early RA patients to identify metabolites and metabolic pathways that were significantly associated with sustained, drug-free remission (DAS28 < 2.6) after tocilizumab- or methotrexate-based therapy [24]. Likewise, Sasaki et al. identified 15 and 20 metabolites in plasma and urine, respectively, that were differentially abundant between active RA (DAS28-ESR ≥ 3.2) and inactive RA (DAS28-ESR < 3.2) [25]. These findings suggest that a wider application of global metabolomic profiling—coupled with advanced analytics [30]—can lead to the discovery of novel and predictive biomarkers that complement current standard laboratory tests for assessing disease activity in RA.

To date, a global metabolomic profiling analysis to demonstrate the predictive value of blood biochemicals in estimating disease activity scores for patients with RA has remained elusive. In this study, on 128 plasma metabolomic profiles from 64 RA patients, we utilize a multi-approach analysis to uncover metabolites that reflect and predict RA disease activity. First, we identify metabolites that stratify patients of “higher” (DAS28-CRP ≥ 3.2) and “lower” (DAS28-CRP < 3.2) disease activity groups. Next, we pinpoint specific metabolites that significantly associate with DAS28-CRP. Interestingly, a few of the metabolites identified through these two approaches were able to differentiate between the two groups of patients divided according to their C-reactive protein (CRP) levels in blood (“high-CRP,” CRP > 3.0 mg/L; “low-CRP,” CRP ≤ 3.0 mg/L); these metabolites may possibly reflect metabolic perturbations affected by worsening inflammatory activity. Finally, we utilize a machine learning technique to predict DAS28-CRP with plasma metabolite abundances. Importantly, we find that the feature selection step led to improved performance in predicting quantitative disease activity, and this translated reasonably well to a validation cohort. Taken together, our findings described herein support a key role for high-throughput metabolomic technologies in identifying blood-borne metabolic signatures of RA disease activity and lay the groundwork for monitoring disease progression and systemic inflammation using blood samples alone.

Materials and methods

Study population, subject enrollment, sample collection, and demographic characteristics

The study population consisted of consecutive patients with RA attending the outpatient practice of the Division of Rheumatology at Mayo Clinic in Rochester, MN. Eligibility required patients to be adults 18 years of age or older with a clinical diagnosis of RA by a rheumatologist, fulfilling the American College of Rheumatology/European League Against Rheumatism 2010 revised classification criteria for RA [2]. Patients were excluded if they did not comprehend English, were unable to provide written informed consent, or were members of a vulnerable population (e.g., incarcerated subjects). This led to a total of 76 patients fulfilling the eligibility criteria, who were partitioned into two groups (Table 1): for the discovery cohort of this study, 64 patients with available blood samples from at least two outpatient visits 6–12 months apart were included (128 total samples); for the validation cohort, 12 patients whose blood samples were available from only a single outpatient visit were included (12 total samples). Demographic and clinical data, including the numbers of tender and swollen joints, patient and evaluator global assessments, CRP (mg/L), body mass index (BMI, kg/m2), smoking status, and results for rheumatoid factor (RF, IU/mL) and anti-cyclic citrullinated peptide antibodies (anti-CCP), were collected from the electronic medical records. The patient samples (140 in total) in the study had established disease with a mean age of 63.54 (range 32–86), and 69.7% (53 of 76) were female. Disease activity varied from remission to high disease activity, with a DAS28-CRP mean of 3.0 (range 1.2–7.0), see Additional file 1 for the distribution of DAS28-CRPs corresponding to all study participants.

Table 1 Demographic characteristics of study participants

Metabolomic profiling

Untargeted metabolomic profiling of plasma samples from both discovery and validation cohorts through ultra-high-performance liquid chromatography-tandem mass spectrometry (UPLC-MS/MS) was performed by Metabolon Inc.’s (Durham, NC, USA) Discovery HD4™ platform. Detailed descriptions of all methods regarding metabolomic profiling are available in Additional file 4.

Analysis Workflow

Figure 1 provides a summary of the analytic strategy used on the 128 plasma samples of the discovery cohort to identify the associations between metabolites and RA disease activity. The analysis workflow consists of two complementary approaches: using mixed-effects logistic regression, the first approach identifies metabolites that are differentially abundant between the higher and lower disease activity groups, which were determined by DAS28-CRP scores [18,19,20, 31] (Fig. 1A); the second approach uses mixed-effects linear regression to model the relationship between DAS28-CRP and metabolite abundances, allowing the detection of key biochemical features that associate with quantitative disease activity (Fig. 1B). To test the predictive accuracy of these selected features when incorporated into a generalized linear model, an additional cohort of twelve plasma metabolomic profiles (from twelve RA patients obtained at single time points) was collected as an independent validation set.

Fig. 1
figure 1

A multi-approach discovery strategy to identify metabolites indicative of RA disease activity. (A) Differentially abundant metabolites between the higher and lower disease activity groups were identified using a mixed-effects logistic regression model adjusted for patient age and sex, as well as for patient ID to control for having multiple samples from the same patient. (B) A selection scheme to identify metabolites associated with DAS28-CRP. Metabolites were selected with mixed-effects linear regression. To further demonstrate their association with DAS28-CRP, these metabolites were used to construct a generalized linear model for predicting DAS28-CRP. The predictive performance of the model was evaluated on the discovery cohort (using a cross-validation technique) and on a validation cohort

Pre-processing of metabolomic profiling data

Statistical analyses on untargeted metabolomic data were performed using scaled imputed data provided by Metabolon, Inc. Briefly, the raw data were normalized to account for inter-day variation, which is a result of UPLC-MS/MS runs over multiple days; then, the peak intensities were rescaled to set each metabolite’s median equal to 1. Missing values were then imputed with the minimum observed value of the metabolite across all samples, finally yielding the scaled imputed data. In addition, metabolites with missing values in over 20% of the entire samples were removed, resulting in 686 metabolites remaining for further analysis. R (v3.6.1), lme4 package (v1.1.21) [32], Python3 (v3.7.5), and sklearn (v0.22.2) were used to perform all data pre-processing and statistical analyses.

Delineation of RA disease activity groups

Following previous reports [18,19,20, 31], samples from RA patients were divided into two disease activity groups based upon DAS28-CRP: “lower” (DAS28-CRP ≤ 3.2, n = 76) and “higher” (DAS28-CRP > 3.2, n = 52). These pre-defined two disease activity groups were used as the nominal response variable in a mixed-effects logistic regression model to identify differentially abundant metabolites between the two groups. The demographic characteristics of samples (n = 128) divided into lower and higher disease activity are summarized in Additional file 5.

Identification of differentially abundant metabolites while controlling for confounding factors

The following patient characteristics were examined to identify the potential confounding factors in the association between plasma metabolites and disease activity (i.e., higher or lower disease activity): age, sex, BMI, smoking history, and treatment use (for methotrexate, prednisone, non-methotrexate csDMARDs, TNFi-bDMARDs, and non-TNFi-bDMARDS). Based upon Fisher’s exact test, patient age (age ≤ 60, age > 60) and sex (male, female) were observed to have statistically significant associations with the two disease activity groups; the P-value for age and sex was P = 0.01 (odds ratio [OR] = 2.74, 95% confidence interval [CI] = 1.15–6.73) and P = 0.02 (OR = 0.37, 95% CI = 0.14–0.88), respectively. On the other hand, no statistically significant associations were observed between these two disease activity groups and BMI (BMI ≤ 30, BMI > 30; P = 0.32), disease duration (duration ≤ 9 years, duration > 9 years; P = 0.14), smoking history (smoked at least once, never smoked; P = 0.36), or treatment use (user, non-user) for methotrexate (P = 0.83); prednisone (P = 0.58); TNFi-bDMARDs, i.e., adalimumab, certolizumab, etanercept, and infliximab (P = 0.18); non-TNFi-bDMARDs, i.e., abatacept, rituximab, and tocilizumab (P = 0.76); or other non-methotrexate csDMARDs, i.e., azathioprine, hydroxychloroquine, leflunomide, and sulfasalazine (P = 0.71). In addition, no significant changes in treatment use were observed between the two visits; P-values of the associations between treatment use and time point based upon McNemar’s chi-squared test for paired nominal data were as follows: methotrexate (P = 1), prednisone (P = 1), TNFi-bDMARDs (P = 0.75), non-TNFi-bDMARDs (P = 1), and non-methotrexate csDMARDs (P = 0.07). Therefore, the mixed-effects logistic regression model was adjusted for age and sex as fixed effects, but not for all other aforementioned covariates. In accordance with these results, age and sex have been previously reported to be connected to RA disease activity [33,34,35]. Herein, patient ID was considered as a random effect in the model to account for intra-subject variance due to having repeated measurements from a single patient. By controlling for patient ID (which are unique to each patient) as a random effect, we are acknowledging the non-independence in our data, that is, sampling that has taken place from within a patient. Leveraging multiple samples from the same patient allows us to compensate for the small number of samples in higher disease activity (DAS28-CRP > 3.2) in each visit (visit 1 and visit 2 having 25 and 27 samples, respectively) by maximizing the degree of freedom for the quantitative disease activity measure, and thereby to boost statistical power. Importantly, no significant difference was observed in DAS28-CRP between visit 1 and visit 2 (P = 0.98, Wilcoxon signed-rank test). Metabolites whose corresponding coefficients of the regression model were of P-value < 0.05 were considered as differentially abundant, that is, having a statistically significant association with the disease activity group.

Selection of metabolites associated with DAS28-CRP

Selection of metabolites associated with DAS28-CRP was performed with a mixed-effects linear regression model (DAS28-CRP as the continuous response variable), which controls for fixed effects (scaled metabolite abundances, patients’ age and sex) and for random effects (patient ID). Satterthwaite’s degrees of freedom method supported by lmerTest (v3.1.1) [36] was applied to test for the statistical significance (P-value) of associations between metabolites and DAS28-CRP. P-values were retrieved from the corresponding regression coefficients of the predictor variables.

Evaluation of predictive performance of DAS28-CRP-associated metabolites

A generalized linear model (GLM) was used to estimate DAS28-CRP scores using the aforementioned significantly associated metabolites as predictor variables. The predictive performance of the parameterized model was evaluated by two different techniques: First, a modified leave-one-out cross-validation approach was applied to the 128 samples of the training group (discovery cohort). More specifically, in each cross-validation loop, both samples from the same patient were allocated as an internal validation set, while all remaining samples (126 samples from 63 patients) were used to select metabolites significantly associated with DAS28-CRP (P < 0.05). These selected biochemical features were then included in a GLM for predicting DAS28-CRP scores of the remaining two samples (of the internal validation group) from their metabolite abundances. The second approach considers testing a GLM, which was composed of the DAS28-CRP-associated metabolites identified from all 128 samples of the training group, on the independent validation group of 12 plasma samples (validation cohort). For both techniques, model performance was reported using mean absolute error (MAE) and standard deviation (SD).

Identification of metabolites associated with treatment use

A marginal, mixed-effects linear regression model was used to relate metabolite abundance with treatment use. Scaled metabolite abundance, treatment use, and patient ID were set as the response variable, predictor variable (fixed effect), and random effect, respectively. The use of the following treatments was assessed individually: methotrexate, prednisone, non-methotrexate csDMARDs, TNFi-bDMARDs, and non-TNFi-bDMARDs (names of individual drugs in each treatment group are provided in the footnote of Table 1). P-values were retrieved from the corresponding regression coefficient of the predictor variable (i.e., use or non-use), and a significance of P < 0.05 was reported as statistically significant.

Identification of differentially abundant metabolites between the two CRP groups

Metabolites that are significantly associated with disease activity groups and DAS28-CRP scores were further investigated to find those associated with patient groups delineated by CRP levels. First, all samples were divided into two groups as follows: “high-CRP” (CRP > 3.0 mg/L, n = 52) and “low-CRP” (CRP ≤ 3.0 mg/L, n = 76). Next, a marginal, mixed-effects linear regression model was used to define the abundance of a metabolite based upon the following fixed effects: CRP group, sex, age, smoking history, and treatment with prednisone, methotrexate, non-methotrexate csDMARDs, TNFi-bDMARDs, or non-TNFi-bDMARDs. Additionally, patient ID was treated as a random effect. Any covariates whose association with metabolite abundance was statistically significant (i.e., P-value of the corresponding regression coefficient < 0.05) were then included in an adjusted mixed model for metabolite abundance. Finally, metabolites were considered as differentially abundant between the two CRP groups if the association between metabolite abundance and CRP group was still found to be significant in the adjusted model (P < 0.05).

Results

Differentially abundant metabolites between higher and lower disease activity groups

As shown in our analysis workflow (Fig. 1), we first sought metabolites that were significantly different in abundance between two major disease activity groups. For this, we divided the 128 metabolomic profiles into two major categories (“higher” vs. “lower”) based upon the reported disease activity of the corresponding patient at the time of sample collection (“Materials and methods” section). Using a mixed-effects logistic regression model (“Materials and methods” section), we identified 33 metabolites as differentially abundant between the higher (n = 52) and lower (n = 76) DAS28-CRP groups (Fig. 2). Most of these metabolites (31 of 33) were observed to have significantly increased abundances in lower disease activity, whereas the remaining two (glucuronate and hypoxanthine) were found to be significantly increased in higher disease activity. Notably, of the 31 metabolites increased in lower disease activity, seven metabolites (3-hydroxydecanoylcarnitine, dihomo-linoleoylcarnitine (C20:2), eicosenoylcarnitine (C20:1), linoleoylcarnitine (C18:3), linoleoylcarnitine (C18:2), stearoylcarnitine (C18), palmitoylcarnitine (C16)) are a part of acylcarnitine metabolism and represent a 3.6-fold enrichment in metabolites involved in this particular pathway (P = 1.9 × 10−3, hypergeometric test). It is important to note that the differences seen are relatively small in terms of fold change, with most of the metabolites varying by 1.1–1.3-fold. Despite these subtle differences within RA patients of varying disease activities, we were still able to obtain statistically significant signal even after considering and controlling for all known potentially confounding factors (which often leads to a reduction in statistical power), while adhering to our cutoffs for statistical significance (P < 0.05).

Fig. 2
figure 2

Plasma metabolites differentiating between higher and lower disease activity groups in RA. A total of 2 and 31 metabolites were found to be significantly increased in higher (DAS28-CRP > 3.2, n = 52) and lower (DAS28-CRP ≤ 3.2, n = 76) disease activity groups, respectively. Each point corresponds to a metabolite (686 total). Differentially abundant metabolites were found using a mixed-effects logistic regression model on the discovery cohort (128 samples), for which age and sex were adjusted. Metabolites with a P-value < 0.05 (based upon the corresponding coefficient of the regression model) were considered as significantly different between the groups. P-values and fold changes for all metabolites are listed in Additional file 6. Metabolites in bold have been previously described in the literature for their associations with RA

N2-acetyl,N6-methyllysine (|log2(FC)| = 1.11, P = 1.26 × 10−2) and trigonelline (N′-methylnicotinate) (|log2(FC)| = 0.74, P = 2.09 × 10−2), which were both found to have increased abundance in lower disease activity, were the top two metabolites having the largest fold changes between the two groups. Although the direct relevance of N2-acetyl,N6-methyllysine to RA is currently not well understood, N2-acetyl,N6-methyllysine is part of the lysine metabolism pathway, which has been reported to be associated with RA in the following studies: (i) according to Teitsma et al., serum metabolites associated with lysine degradation were observed to have a higher concentration in early RA patients who achieved sustained drug-free remission (after tocilizumab- or methotrexate-based treatment) compared to those who never achieved a drug-free status [24]; (ii) Yang et al. reported that metabolic products of lysine degradation (carnitine and pipecolic acid) were significantly increased in RA patients than in normal subjects [37].

In regard to trigonelline, which is a product of niacin (vitamin B3) metabolism, this alkaloid has been suggested to have therapeutic potential for diabetes and central nervous system disease [38], and also reported to demonstrate anti-inflammatory properties in mice [39]. In accordance with our results showing a decreased abundance of trigonelline in higher disease activity, trigonelline could be of interest in future studies on inflammatory responses in RA.

Biliverdin (|log2(FC)| = 0.48, P = 1.38 × 10−2) and bilirubin (E,E) (|log2(FC)| = 0.43, P = 1.18 × 10−2), which are known metabolic products of the heme catabolic pathway, were also observed to have significantly increased abundances in lower disease activity. In particular, biliverdin has been shown to (i) inhibit the activation of pro-inflammatory transcription factors, including NFkB both in vitro and in vivo [40,41,42,43,44]; (ii) inhibit the proliferation of primary T cells stimulated with anti-CD3 and anti-CD28 monoclonal antibodies by inhibiting NFAT/NF-kB activation in a mouse model of heart transplantation [45]; and (iii) improve corneal inflammation mediated by heme-oxygenase 2 (HO-2) deficiency in a transgenic mouse model [41]. Moreover, bilirubin, which is derived from the reduction of biliverdin by biliverdin reductase, has been reported as a potential biomarker for RA in line with our findings. For example, Peng et al. observed a decreased concentration of serum bilirubin in RA patients compared to healthy controls, as well as in RA patients with worsening disease activity [46]. Additionally, Fischman et al. found that total bilirubin levels are inversely related to RA disease activity even after adjusting for multiple confounders (e.g., age, sex, race) and discussed the possibility of bilirubin (a known anti-oxidant) having a physiological anti-inflammatory effect [47]. This point is further elaborated upon by Jangi et al. [48], who have described in detail the immunosuppressive properties of unconjugated bilirubin in RA and other inflammatory disorders. The full list of differentially abundant metabolites and their associated pathways is shown in Additional file 6.

Metabolic feature selection improves DAS28-CRP prediction accuracy

Having uncovered metabolites demonstrating altered abundance between two major disease activity groups, we next asked whether quantitative disease activity can be predicted with plasma metabolomes. As untargeted metabolomic profiling can yield a considerable amount of noise and random fluctuations in observed signals [49], it is necessary to first select informative metabolic features that reliably capture relevant aspects of the phenotype of interest [50]. For this, we used mixed-effects linear regression models to select metabolites significantly associated with DAS28-CRP. Afterwards, the abundances of the selected metabolic features were incorporated into a generalized linear model (GLM) to predict DAS28-CRP. For comparison purposes, a GLM was constructed without metabolic feature selection and thereby taking into consideration all features of a metabolomic profile. Details regarding GLM construction and performance evaluation are provided in the “Materials and methods” section.

When applying a modified leave-one-out cross-validation technique to the training group samples (n = 128), we found that the GLM incorporating metabolites that were significantly associated with DAS28-CRP outperformed the model without feature selection (i.e., using all metabolites). As shown in Fig. 3, the distribution of absolute errors between the observed and predicted DAS28-CRP scores was smaller (with respect to the cumulative area under the error curve) for the GLM with feature selection than that without feature selection. To this point, the prediction MAE (± SD) of the GLM with and without feature selection was 1.51 (±1.89) and 2.02 (±2.52), respectively.

Fig. 3
figure 3

Evaluation of DAS28-CRP predictive performance in cross-validation. A modified leave-one-out cross-validation approach was used on the samples of the training group (128 samples) to test the performance of a generalized linear model (GLM) in predicting DAS28-CRP scores from metabolite abundances. Distributions of absolute errors from models with and without a feature selection scheme were compared to identify the more robust model. The GLM with the feature selection scheme performed better (MAE ± SD: 1.51 ± 1.89) than the model without feature selection (MAE ± SD: 2.02 ± 2.52)

Having confirmed that feature selection can lead to a more accurate prediction model in cross-validation, we applied the same scheme to all metabolome samples of the discovery cohort to obtain a final set of metabolites associated with DAS28-CRP (P < 0.05). After adjusting for potential confounding factors (“Materials and methods” section), this resulted in a collection of 51 plasma metabolites (Table 2). These metabolites were used to construct a final GLM, whose predictive accuracy was tested on an independent validation cohort (n = 12) of plasma metabolomic profiles from twelve RA patients (importantly, this additional cohort was not drawn from the same population distribution from which the features were derived). On this previously unseen cohort, the GLM constructed with only the 51 selected metabolites performed considerably better than the model without the feature selection scheme by over twofold (Fig. 4A); the prediction MAE of the GLM with and without feature selection was 0.97 (±0.47) and 2.01 (±2.18), respectively. Likewise, when the actual and predicted DAS28-CRPs were plotted together for both GLMs (Fig. 4B), we found that the model with the selection scheme performed more favorably. More specifically, a stronger correlation between the actual and predicted disease activity scores was observed in the model with feature selection (Spearman’s ρ = 0.69, P = 1.40 × 10−2, 95% CI: [0.18, 0.90]) compared to the model without (Spearman’s ρ = 0.18, P = 5.72 × 10−2 , 95% CI: [−0.44, 0.68]).

Table 2 Plasma metabolites significantly associated with DAS28-CRP
Fig. 4
figure 4

GLM with feature selection provides improved DAS28-CRP prediction accuracy in an independent validation group (12 samples). A Performance of GLMs in predicting quantitative disease activity was evaluated on samples of an independent validation group. Distributions of absolute errors from models with and without a feature selection scheme were compared to identify the more robust model. B Selection of metabolic features prior to model training resulted in higher predictive performance, as evidenced by the stronger correlation between observed and predicted DAS28-CRPs. Three samples predicted to have negative DAS28-CRP values are omitted from the scatter plot. The dashed violet line indicates “y = x,” i.e., an exact match between the observed and predicted values. 95% confidence interval for ρ with feature selection [0.18, 0.90]; without feature selection [−0.44, 0.68]

Commonly identified metabolites from two different analytic approaches

To summarize the findings above, we found that, from the 686 total detectable metabolites in a metabolomic profile, 33 (4.8%) were differentially abundant between higher and lower disease activity, and 51 (7.4%) were significantly associated with DAS28-CRP (Fig. 5). These separate findings amounted to a total of 67 unique metabolites, among which 40 were found to have no association with the use of prednisone, methotrexate, other non-methotrexate csDMARDs, TNFi-bDMARDs, or non-TNFi-bDMARDs (“Materials and methods” section). Notably, eight metabolites (6-bromotryptophan, bilirubin (E,E), biliverdin, glucuronate, N-acetyltryptophan, N-acetyltyrosine, serine, and trigonelline) were not only consistently detected across both analytic approaches, but also found to have no association with any treatment use; these results strongly suggest key metabolic pathways and modules potentially contributing to, or serving as indicators of, RA pathogenesis independent of confounding treatment effects. Consistent with this idea, additional studies into the metabolites found in this study (the majority of which have yet to be linked to RA) may be able to provide new insight into the perturbed physiological metabolic processes—which are then in turn reflected in blood—underlying disease progression in RA.

Fig. 5
figure 5

Venn diagram of all plasma metabolites identified through the multi-approach discovery strategy. A total of 67 unique metabolites were identified, among which 40 were found to have no association with the use of treatment. Notably, eight metabolites (6-bromotryptophan, bilirubin (E,E), biliverdin, glucuronate, N-acetyltryptophan, N-acetyltyrosine, serine, and trigonelline) in bold were not only consistently detected across both analytic approaches, but also found to have no association with any treatment use. Colored circles indicate metabolites whose abundances associate with treatment use. Metabolites with red triangles were found to have increasing abundances with worsening disease activity, whereas metabolites with blue triangles were found to have decreasing abundances with worsening disease activity

Metabolites associated with CRP patient groups

Elevated levels of C-reactive protein (CRP) in the blood are well known to often indicate increased inflammatory conditions, which may be caused by a wide variety of acute (e.g., infections) and chronic disorders (e.g., rheumatoid arthritis, inflammatory bowel disease). In RA patients, CRP levels have been observed to increase after acute mental stress tasks [51] and also to be linked to risk of cardiovascular disease [52]. Furthermore, several serum metabolites were found to reflect inflammatory activity in patients with early arthritis [53].

We further investigated the aforementioned 67 plasma metabolites to see whether any were differentially abundant between the two CRP patient groups, i.e., “high-CRP” (CRP > 3.0 mg/L, n = 52) and “low-CRP” (CRP ≤ 3.0 mg/L, n = 76) (“Materials and methods” section). While controlling for potential confounding variables, we identified eight total metabolites that were significantly associated with CRP patient group. More specifically, the abundances of mannose, beta-hydroxyisovalerate, (14 or 15)-methylpalmitate (a17:0 or i17:0), erucate (22:1n9), 10-undecenoate (11:1n1), N-acetylcitrulline were higher in high-CRP, while those of serine and linoleoylcarnitine (C18:3) were lower in high-CRP (Fig. 6). Application of these plasma metabolites, which were found to be connected to both RA disease activity and circulating CRP levels, may lead to the development of new clinical laboratory tests to further enable precision medicine for RA patients.

Fig. 6
figure 6

Metabolites differentially abundant between the two CRP patient groups. Among the 67 total metabolites identified through our multi-approach analysis on the discovery cohort (n = 128), eight metabolites were identified to have significant associations with the CRP group while controlling for confounding variables (regression coefficient for CRP, P < 0.05). A Metabolites with higher abundances in the high-CRP group: mannose, beta-hydroxyisovalerate, (14 or 15)-methylpalmitate (a17:0 or i17:0), erucate (22:1n9), 10-undecenoate (11:1n1), and N-acetylcitrulline. B Metabolites with higher abundances in the low-CRP group: serine and linoleoylcarnitine (C18:3)

Plasma metabolites associate with clinical improvement in RA

Based upon the European League Against Rheumatism (EULAR) response criteria for DAS28-CRP [54], we found that sixteen of the 64 patients in the discovery cohort showed moderate or good improvement in disease activity from visit 1 to visit 2, while the remaining 48 patients did not show clinical improvement at the time of their second visit. This discovery provided an entry point for the following analysis: For each of these two patient groups, i.e., “improved” (n = 16) and “non-improved” (n = 48) patients, we aimed to identify metabolites whose abundances significantly changed from visit 1 to visit 2, while controlling for the same confounding factors (mixed-effects regression model, P < 0.05). As a result, we identified eleven metabolites whose abundances significantly changed in the improved patient group (Additional file 7), while nineteen metabolites showed significant changes in the non-improved patient group (Additional file 8). The following three metabolites, which were discovered in our previous analyses on the 128 plasma metabolome samples of the discovery cohort, were detected once again: erucate (22:1n9), a metabolite identified to be associated with both the DAS28-CRP and CRP patient groups, was identified to be significantly different between visit 1 and visit 2 in patients who did not show clinical improvement (non-improved); 3-amino-2-piperidone, a metabolite identified to be differentially abundant between higher and lower disease activities in our study, was identified to be significantly different between visit 1 and visit 2 in patients who showed clinical improvement (improved); and gamma-glutamylmethionine, a metabolite identified to be differentially abundant between higher and lower disease activities, was identified to be significantly different between visit 1 and visit 2 in the non-improved group. These results allow us to expand our future direction to investigate metabolites associated with clinical improvement in patients with RA.

Discussion

Dysfunctions in cellular and tissue metabolism are involved in a broad range of autoimmune disorders [55,56,57,58], including RA [59,60,61]. These metabolic implications highlight the importance of investigating which biochemical functions and metabolic states are altered during the onset and progression of the disease. To this end, metabolomics platforms (and the accompanying wealth of data) can present unique opportunities for discovering novel disease “molecular signatures” [62], which can be interpreted through the lens of annotated biochemical relationships. Moreover, high-throughput profiling can enable the identification of circulating pro-inflammatory (disease-triggering) and anti-inflammatory (disease-protective) metabolites in RA, as elaborated upon by Coras et al. [63]; such discoveries may facilitate the design of either dietary or gut microbiome-based intervention strategies to improve wellness or alter the course of disease for RA patients. In this study, by performing a global metabolomic profiling analysis on 128 plasma samples obtained from 64 patients with RA, we identified biochemical signatures associated with, and predictive of, disease activity. Mainly, through a combination of statistical approaches for metabolic signature discovery, we identified several metabolites that (i) differ significantly between the lower and higher disease activity groups and (ii) significantly associate with DAS28-CRP. Of note, our study is the first to leverage biochemical features from a plasma metabolomic profile to predict quantitative disease activity.

Interestingly, we identified eight metabolites (6-bromotryptophan, bilirubin (E,E), biliverdin, glucuronate, N-acetyltryptophan, N-acetyltyrosine, serine, and trigonelline) that were commonly found across different statistical approaches, possibly capturing representative metabolite signals of RA progression. We discussed above the reported roles of bilirubin and biliverdin in RA. Moreover, these two metabolites were previously reported for their cytoprotective and anti-inflammatory effects [64,65,66,67,68] and even suggested as an “RA protective factor” by Fischman et al. [47]. Interestingly, high concentrations of bilirubin and biliverdin were reported in other inflammatory disorders, such as atherosclerotic diseases [69] and autoimmune encephalitis [67]. In regard to the other six metabolites, clear and definitive connections with RA have not yet been established. However, if our results on bilirubin and biliverdin were to serve as benchmarks for reliably identifying plasma metabolites important to RA disease activity, then these remaining metabolites may be deemed as leading candidates for future investigations.

Glucuronate was found to show elevated abundance in higher disease activity than in lower disease activity. This glucose derivative is involved in the detoxification of xenobiotics via glucuronidation in the mammalian liver. Interestingly, this process can be reversed by gut bacteria-harboring β-D-glucuronidases [70] and thereby releasing (potentially toxic) exogenous compounds into the gut lumen and subsequently into circulation [71, 72]. In that respect, examining a possible role involving dysbiosis in the gut microbiota—combined with metabolomic approaches to infer relationships between gut microbes and blood metabolites in RA, as shown by Chen et al. [73]—may help elucidate a microbial-based mechanism explaining the observed alterations in plasma glucuronate.

Serine was seen to decrease with worsening disease activity. In line with our results, albeit in an RA mouse model with collagen-induced arthritis (CIA), plasma levels of serine and other free amino acids were found to have significantly decreased in the CIA group compared to control mice [74]. In another study wherein synovial fluid of RA patients was examined for citrullinated proteins (which is widely known to result in a rise in anti-citrullinated protein antibodies in RA), Tilvawala et al. found increased citrullination in a wide array of serine protease inhibitors (serpins) and serine proteases [75]; in the same study, the investigators demonstrated in vitro that citrullinating serine protease inhibitors nearly abolishes their inhibitory activity towards their target proteases. Although we have yet to uncover whether a decrease in plasma serine levels (with worsening disease activity) is linked to the citrullination of serine proteases in synovial fluid, we speculate that changes in serine may reflect dysregulated protein degradation during systemic inflammatory activity and joint destruction in RA.

Of note, the metabolites reported in our study display considerable differences with those previously found to be associated with RA. For example, Chu et al. reported multiple lists of circulating plasma metabolites relevant to RA risk in the NHS (Nurses’ Health Study)/NHSII [76]. Among these lists, our study found only one (i.e., N-acetyltryptophan) out of seven metabolites associated with the incidence of RA, and one (i.e., 4-guanidinobutanoate) out of eight metabolites associated with the incidence of RA specifically within 5 years of diagnosis. Likewise, our study identified one metabolite (i.e., glucose) among eleven previously reported to differentiate early RA patients from healthy controls [53]. Such discrepancies could stem from differences in a range of factors between these other studies and ours, such as clinical and demographic characteristics, sample size, disease duration, and analytical approach (e.g., metabolomic instruments, statistical tests). The results of our pilot study herein motivate future investigations—composed of larger sample sizes and of a wider variety of cohort characteristics—that may eventually provide findings more consistently aligned with previous work by others.

We note a few limitations of this study: First, we acknowledge that our study includes a relatively small number of samples within each disease activity group of the discovery cohort and the validation cohort. Nevertheless, we were able to detect statistically significant metabolites in all analytical strategies, demonstrating that our data provides reasonably sufficient statistical power. Certainly, a much larger cohort would have been ideal; however, this is a small pilot study on our stored plasma samples, and obtaining an additional cohort is outside of its scope. Encouragingly, despite the low sample size of the validation cohort, we were able to successfully show that feature selection is a necessary step in the model-training process, and we expect this finding to translate well to larger cohorts in our future studies. Nevertheless, in order to more meticulously examine the role of blood metabolites in RA, future investigations will warrant a higher number of samples and more detailed subject characteristics. Second, to define RA disease activity, we solely used the DAS28-CRP scoring system, which is dependent upon acute-phase responses that may not accurately reflect patients who have an inflammation-free state [24]. Our future plans include performing our analytical pipeline with other RA disease activity metrics (e.g., clinical disease activity index (CDAI), simple disease activity index (SDAI)) to test the robustness of our findings. Third, all of our multivariate analyses followed adjustment for patient age and sex only. Other potential confounders that may affect the concentration of blood metabolites, such as diet, exercise habits, lifestyle factors, time of the day of sample collection, and gut microbiome, were not considered as predictor variables in our analyses. Fourth, comorbidity can certainly be a significant confounding factor. Alternatively, comorbidities could theoretically be contributors or mediators to inflammatory disease activity in patients with RA. At this stage, it would have been premature to adjust statistically for the effects of particular comorbidities or for the presence of multi-morbidity (i.e., multiple chronic conditions) before carefully investigating their potential interaction with plasma metabolites and RA disease activity. This may lead to “over-adjustment” and falsely concluding that certain metabolites are not significant when they may, in fact, be very important. Future studies will be necessary to explore the potential interactions between comorbidities, RA disease activity, and plasma metabolite levels in RA. Last, despite the similarities in our findings with previous investigations (as noted above), many of our results are reported for the first time and remain to be validated by others. Possible causes of discrepancies with the work of others include the comparatively small number of samples in this pilot study, technical and biological sources of random noise, the uniqueness of our recruited patient cohort, variabilities in detection protocols and instrument sensitivity, and the use of alternate statistical techniques and potential over-fitting. Future efforts, by us and others, are likely to elucidate truly robust signals and collectively strengthen the confidence in our novel findings.

Despite the aforementioned limitations, our study establishes the far-reaching utility of using cutting-edge technological and analytical approaches for plasma metabolomic profiling and justifies analogous investigations at larger scales. The identified metabolites could be a reflection of the perturbed metabolic processes concurrent with worsening disease activity, and our findings will inspire future studies into how inflammation and pain in RA are coupled to physiological metabolism. Moreover, our identified sets of signature metabolites offer a promising glimpse into the biomolecular marker panels for diagnosing disease activity of RA patients solely through blood (thereby complementing current diagnostic approaches), with the overall aim to make such assessments faster, cheaper, and less invasive. In turn, studies such as ours are expected to contribute towards fully realizing the potential of virtual and digital healthcare by foregoing the need for patients to physically arrive at the clinic to meet their primary care provider in person.

As the gut microbiome has been recognized to be implicated in RA [73, 77,78,79]—possibly through complex mechanisms underlying microbe-microbe and host-microbe biochemical cross-talk [80]—integrated profiling of both stool metagenome and blood metabolome would provide an in-depth, comprehensive view of functional dysbiosis during RA onset and progression. Interestingly, a recent study showed that blood metabolites can be predictive of gut microbiome alpha-diversity [81]. Such investigations into integrating across multiple data types of the same phenotype can help to amplify the primary biological signal of interest relative to noise, as well as provide actionable insights. In conclusion, the results reported herein are poised to eventually improve disease management and outcomes of patients with RA and other rheumatic diseases, as well as to provide novel means of monitoring health and wellness [82].

Conclusions

This is the first study that uses global metabolomic profiling on blood to demonstrate the predictive value of circulating biochemicals in estimating quantitative disease activity in RA. We uncovered RA plasma metabolites that stratify patients of the “higher” and “lower” disease activity groups, significantly associate with DAS28-CRP, and differentiate between the two groups of patients divided by their blood C-reactive protein (CRP) levels. In all, our study supports the key role of high-throughput metabolomic technologies in identifying blood-borne biochemical signatures and metabolic pathways reflective of RA progression and systemic inflammation.

Availability of data and materials

All raw metabolomic datasets, as well as source codes used to reproduce all results in this study, are available at https://github.com/jaeyunsung/RA_plasma_metabolomics_2021.

Abbreviations

RA:

Rheumatoid arthritis

DMARDs:

Disease-modifying anti-rheumatic drugs

bDMARDs:

Biologic disease-modifying anti-rheumatic drugs

csDMARDs:

Conventional synthetic disease-modifying anti-rheumatic drugs

UPLC-MS/MS:

Ultra-high-performance liquid chromatography-tandem mass spectrometry

DAS28-CRP:

Disease Activity Score-28 using C-reactive protein

BMI:

Body mass index

GLM:

Generalized linear model

MAE:

Mean absolute error

SD:

Standard deviation

ESR:

Erythrocyte sedimentation rate

CRP:

C-reactive protein

Anti-CCP:

Anti-cyclic citrullinated peptide antibodies

HMDB:

Human metabolome database

CIA:

Collagen-induced arthritis

NHS:

Nurses’ Health Study

CDAI:

Clinical disease activity index

SDAI:

Simple disease activity index

References

  1. Guo Q, Wang Y, Xu D, Nossent J, Pavlos NJ, Xu J. Rheumatoid arthritis: pathological mechanisms and modern pharmacologic therapies. Bone Res. 2018;6(1):1–14.

    Article  Google Scholar 

  2. Aletaha D, Neogi T, Silman AJ, Funovits J, Felson DT, Bingham CO III, et al. 2010 Rheumatoid Arthritis Classification Criteria: an American College of Rheumatology/European League Against Rheumatism collaborative initiative. Arthritis Rheum. 2010;62(9):2569–81. https://doi.org/10.1002/art.27584.

    Article  PubMed  Google Scholar 

  3. Smolen JS, Aletaha D, Barton A, Burmester GR, Emery P, Firestein GS, et al. Rheumatoid arthritis. Nat Rev Dis Primers. 2018;4(1):18001. https://doi.org/10.1038/nrdp.2018.1.

    Article  PubMed  Google Scholar 

  4. Firestein GS. Evolving concepts of rheumatoid arthritis. Nature. 2003;423(6937):356–61. https://doi.org/10.1038/nature01661.

    Article  CAS  PubMed  Google Scholar 

  5. Aletaha D, Smolen JS. Diagnosis and management of rheumatoid arthritis: a review. JAMA. 2018;320(13):1360–72. https://doi.org/10.1001/jama.2018.13103.

    Article  Google Scholar 

  6. MacGregor AJ, Snieder H, Rigby AS, Koskenvuo M, Kaprio J, Aho K, et al. Characterizing the quantitative genetic contribution to rheumatoid arthritis using data from twins. Arthritis Rheum. 2000;43(1):30–7. https://doi.org/10.1002/1529-0131(200001)43:1<30::AID-ANR5>3.0.CO;2-B.

    Article  CAS  PubMed  Google Scholar 

  7. Crowson CS, Matteson EL, Myasoedova E, Michet CJ, Ernste FC, Warrington KJ, et al. The lifetime risk of adult-onset rheumatoid arthritis and other inflammatory autoimmune rheumatic diseases. Arthritis Rheum. 2011;63(3):633–9. https://doi.org/10.1002/art.30155.

    Article  PubMed  PubMed Central  Google Scholar 

  8. Ngo ST, Steyn FJ, McCombe PA. Gender differences in autoimmune disease. Front Neuroendocrinol. 2014;35(3):347–69. https://doi.org/10.1016/j.yfrne.2014.04.004.

    Article  CAS  PubMed  Google Scholar 

  9. Sugiyama D, Nishimura K, Tamaki K, Tsuji G, Nakazawa T, Morinobu A, et al. Impact of smoking as a risk factor for developing rheumatoid arthritis: a meta-analysis of observational studies. Ann Rheum Dis. 2010;69(01):70–81. https://doi.org/10.1136/ard.2008.096487.

    Article  CAS  PubMed  Google Scholar 

  10. Khanna S, Jaiswal KS, Gupta B. Managing rheumatoid arthritis with dietary interventions. Fron Nutr. 2017;4:52. https://doi.org/10.3389/fnut.2017.00052.

    Article  CAS  Google Scholar 

  11. Julià A, Alonso A, Marsal S. Metabolomics in rheumatic diseases. Int J Clin Rheumatol. 2014;9(4):353–69. https://doi.org/10.2217/ijr.14.25.

    Article  CAS  Google Scholar 

  12. Aletaha D, Smolen J. The rheumatoid arthritis patient in the clinic: comparing more than 1300 consecutive DMARD courses. Rheumatology. 2002;41(12):1367–74. https://doi.org/10.1093/rheumatology/41.12.1367.

    Article  CAS  PubMed  Google Scholar 

  13. Yap H-Y, Tee SZ-Y, Wong MM-T, Chow S-K, Peh S-C, Teow S-Y. Pathogenic role of immune cells in rheumatoid arthritis: implications in clinical treatment and biomarker development. Cells. 2018;7(10):161. https://doi.org/10.3390/cells7100161.

    Article  CAS  PubMed Central  Google Scholar 

  14. Haringman JJ, Gerlag DM, Zwinderman AH, Smeets TJ, Kraan MC, Baeten D, et al. Synovial tissue macrophages: a sensitive biomarker for response to treatment in patients with rheumatoid arthritis. Ann Rheum Dis. 2005;64(6):834–8. https://doi.org/10.1136/ard.2004.029751.

    Article  CAS  PubMed  Google Scholar 

  15. Chung S-J, Kwon Y-J, Park M-C, Park Y-B, Lee S-K. The correlation between increased serum concentrations of interleukin-6 family cytokines and disease activity in rheumatoid arthritis patients. Yonsei Med J. 2011;52(1):113–20. https://doi.org/10.3349/ymj.2011.52.1.113.

    Article  CAS  PubMed  Google Scholar 

  16. Arts EE, Fransen J, Den Broeder AA, van Riel PL, Popa CD. Low disease activity (DAS28≤3.2) reduces the risk of first cardiovascular event in rheumatoid arthritis: a time-dependent Cox regression analysis in a large cohort study. Ann Rheum Dis. 2017;76(10):1693–9. https://doi.org/10.1136/annrheumdis-2016-210997.

    Article  CAS  PubMed  Google Scholar 

  17. Van Gestel A, Prevoo M, Van’t Hof M, Van Rijswijk M, Van de Putte L, Van Riel P. Development and validation of the European League Against Rheumatism response criteria for rheumatoid arthritis: comparison with the preliminary American College of Rheumatology and the World Health Organization/International League Against Rheumatism criteria. Arthritis Rheum. 1996;39(1):34–40.

    Article  Google Scholar 

  18. Prevoo M, Van’T Hof MA, Kuper H, Van Leeuwen M, Van De Putte L, Van Riel P. Modified disease activity scores that include twenty-eight-joint counts development and validation in a prospective longitudinal study of patients with rheumatoid arthritis. Arthritis Rheum. 1995;38(1):44–8. https://doi.org/10.1002/art.1780380107.

    Article  CAS  PubMed  Google Scholar 

  19. Fransen J, Van Riel P. The Disease Activity Score and the EULAR response criteria. Clin Exp Rheumatol. 2005;23(5):S93.

    CAS  PubMed  Google Scholar 

  20. Anderson J, Caplan L, Yazdany J, Robbins ML, Neogi T, Michaud K, et al. Rheumatoid arthritis disease activity measures: American College of Rheumatology recommendations for use in clinical practice. Arthritis Care Res. 2012;64(5):640–7. https://doi.org/10.1002/acr.21649.

    Article  Google Scholar 

  21. Li J, Che N, Xu L, Zhang Q, Wang Q, Tan W, et al. LC-MS-based serum metabolomics reveals a distinctive signature in patients with rheumatoid arthritis. Clin Rheumatol. 2018;37(6):1493–502. https://doi.org/10.1007/s10067-018-4021-6.

    Article  PubMed  Google Scholar 

  22. Patti GJ, Yanes O, Siuzdak G. Metabolomics: the apogee of the omics trilogy. Nat Rev Mol Cell Biol. 2012;13(4):263–9. https://doi.org/10.1038/nrm3314.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Clish CB. Metabolomics: an emerging but powerful tool for precision medicine. Mol Case Stud. 2015;1(1):a000588. https://doi.org/10.1101/mcs.a000588.

    Article  Google Scholar 

  24. Teitsma XM, Yang W, Jacobs JW, Pethö-Schramm A, Borm ME, Harms AC, et al. Baseline metabolic profiles of early rheumatoid arthritis patients achieving sustained drug-free remission after initiating treat-to-target tocilizumab, methotrexate, or the combination: insights from systems biology. Arthritis Res Ther. 2018;20(1):230. https://doi.org/10.1186/s13075-018-1729-2.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Sasaki C, Hiraishi T, Oku T, Okuma K, Suzumura K, Hashimoto M, Ito H, Aramori I, Hirayama Y. Metabolomic approach to the exploration of biomarkers associated with disease activity in rheumatoid arthritis. Plos One 2019;14(7):e0219400. https://doi.org/10.1371/journal.pone.0219400.

  26. Hambardzumyan K, Bolce RJ, Wallman JK, Van Vollenhoven RF, Saevarsdottir S. Serum biomarkers for prediction of response to methotrexate monotherapy in early rheumatoid arthritis: results from the SWEFOT trial. J Rheumatol. 2019;46(6):555–63. https://doi.org/10.3899/jrheum.180537.

    Article  CAS  PubMed  Google Scholar 

  27. Tasaki S, Suzuki K, Kassai Y, Takeshita M, Murota A, Kondo Y, et al. Multi-omics monitoring of drug response in rheumatoid arthritis in pursuit of molecular remission. Nat Commun. 2018;9(1):1–12.

    Article  CAS  Google Scholar 

  28. Kapoor SR, Filer A, Fitzpatrick MA, Fisher BA, Taylor PC, Buckley CD, et al. Metabolic profiling predicts response to anti–tumor necrosis factor α therapy in patients with rheumatoid arthritis. Arthritis Rheum. 2013;65(6):1448–56. https://doi.org/10.1002/art.37921.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Cuppen BV, Fu J, Van Wietmarschen HA, Harms AC, Koval S, Marijnissen AC, et al. Exploring the inflammatory metabolomic profile to predict response to TNF-α inhibitors in rheumatoid arthritis. Plos One. 2016;11(9):e0163087. https://doi.org/10.1371/journal.pone.0163087.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  30. Liebal UW, Phan AN, Sudhakar M, Raman K, Blank LM. Machine learning applications for mass spectrometry-based metabolomics. Metabolites. 2020;10(6):243. https://doi.org/10.3390/metabo10060243.

    Article  CAS  PubMed Central  Google Scholar 

  31. Welsing PM, van Riel PL. The Nijmegen inception cohort of early rheumatoid arthritis. J Rheumatol Suppl. 2004;69:14–21.

    PubMed  Google Scholar 

  32. Bates D, Mächler M, Bolker B, Walker S. Fitting linear mixed-effects models using lme4. arXiv preprint arXiv:14065823 2014.

    Google Scholar 

  33. Siemons L, ten Klooster PM, Vonkeman HE, van Riel PL, Glas CA, van de Laar MA. How age and sex affect the erythrocyte sedimentation rate and C-reactive protein in early rheumatoid arthritis. BMC Musculoskel Disord. 2014;15(1):368. https://doi.org/10.1186/1471-2474-15-368.

    Article  CAS  Google Scholar 

  34. Radovits B, Fransen J, Van Riel P, Laan R. Influence of age and gender on the 28-joint Disease Activity Score (DAS28) in rheumatoid arthritis. Ann Rheum Dis. 2008;67(8):1127–31. https://doi.org/10.1136/ard.2007.079913.

    Article  CAS  PubMed  Google Scholar 

  35. Ranganath VK, Elashoff DA, Khanna D, Park G, Peter JB, Paulus HE. Age adjustment corrects for apparent differences in erythrocyte sedimentation rate and C-reactive protein values at the onset of seropositive rheumatoid arthritis in younger and older patients. J Rheumatol. 2005;32(6):1040–2.

    PubMed  Google Scholar 

  36. Kuznetsova A, Brockhoff PB, Christensen RHB. Package ‘lmertest’. https://cran.r-project.org/web/packages/lmerTest.

  37. Yang XY, Di Zheng K, Lin K, Zheng G, Zou H, Wang JM, Lin YY, Chuka CM, Ge RS, Zhai W. Energy metabolism disorder as a contributing factor of rheumatoid arthritis: a comparative proteomic and metabolomic study. Plos One 2015;10(7):e0132695. https://doi.org/10.1371/journal.pone.0132695.

  38. Zhou J, Chan L, Zhou S. Trigonelline: a plant alkaloid with therapeutic potential for diabetes and central nervous system disease. Curr Med Chem. 2012;19(21):3523–31. https://doi.org/10.2174/092986712801323171.

    Article  CAS  PubMed  Google Scholar 

  39. Zhou J-Y, Du X-H, Zhang Z, Qian G-S. Trigonelline inhibits inflammation and protects β cells to prevent fetal growth restriction during pregnancy in a mouse model of diabetes. Pharmacology. 2017;100(5-6):209–17. https://doi.org/10.1159/000479088.

    Article  CAS  PubMed  Google Scholar 

  40. Foresti R, Green CJ, Motterlini R. Generation of bile pigments by haem oxygenase: a refined cellular strategy in response to stressful insults. Biochemical Society Symposia. 2004;71:177–192. https://doi.org/10.1042/bss0710177.

  41. Bellner L, Vitto M, Patil KA, Dunn MW, Regan R, Laniado-Schwartzman M. Exacerbated corneal inflammation and neovascularization in the HO-2 null mice is ameliorated by biliverdin. Exp Eye Res. 2008;87(3):268–78. https://doi.org/10.1016/j.exer.2008.06.007.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Nakao A, Otterbein LE, Overhaus M, Sarady JK, Tsung A, Kimizuka K, et al. Biliverdin protects the functional integrity of a transplanted syngeneic small bowel. Gastroenterology. 2004;127(2):595–606. https://doi.org/10.1053/j.gastro.2004.05.059.

    Article  CAS  PubMed  Google Scholar 

  43. Sarady-Andrews JK, Liu F, Gallo D, Nakao A, Overhaus M, Ollinger R, et al. Biliverdin administration protects against endotoxin-induced acute lung injury in rats. Am J Physio Lung Cell Mol Physiol. 2005;289(6):L1131–7.

    Article  CAS  Google Scholar 

  44. Gibbs PE, Maines MD. Biliverdin inhibits activation of NF-κB: reversal of inhibition by human biliverdin reductase. Int J Cancer. 2007;121(11):2567–74. https://doi.org/10.1002/ijc.22978.

    Article  CAS  PubMed  Google Scholar 

  45. Öllinger R, Wang H, Yamashita K, Wegiel B, Thomas M, Margreiter R, et al. Therapeutic applications of bilirubin and biliverdin in transplantation. Antioxid Redox Signal. 2007;9(12):2175–86. https://doi.org/10.1089/ars.2007.1807.

    Article  PubMed  Google Scholar 

  46. Peng Y-F, Wang J-L, Pan G-G. The correlation of serum bilirubin levels with disease activity in patients with rheumatoid arthritis. Clin Chim Acta. 2017;469:187–90. https://doi.org/10.1016/j.cca.2017.04.006.

    Article  CAS  PubMed  Google Scholar 

  47. Fischman D, Valluri A, Gorrepati VS, Murphy ME, Peters I, Cheriyath P. Bilirubin as a protective factor for rheumatoid arthritis: an NHANES study of 2003-2006 data. J Clin Med Res. 2010;2(6):256–60. https://doi.org/10.4021/jocmr444w.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Jangi S, Otterbein L, Robson S. The molecular basis for the immunomodulatory activities of unconjugated bilirubin. Int J Biochem Cell Biol. 2013;45(12):2843–51. https://doi.org/10.1016/j.biocel.2013.09.014.

    Article  CAS  PubMed  Google Scholar 

  49. Enot DP, Lin W, Beckmann M, Parker D, Overy DP, Draper J. Preprocessing, classification modeling and feature selection using flow injection electrospray mass spectrometry metabolite fingerprint data. Nat Protoc. 2008;3(3):446–70. https://doi.org/10.1038/nprot.2007.511.

    Article  CAS  PubMed  Google Scholar 

  50. Grissa D, Pétéra M, Brandolini M, Napoli A, Comte B, Pujos-Guillot E. Feature selection methods for early predictive biomarker discovery using untargeted metabolomic data. Front Mol Biosci. 2016;3:30.

    Article  Google Scholar 

  51. van Zanten JV, Ring C, Carroll D, Kitas GD. Increased C reactive protein in response to acute stress in patients with rheumatoid arthritis. Ann Rheum Dis. 2005;64(9):1299–304. https://doi.org/10.1136/ard.2004.032151.

    Article  Google Scholar 

  52. Graf J, Scherzer R, Grunfeld C, Imboden J. Levels of C-reactive protein associated with high and very high cardiovascular risk are prevalent in patients with rheumatoid arthritis. Plos One. 2009;4(7):e6242. https://doi.org/10.1371/journal.pone.0006242.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Young SP, Kapoor SR, Viant MR, Byrne JJ, Filer A, Buckley CD, et al. The impact of inflammation on metabolomic profiles in patients with arthritis. Arthritis Rheum. 2013;65(8):2015–23. https://doi.org/10.1002/art.38021.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Wells G, Becker J, Teng J, Dougados M, Schiff M, Smolen J, et al. Validation of the 28-joint Disease Activity Score (DAS28) and European League Against Rheumatism response criteria based on C-reactive protein against disease progression in patients with rheumatoid arthritis, and comparison with the DAS28 based on erythrocyte sedimentation rate. Ann Rheum Dis. 2009;68(6):954–60. https://doi.org/10.1136/ard.2007.084459.

    Article  CAS  PubMed  Google Scholar 

  55. Patel CH, Leone RD, Horton MR, Powell JD. Targeting metabolism to regulate immune responses in autoimmunity and cancer. Nat Rev Drug Discov. 2019;18(9):669–88. https://doi.org/10.1038/s41573-019-0032-5.

    Article  CAS  PubMed  Google Scholar 

  56. Bettencourt IA, Powell JD. Targeting metabolism as a novel therapeutic approach to autoimmunity, inflammation, and transplantation. J Immunol. 2017;198(3):999–1005. https://doi.org/10.4049/jimmunol.1601318.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  57. Spies CM, Straub RH, Buttgereit F. Energy metabolism and rheumatic diseases: from cell to organism. Arthritis Res Ther. 2012;14(3):216. https://doi.org/10.1186/ar3885.

    Article  PubMed  PubMed Central  Google Scholar 

  58. Yang Z, Matteson EL, Goronzy JJ, Weyand CM. T-cell metabolism in autoimmune disease. Arthritis Res Ther. 2015;17(1):29. https://doi.org/10.1186/s13075-015-0542-4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  59. Fearon U, Hanlon M, Wade S, Fletcher J. Altered metabolic pathways regulate synovial inflammation in rheumatoid arthritis. Clin Exp Immunol. 2019;197(2):170–80. https://doi.org/10.1111/cei.13228.

    Article  CAS  PubMed  Google Scholar 

  60. Chimenti M, Triggianese P, Conigliaro P, Candi E, Melino G, Perricone R. The interplay between inflammation and metabolism in rheumatoid arthritis. Cell Death Dis. 2015;6(9):e1887. https://doi.org/10.1038/cddis.2015.246.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Sanchez-Lopez E, Cheng A, Guma M. Can metabolic pathways be therapeutic targets in rheumatoid arthritis? J Clin Med. 2019;8(5):753. https://doi.org/10.3390/jcm8050753.

    Article  CAS  PubMed Central  Google Scholar 

  62. Sung J, Wang Y, Chandrasekaran S, Witten DM, Price ND. Molecular signatures from omics data: from chaos to consensus. Biotechnol J. 2012;7(8):946–57. https://doi.org/10.1002/biot.201100305.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  63. Coras R, D Murillo-Saich J, Guma M. Circulating pro-and anti-inflammatory metabolites and its potential role in rheumatoid arthritis pathogenesis. Cells 2020, 9(4):827.

  64. Wegiel B, Otterbein L. Go green: the anti-inflammatory effects of biliverdin reductase. Front Pharmacol. 2012;3:47.

    Article  Google Scholar 

  65. Wegiel B, Gallo D, Csizmadia E, Roger T, Kaczmarek E, Harris C, et al. Biliverdin inhibits Toll-like receptor-4 (TLR4) expression through nitric oxide-dependent nuclear translocation of biliverdin reductase. Proc Natl Acad Sci. 2011;108(46):18849–54. https://doi.org/10.1073/pnas.1108571108.

    Article  CAS  PubMed  Google Scholar 

  66. Vogel ME, Zucker SD. Bilirubin acts as an endogenous regulator of inflammation by disrupting adhesion molecule-mediated leukocyte migration. Inflamm Cell Signal 2016;3(1):e1178.

  67. Liu Y, Li P, Lu J, Xiong W, Oger J, Tetzlaff W, et al. Bilirubin possesses powerful immunomodulatory activity and suppresses experimental autoimmune encephalomyelitis. J Immunol. 2008;181(3):1887–97. https://doi.org/10.4049/jimmunol.181.3.1887.

    Article  CAS  PubMed  Google Scholar 

  68. Tang W, Wu Y, Wu B. Heme catabolic pathway in inflammation and immune disorders. Front Pharmacol. 2019;10:825.

    Article  Google Scholar 

  69. Öllinger R, Yamashita K, Bilban M, Erat A, Kogler P, Thomas M, et al. Bilirubin and biliverdin treatment of atherosclerotic diseases. Cell Cycle. 2007;6(1):39–43. https://doi.org/10.4161/cc.6.1.3700.

    Article  PubMed  Google Scholar 

  70. Pellock SJ, Redinbo MR. Glucuronides in the gut: sugar-driven symbioses between microbe and host. J Biol Chem. 2017;292(21):8569–76. https://doi.org/10.1074/jbc.R116.767434.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Takasuna K, Hagiwara T, Hirohashi M, Kato M, Nomura M, Nagai E, et al. Involvement of β-glucuronidase in intestinal microflora in the intestinal toxicity of the antitumor camptothecin derivative irinotecan hydrochloride (CPT-11) in rats. Cancer Res. 1996;56(16):3752–7.

    CAS  PubMed  Google Scholar 

  72. LoGuidice A, Wallace BD, Bendel L, Redinbo MR, Boelsterli UA. Pharmacologic targeting of bacterial β-glucuronidase alleviates nonsteroidal anti-inflammatory drug-induced enteropathy in mice. J Pharmacol Exp Ther. 2012;341(2):447–54. https://doi.org/10.1124/jpet.111.191122.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Chen J, Wright K, Davis JM, Jeraldo P, Marietta EV, Murray J, et al. An expansion of rare lineage intestinal microbes characterizes rheumatoid arthritis. Genome Med. 2016;8(1):1–14.

    Article  Google Scholar 

  74. He M, Harms AC, van Wijk E, Wang M, Berger R, Koval S, et al. Role of amino acids in rheumatoid arthritis studied by metabolomics. Int J Rheum Dis. 2019;22(1):38–46. https://doi.org/10.1111/1756-185X.13062.

    Article  CAS  PubMed  Google Scholar 

  75. Tilvawala R, Nguyen SH, Maurais AJ, Nemmara VV, Nagar M, Salinger AJ, et al. The rheumatoid arthritis-associated citrullinome. Cell Chem Biol. 2018;25(6):691–704. e696.

    Article  CAS  Google Scholar 

  76. Chu SH, Cui J, Sparks JA, Lu B, Tedeschi SK, Speyer CB, et al. Circulating plasma metabolites and risk of rheumatoid arthritis in the Nurses’ Health Study. Rheumatology. 2020;59(11):3369–79. https://doi.org/10.1093/rheumatology/keaa125.

    Article  CAS  PubMed  Google Scholar 

  77. Marietta EV, Murray JA, Luckey DH, Jeraldo PR, Lamba A, Patel R, et al. Suppression of inflammatory arthritis by human gut-derived Prevotella histicola in humanized mice. Arthritis Rheumatol. 2016;68(12):2878–88. https://doi.org/10.1002/art.39785.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  78. Liu X, Zou Q, Zeng B, Fang Y, Wei H. Analysis of fecal Lactobacillus community structure in patients with early rheumatoid arthritis. Curr Microbiol. 2013;67(2):170–6. https://doi.org/10.1007/s00284-013-0338-1.

    Article  CAS  PubMed  Google Scholar 

  79. Zhang X, Zhang D, Jia H, Feng Q, Wang D, Liang D, et al. The oral and gut microbiomes are perturbed in rheumatoid arthritis and partly normalized after treatment. Nat Med. 2015;21(8):895–905. https://doi.org/10.1038/nm.3914.

    Article  CAS  PubMed  Google Scholar 

  80. Sung J, Kim S, Cabatbat JJT, Jang S, Jin Y-S, Jung GY, et al. Global metabolic interaction network of the human gut microbiota for context-specific community-scale analysis. Nat Commun. 2017;8(1):1–12.

    Article  CAS  Google Scholar 

  81. Wilmanski T, Rappaport N, Earls JC, Magis AT, Manor O, Lovejoy J, et al. Blood metabolome predicts gut microbiome α-diversity in humans. Nat Biotechnol. 2019;37(10):1217–28. https://doi.org/10.1038/s41587-019-0233-9.

    Article  CAS  PubMed  Google Scholar 

  82. Gupta VK, Kim M, Bakshi U, Cunningham KY, Davis JM, Lazaridis KN, Nelson H, Chia N, Sung J. A predictive index for health status using species-level gut microbiome profiling. Nat Commun 2020, 11(1):1-16, 4635, doi: https://doi.org/10.1038/s41467-020-18476-8.

Download references

Acknowledgements

First and foremost, we thank our dear patients who volunteered for this study. We also thank the Mayo Clinic Division of Rheumatology study coordinators (Jennifer Sletten and Kathleen McCarthy-Fruin) for their help in making this work possible.

Funding

This work was supported in part by the Mayo Clinic Center for Individualized Medicine (to B.H., V.K.G., and J.S.) and Mark E. and Mary A. Davis to Mayo Clinic Center for Individualized Medicine (J.M.D. and J.S.).

Author information

Authors and Affiliations

Authors

Contributions

B.H., J.M.D., and J.S. conceived the problem and designed all the analytical methodologies. B.H. performed all computational analyses. All authors analyzed the data. B.H., J.M.D., and J.S. wrote the manuscript, with editorial contributions from other authors. J.M.D. is the principal investigator of the Mayo Clinic Rheumatology Biobank, from which blood samples were collected from patients with rheumatoid arthritis. All authors reviewed and approved the final manuscript.

Corresponding author

Correspondence to Jaeyun Sung.

Ethics declarations

Ethics approval and consent to participate

This study was approved by the Mayo Clinic Institutional Review Board (no. 14-000616 and no. 14-000680) in accordance with the Declaration of Helsinki. All methods and procedures were performed in accordance with the Mayo Clinic Institutional Review Board guidelines and regulations.

Consent for publication

All patients provided written informed consent.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Supplementary Figure 1

. Histogram of DAS28-CRPs corresponding to the 128 total samples of the discovery cohort.

Additional file 2: Supplementary Table 1

. Subject characteristics of 64 patients (128 samples) from the discovery cohort.

Additional file 3: Supplementary Table 2

. Subject characteristics of 12 patients (12 samples) from the validation cohort.

Additional file 4.

Supplementary Information on detailed methods regarding metabolomic profiling (Metabolon).

Additional file 5: Supplementary Table 3

. Demographic characteristics of discovery cohort samples (n = 128) divided into lower (DAS28-CRP ≤ 3.2) and higher disease activity (DAS28-CRP > 3.2).

Additional file 6: Supplementary Table 4

. Differentially abundant metabolites between lower (DAS28-CRP ≤ 3.2) and higher (DAS28-CRP > 3.2) disease activity groups.

Additional file 7: Supplementary Table 5

. Metabolites displaying significant changes in abundances in patients with clinical improvement (n = 16).

Additional file 8: Supplementary Table 6

. Metabolites displaying significant changes in abundances in patients without clinical improvement (n = 48).

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Hur, B., Gupta, V.K., Huang, H. et al. Plasma metabolomic profiling in patients with rheumatoid arthritis identifies biochemical features predictive of quantitative disease activity. Arthritis Res Ther 23, 164 (2021). https://doi.org/10.1186/s13075-021-02537-4

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s13075-021-02537-4

Keywords