In this study, we enrolled 93 RA patients and 16 healthy subjects ( Table 1 ). Patients with RA fulfilled the 2010 American College of Rheumatology and European League Against Rheumatism classification criteria (22). RA patients were classified based on RF levels into RF-negative (n = 16), RF-low (n = 24), and RF-high (n = 53) groups. Patients with RF ≤20 IU/ml were considered RF-negative, those with RF levels between 20 and 60 IU/ml were considered RF-low, and those with RF >60 IU/ml were considered RF-high. Subjects who had recently used antibiotics, probiotics, or prebiotics were excluded from the study. The demographic and clinical characteristics of the study subjects are summarized in Table 1 . The age of healthy subjects (37–62) and RA patients including RF-negative (41–57), RF-low (40–59), and RF-high (40–60) was not significantly different ( Table 1 ). The ratio of female among the subjects is 94% from healthy control subjects, 100% from RF-negative, 92% from RF-low, and 85% from RF-high. Additionally, disease activity score in 28 joints (DAS28), erythrocyte sedimentation rate (ESR), C-reactive protein levels, and current use of methotrexate or a biological disease-modifying anti-rheumatic drug were not statistically different between RF-negative, RF-low, and RF-high RA patients. The study design was approved by the Institutional Review Board of Seoul St. Mary’s Hospital, The Catholic University of Korea (approval ID: KC17TNSI0570). Written informed consent was obtained from all study participants.

Individual human fecal samples were sent immediately with ice to the research site after collection in a plastic container and were stored at −70°C within 12 h of arrival. Isolation of the bacterial genomic DNA from the fecal samples was performed using a QIAamp DNA Stool Mini Kit (Qiagen, Hilden, Germany) following the instructions of the manufacturer, followed by bead beating on a TissueLyser system (Qiagen). Quantification of the bacterial genomic DNA was carried out using a Qubit 3.0 Fluorometer (Thermo Fisher Scientific, Waltham, MA, USA). For sequencing, 16S rRNA gene amplification and index PCR were conducted following the Illumina 16S Metagenomic Sequencing Library preparation guide (Illumina, San Diego, CA, USA). For amplification of V3 and V4 regions, the 16S sequence was amplified using forward primer (5′-TCGTCGGCAGCGTCAGATGTGTATAAGAGACAGCCTACGGGNGGCWGCAG-3′) and reverse primer (5′-GTCTCGTGGGCTCGGAGATGTGTATAAGAGACAGGACTACHVGGGTATCTAATCC-3′). For indexing, Nextera XT Index 1 and 2 Primers from a Nextera XT Index kit (Illumina) were used. Each PCR product was purified using AMPure XP beads (Beckman Coulter, Pasadena, CA, USA). DNA sequencing was conducted in using the paired-end method (300 bpx2) with an Illumina MiSeq instrument according to the Illumina protocol. Raw 16S rRNA sequences were bioinformatically analyzed using QIIME 2 version 2019.4 (23), as previously described (24).

Transgalactooligosaccharide (TOS) antibiotic mupirocin agar medium, a selective medium for Bifidobacterium, was used to separate Bifidobacterium from 100 mg of human fecal samples. After plating fecal samples diluted with phosphate-buffered saline (PBS) and culturing for 30 h at 37°C in an anaerobic jar under anaerobic condition, colonies were randomly taken from a plate from each sample and inoculated with cysteine De Man, Rogosa, and Sharpe (cMRS) broth and cultured for 24 h. The above process was performed in a glove box for blocking oxygen. Identification of bacteria strain was done with 16S rRNA sequencing service of Macrogen Inc. (Seoul, Korea).

Bacterial samples (50µl) were inoculated into 5 L of cMRS liquid medium and incubated for 20 h at 37°C. After incubation, the bacteria were collected via centrifugation (2236HR high-speed centrifuge; Gyrozen, Gimpo, Korea) at 20°C and 7,000 rpm. After two washes with PBS, the bacterial samples were dried at 36°C and centrifuged at 2,000 rpm for 24 h using a Scanvac Speed Vacuum Concentrator (LaboGene Aps, Lillerød, Denmark). All bacteria used in this study were isolated from healthy subjects and RA patients. Donor characteristics are summarized in Table 1 .

Seven-week-old male DBA/1J mice (Orient Bio, Gyeonggi-do, Korea) and 7-week-old male NOD-scid IL2rγ^null (NSG) mice (Jackson Laboratory, Bar Harbor, ME, USA) were maintained under specific pathogen-free conditions in an animal facility with controlled humidity (55 ± 5%), light (12 h/12 h light/dark), and temperature (22 ± 1°C). The air in the facility passed through a HEPA filter system designed to exclude bacteria and viruses. Animals were fed mice chow (normal chow and 60 kcal chow) and water ad libitum. All experimental procedures were approved by the Institutional Animal Care and Use Committee at the School of Medicine and the Animal Research Ethics Committee of the Catholic University of Korea and were conducted in accordance with the Laboratory Animals Welfare Act according to the Guide for the Care and Use of Laboratory Animals. All experimental procedures were evaluated and conducted in accordance with the protocols approved by the Animal Research Ethics Committee at the Catholic University of Korea (Permit Number: CUMC 2019-0242-01). All procedures performed in this study followed the ethical guidelines for animal use.

Collagen-induced arthritis (CIA) was established in male DBA/1J mice. Briefly, mice were inoculated with 100 µg of chicken type II collagen (CII; Chondrex Inc., Redmond, WA, USA) liquified overnight in 0.1 N acetic acid (4 mg/ml) in complete Freund’s adjuvant. Immunizations were performed intradermally via injection into the bottom of the tail. At 2 weeks following immunization, the mice were boosted with 100 µg of CII in incomplete Freund’s adjuvant (Chondrex Inc.). In addition, the obese CIA mouse group was fed a high-fat diet (60 kcal of fat) at primary immunization. A HFD was maintained until the sacrifice.

To induce arthritis in humanized mice, 5 × 10⁵ peripheral blood mononuclear cells (PBMCs) of RA patients were injected into the tail vein. At 2 weeks after, orbital blood was collected to confirm engraftment by flow cytometry. Orbital blood was immunostained with CD4 T cells and CD8 T cells, and the cells were stained with Alexa Flour 700-conjugated anti-human CD4 antibody (BD Biosciences, San Diego, CA, USA) and allophycocyanin (APC)-Cy7-conjugated anti-human CD8 (BD Biosciences). Thereafter, 1 week after confirmation of engraftment, the mice were sensitized with 50 µg of lipopolysaccharide (Sigma Aldrich, St. Louis, MO, USA). Each in vivo experiment was repeated a total of three times.

Arthritis severity was determined using the mean arthritis index, which ranges from 0 to 4, as follows: 0, no evidence of erythema or swelling; 1, erythema and mild swelling confined to the midfoot (tarsals) or ankle joint; 2, erythema and mild swelling extending from the ankle to the midfoot; 3, erythema and moderate swelling extending from the ankle to the metatarsal joint; and 4, erythema and severe swelling encompassing the ankle, foot, and digits. The severity of arthritis was determined as the sum of scores from all legs assessed by two independent observers blinded to the experimental groups.

Mouse joint tissues were fixed in 4% paraformaldehyde (Sigma-Aldrich, St. Louis, MO, USA), decalcified in a histological decalcifying solution (Calci-Clear Rapid; National Diagnostics, Atlanta, GA, USA), and embedded in paraffin wax for histological analyses. Sections (7 μm) were prepared and stained with hematoxylin (YD Diagnostics, Yongin, Korea), eosin (Muto Pure Chemicals Co., Ltd., Tokyo, Japan), and safranin O (Sigma-Aldrich). Cartilage damage was marked as described previously (25).

Joint tissues were incubated overnight at 4°C with primary antibodies against TNFα (Abcam, Cambridge, UK), IL-17 (Abcam), IL-6 (Abcam), and IL-1β (Novus Biologicals, Littleton, CO, USA). Subsequently, samples were incubated with a biotinylated streptavidin–peroxidase complex for 1 h, and the signals were developed using chromogen 3,3′-diaminobenzidine (Thermo Scientific, Rockford, IL, USA). The sections were examined under a photomicroscope (Olympus, Tokyo, Japan). The number of positive cells in high-power digital images (magnification, ×400) was counted using Adobe Photoshop software (Adobe, San Jose, CA, USA). Stained cells were counted independently by three observers, and the mean values were evaluated.

Serum was gathered 7 weeks after the first immunization for concentration of the IgA, IgM, IgG2a, and collagen-specific IgG2a levels, and the serum was stored at −70°C until further use. The levels of the IgA, IgM, IgG2a, and collagen-specific IgG2a antibodies in the serum were measured using ELISA (Bethyl Laboratories, Montgomery, TX, USA).

Cytokine expression in mice was analyzed via intracellular staining with the following antibodies: fluorescein isothiocyanate-conjugated anti-IL-17, APC-conjugated anti- interferon-gamma (IFN-γ), and PerCP-conjugated anti-CD4. Before staining, cells were stimulated for 4 h with phorbol myristate and ionomycin in the presence of GolgiStop (BD Biosciences). To analyze regulatory T (Treg) populations, splenocytes were stained with PerCP-conjugated anti-CD4, APC-conjugated anti-CD25, and PE-conjugated anti-Foxp3 antibodies. T helper and Treg populations were determined using specific antibodies (eBioscience, San Diego, CA, USA). After surface staining for 30 min, the cells were permeabilized with Cytofix/Cytoperm solution (BD Biosciences). Thereafter, the cells were intracellularly stained with fluorescent antibodies. Flow cytometry was performed with the aid of a cytoFLEX Flow Cytometer (Beckman Coulter, Brea, CA, USA), and flow cytometry data were analyzed using FlowJo (Tree Star, Ashland, OR, USA).

PBMCs were isolated from healthy donors. Briefly, blood from donors was mixed with PBS, and the mixture was cautiously layered onto 10 ml of Ficoll Plaque Plus solution (GE Healthcare Life Sciences, Marlborough, MA, USA). Blood samples were centrifuged for 30 min at 2,500 rpm and 20°C to separate blood contents. After centrifugation, the layer containing PBMCs was collected, and the cells were washed with PBS and cultured in RPMI medium containing 10% fetal bovine serum. PBMCs were incubated with plate-bound anti-CD3 (0.5 μg/ml) and then subjected to IL-17 analysis by ELISA for 3 days or subjected to RNA sequencing and real-time polymerase chain reaction analysis for 2 days.

The supernatant was collected 3 days after B. longum RAPO or vehicle treatment with anti-CD3 (0.5 μg/ml). IL-17 was analyzed by sandwich enzyme-linked immunosorbent assay (ELISA; IL-17 DuoSet ELISA; R&D Systems, Lille, France). The absorbance at 450 nm was determined using an ELISA microplate reader (Molecular Devices).

We preprocessed the raw reads from the sequencer to remove low-quality and adapter sequence before analysis and aligned the processed reads to Homo sapiens (hg38) using HISAT v2.1.0 (26). HISAT utilizes two types of indexes for alignment (a global, whole-genome index and tens of thousands of small local indexes). These two types of indexes are constructed using the same BWT (Burrows–Wheeler transform) and a graph FM index (GFM) as Bowtie2. Because of its use of these efficient data structures and algorithms, HISAT generates spliced alignments several times faster than Bowtie and BWT that are widely used. The reference genome sequence of Homo sapiens (hg38) and the annotation data were downloaded from the UCSC table browser (http://genome.uscs.edu). Transcript assembly and abundance estimation used StringTie (27, 28). After alignment, StringTie v1.3.4d was used to assemble aligned reads into transcripts and to estimate their abundance. It provides the relative abundance estimates as read count values of transcript and gene expressed in each sample. Also, transcript assembly of known transcripts, novel transcripts, and alternative splicing transcripts was processed by StringTie v2.1.3b.

The relative abundances of gene were measured in read count using StringTie. We performed statistical analysis to find differentially expressed genes using the estimates of abundances for each gene in the samples. Genes with one more than zeroed read count values in the samples were excluded. To facilitate log2 transformation, 1 was added to each read count value of filtered genes. Filtered data were log2-transformed and subjected to TMM normalization. Statistical significance of the differential expression data was determined using exactTest using edgeR and fold change in which the null hypothesis was that no difference exists among groups. False discovery rate (FDR) was controlled by adjusting p-value using the Benjamini–Hochberg algorithm. For the DEG set, hierarchical clustering analysis was performed using complete linkage and Euclidean distance as a measure of similarity. Gene enrichment and functional annotation analysis and pathway analysis for significant gene list were performed based on gProfiler (https://biit.cs.ut.ee/gprofiler/gost) and KEGG pathway (http://www.genome.jp/kegg/pathway.html).

Total RNA was extracted using TRI reagent (Molecular Research Center, USA) and cDNA was synthesized with the Dyne First Strand cDNA Synthesis Kit (Dyne Bio, South Korea) according to the protocol of the manufacturer. mRNA was quantified using the StepOnePlus™ Real-Time PCR systems (Applied Biosystems, USA) with SensiFAST SYBR Hi-ROX (Bioline, USA). mRNA levels were normalized to that of β-actin. The following primers were used: β-actin 5′-GGACTT CGAGCAAGAGATGG-3′ and 5′-TGTGTTGGGGTACAGG TCTTTG-3′; IL-17A 5′-CAACCGATCCACCTCACC TT-3′ and 5′-GGCACTTTGCCTCCCAGAT-3′; RORC 5′-AGTCGGAAGGCAAGATCAGA-3′ and 5′-CAAGAGAGGTTCTGG GCAAG-3′; IL-21 5′-TGTGAATGACTTGGACCCTGAA-3′ and 5′-AAACAG GAAATAGCTGACCACTCA-3′; IRF-4 5′-CCTGCAAGCTCTTTGACACA-3′ and 5′-GAGTCACCTGGAATCTTGGC-3′; and IL-23R 5′-AGGTACTGGCAGCCTTGGAGTT-3′ and 5′-CCCTGTAGAGATGGAAGCAACTG-3′.

Non-parametric statistical tests (Mann–Whitney test and Kruskal–Wallis test) were used to compare the alpha diversity of each group and the relative abundance of each taxon. Pairwise permutational analysis of variance (PERMANOVA) was used to compare beta diversity between groups. Spearman correlation analysis was used to analyze the relationship between the relative abundance of each taxon and RF levels. Data were analyzed using GraphPad Prism (V.5 for Windows; GraphPad Software Inc., La Jolla, CA, USA) and are presented as means ± standard deviations. Comparisons of numerical data between two groups were performed using Student’s t-test or the Mann–Whitney U-test. Differences in the mean values of more than two groups were assessed using analysis of variance with a post-hoc test. p-values <0.05 (two-tailed) were considered to indicate statistical significance.

To examine the relationship between RF levels and changes in intestinal microbial diversity, we performed gut microbiome profiling through 16S rRNA sequencing of fecal samples from 77 RF-positive patients, 16 RF-negative patients, and 16 healthy volunteers ( Table 1 ). Raw sequences were filtered and merged, and 1,284,582 reads were grouped into operational taxonomic units (OTUs) based on 97% similarity, with a mean of 11,059 across 110 samples. Alpha diversity (Faith’s PD), OTU number, Pielou’s evenness, and Shannon’s diversity index did not differ significantly between the RF-negative, RF-low, and RF-high groups ( Supplementary Figures 1A–D ). However, there were significant differences in beta diversity between RF-negative and RF-low RA patients and healthy subjects ( Supplementary Figures 1F, G ). Notably, RF-negative RA patients had a significantly lower Shannon index compared with healthy subjects ( Supplementary Figure 1A ). RF-low individuals exhibited a significantly lower Shannon index and OTU number compared with healthy volunteers ( Supplementary Figures 1A, B ). Jaccard distance-based beta diversity analysis revealed significant differences in the microbial community between healthy subjects and RF-low RA patients (p = 0.027, PERMANOVA; Supplementary Figure 1F ), as well as between healthy subjects and RF-high RA patients (p = 0.010, PERMANOVA; Supplementary Figure 1F ). Unweighted UniFrac distance-based microbial structure analysis showed significant differences between healthy subjects and RF-negative RA patients (p = 0.045, PERMANOVA; Supplementary Figure 1G ), as well as between healthy subjects and RF-high RA patients (p = 0.003, PERMANOVA; Supplementary Figure 1G ). Other beta diversity matrices did not differ significantly among the groups ( Supplementary Figures 1E, H ).

Although there were no differences in microbial diversity in RA patients according to RF level, the Kruskal–Wallis test revealed significant differences in the composition of individual taxa. We identified 11 differentially abundant taxa at the false discovery rate of 10%; 10 taxa were enriched in healthy volunteers, and one taxon was enriched in RA patients ( Table 2 ). At the phylum level, Actinobacteria were less abundant in RF-high RA patients than in healthy individuals and RF-negative RA patients ( Figures 1A, B ). The relative abundance of Actinobacteria in all RA patients was negatively correlated with RF level (R = −0.29, p = 00051; Figure 1C ). By contrast, the abundance of the phylum Firmicutes in all RA patients was positively correlated with RF levels (R = 0.22, p = 0.0322; Figure 1D ). Within the Actinobacteria phylum, Bifidobacterium were significantly more abundant in RF-negative than in RF-high RA patients ( Figure 2A ). Another genus of the Actinobacteria phylum, Collinsella, was more abundant in healthy volunteers than in RF-negative or RF-high RA patients; however, Collinsella abundance did not differ significantly among RA patients ( Figure 2B ). Similarly, Clostridium species were more abundant in healthy individuals than in the different RA patient groups, although their abundance was similar across RA patient groups ( Figure 2C ). Compared with RF-negative RA patients and healthy individuals, RF-high RA patients exhibited a lower prevalence of the genera Bifidobacterium, Collinsella, and Clostridium ( Supplementary Figures 3A–C ). The prevalence of Veillonellaceae family was the highest among the control and RA patients ( Supplementary Figure 3D ).

Among the genera of the Actinobacteria phylum, the genus Bifidobacterium exhibited a negative correlation with RF levels (R = −0.22, p = 0.0326; Figure 2D ). Furthermore, the abundance of the genus Ruminococcus was positively correlated with RF levels (R = 0.25, p = 0.0161; Figure 2D ).

Bifidobacterium genera can be cultivated despite being anaerobes and are widely known for immune regulation ability. Hence, we decided to focus on Bifidobacterium genus among the differentially abundant genera for further study. To identify the most potent Bifidobacterium strain among the Bifidobacterium isolated from healthy individuals and RA patients, we treated PBMCs (stimulated with anti-CD3) with B. pseudocatenulatum 20T2, B. pseudocatenulatum 20T1, B. pseudocatenulatum 18T6, B. pseudocatenulatum 5T4, B. longum 4L7, B. longum 4L6, B. adolescentis 8T3, B. adolescentis 8T1, B. bifidum 7T2, and B. longum RAPO. Among these strains, B. longum RAPO exhibited the most potent inhibitory effect on IL-17 secretion ( Figure 3 ).

To evaluate the effects of B. longum RAPO in CIA, we orally administrated 1 × 10⁸ CFU/mouse B. longum RAPO alone, methotrexate 3 mg/kg alone, or vehicle to CIA mice at 3 weeks after immunization. Mice developed RA at 21 days after immunization; however, oral administration of B. longum RAPO to CIA mice decreased CIA, as demonstrated by a significantly reduced arthritis score compared with the vehicle-treated mice ( Figure 4A ). In mice that were given oral administration of B. longum RAPO, the severity of arthritis damage was attenuated, and cartilage protection was improved, compared with the control mice. Inflammation score, bone damage, and cartilage damage were also significantly decreased in mice treated with B. longum RAPO ( Figure 4B ). Importantly, the frequencies of CD4⁺ T cells producing IFN-γ and IL-17 were significantly lower in mice treated with B. longum RAPO than in control mice treated with the vehicle. Furthermore, administration with B. longum RAPO resulted in an increased Treg population, compared with the control ( Figure 4C ). Serum IgA, rheumatoid factor IgM, IgG2a, and CII-specific IgG2a levels were reduced by the administration of B. longum RAPO ( Figure 4D ). *, p < 0.05; **, p < 0.01; ***, p < 0.001 (vs. vehicle-treated group).

An obese arthritis model was reported according to a previous report (29). Obese CIA mice had a more serious state of disease and a higher arthritis score compared with CIA. To investigate the effect of B. longum RAPO in obese CIA, obese CIA-induced RA mice were fed the HFD for 7 weeks. There was a significant decrease in the arthritis score and incidence in the B. longum RAPO-administered mice compared with the vehicle-treated obese CIA mice ( Figure 5A ). We then performed effector T-cell staining on splenocytes from B. longum RAPO- and vehicle-treated mice. Flow cytometry analysis of splenocytes showed that the Th17 (IL-17⁺ in CD4⁺ T cells) population was reduced by the administration of B. longum RAPO; on the other hand, the populations of Th1 (IFN-γ⁺ in CD4 T cells) and Th2 (IL-4⁺ in CD4⁺ T cells) did not differ significantly between B. longum RAPO- and vehicle-treated mice with obese CIA ( Figure 5B ). In mice that received B. longum RAPO, the severity of arthritis was attenuated, and cartilage conservation was improved, compared with the vehicle-treated mice ( Figure 5C ). We also analyzed B. longum RAPO effects on inflammatory cytokine staining of joints which showed that B. longum RAPO treatment diminished the levels of IL-1β, IL-6, IL-17, and TNFα ( Figure 5D ). *, p < 0.05; **, p < 0.01; ***, p < 0.001 (vs. vehicle-treated group).

The PBMCs of RA patients were cultured under anti-CD3 0.5 μg/ml conditions with B. longum RAPO for 72 h. Flow cytometry showed that B. longum RAPO treatment significantly suppressed Th1 and Th17 cell proliferation ( Figure 6A ). We next analyzed by RNA sequencing the gene expression profiles of RA PBMC treatment with B. longum RAPO in the presence of anti-CD3 0.5 μg/ml for 48 h. When demonstrated by absent/present classification and using at least two-fold difference in expression as the cutoff, 5,089 genes were differentially expressed in vehicle-treated RA PBMCs versus B. longum RAPO-treated RA PBMCs. Among those in the Th17 pathway, Il17, IRF4, RORC, SIl21, and Il23 were downregulated by B. longum RAPO ( Figure 6B ). Also, Th17-related genes confirmed that Il17, IRF4, RORC, SIl21, and Il23 were downregulated by qPCR ( Figure 6C ). *, p < 0.05; **, p < 0.01; ***, p < 0.001 (vs. vehicle-treated group).

Arthritis score and incidence were examined to investigate whether oral administration of B. longum RAPO can have an effect on human avatar arthritis mice. The engraftment of the patient PBMC was confirmed by analyzing the expression of CD4 and CD8 T cells using human antibodies ( Supplementary Figure 4D ). The B. longum RAPO-treated group had a lower arthritis score and incidence than the vehicle group ( Figure 6D ).

Also, it was confirmed that Th17 cells were effectively inhibited like CIA and obese CIA ( Figure 6E ). Histological analysis was performed using safranin O staining to confirm the degree of cartilage destruction in avatar mice of RA patient. The cartilage damage score was low in the B. longum RAPO-administered group ( Figure 6F ). **, p < 0.01; ***, p < 0.001 (vs. vehicle-treated group).