Functionally distinct disease-associated fibroblast subsets in rheumatoid arthritis

Abstract

Fibroblasts regulate tissue homeostasis, coordinate inflammatory responses, and mediate tissue damage. In rheumatoid arthritis (RA), synovial fibroblasts maintain chronic inflammation which leads to joint destruction. Little is known about fibroblast heterogeneity or if aberrations in fibroblast subsets relate to pathology. Here, we show functional and transcriptional differences between fibroblast subsets from human synovial tissues using bulk transcriptomics of targeted subpopulations and single-cell transcriptomics. We identify seven fibroblast subsets with distinct surface protein phenotypes, and collapse them into three subsets by integrating transcriptomic data. One fibroblast subset, characterized by the expression of proteins podoplanin, THY1 membrane glycoprotein and cadherin-11, but lacking CD34, is threefold expanded in patients with RA relative to patients with osteoarthritis. These fibroblasts localize to the perivascular zone in inflamed synovium, secrete proinflammatory cytokines, are proliferative, and have an in vitro phenotype characteristic of invasive cells. Our strategy may be used as a template to identify pathogenic stromal cellular subsets in other complex diseases.

Fig. 1

Distinct protein and mRNA expression between fibroblast subsets. a Gating strategy for synovial fibroblasts with heterogeneous expression of surface proteins. b Analysis of variance (ANOVA) reveals 436 genes with significant (1% FDR) variation across seven gated populations that are measured and statistically significant in both microarray and RNA-seq datasets. Each column in the heatmap corresponds to the average of multiple samples of a cell sorting gate. Each square beneath a column represents a donor from which this sample was taken. c Principal component analysis (PCA) with 2,986 genes (1% FDR, ANOVA) in microarray data separates the 32 microarray samples into three subsets: CD34^–THY1^–, CD34^–THY1⁺, and CD34⁺. d Pairwise Pearson's correlation of microarrays also suggests three major subsets of fibroblasts

Fig. 2

Unbiased clustering of fibroblasts by single-cell RNA-seq. a Hierarchical clustering is concordant with our surface protein gating. The first track indicates which donor the cell came from. The second track indicates which cell sorting gate the cell falls into, by surface protein expression. The heatmap indicates mRNA expression, scaled across cells to reveal variation between cells. Below the heatmap, the first track shows linear discriminant analysis (LDA) classification by mRNA expression of 968 genes. The last track shows posterior probabilities of LDA classification. b For the heatmap, we selected 23 genes with greatest 1% mean and greatest 1% variance in expression levels. c The number of cells in each subset as determined by surface protein levels (x-axis) and by LDA classification (y-axis) is consistent. RMSE, root mean squared error

Fig. 3

Anatomical localization and proportions of fibroblast subsets. a Anatomical localization of fibroblast subsets and leukocytes in RA and OA synovial tissue. Hoechst 33258: white, CD45: cyan, PDPN: blue, CD34: green, THY1: red. Small arrowheads: CD34^–THY1^– fibroblasts. Big arrowheads: CD34^–THY1⁺ fibroblasts. Arrows: CD34⁺ fibroblasts. Scale = 100 µm. b Expression of fibroblast subset markers in RA and OA synovial tissue. Hoechst 33258: white, CDH11: cyan, PDPN: blue, CD34: green, THY1: red. c Proportions of fibroblast subsets in synovial tissue in OA (n = 26) and RA (n = 16) evaluated by flow cytometry. P values from the Wilcoxon's rank-sum test. d Proportions of CDH11⁺ cells in CD34^–THY1^– fibroblasts, CD34^–THY1⁺ fibroblasts, and CD34⁺ fibroblasts. L lining area, S sublining area

Fig. 4

Correlation between proportion of fibroblast subsets and synovial inflammation. a The Spearman's correlation between the proportion of CD34^–THY1⁺ fibroblasts in total fibroblasts and the proportion of CD45⁺ cells in total live cells in synovial tissue biopsies from RA knee joints. b–d The linear correlation coefficient between the proportion of CD34^–THY1⁺ fibroblasts and Krenn inflammation score (b), ultrasound synovial hypertrophy score (c), and disease duration (years) (d). Proportion of cells were evaluated by flow cytometry. Plots b–d show linear model coefficient estimates and ttest P values

Fig. 5

Functional differences between fibroblast subsets. a Gene expression pattern of key effector molecules to predict the function of fibroblast subsets. b Gene set enrichment analysis with RNA-seq samples and MSigDB Hallmark Gene Signatures. Terms with <5% FDR are displayed. c The proportion of Ki67-positive cells in freshly isolated synovial fibroblast subsets. d Protein expression of MMP14 in the cell lysates of expanded fibroblast subsets in the presence of 1 ng/mL TNF-α. e The number of migrated and invaded cells in transwell matrix invasion assay. Fifty nanogram per mL of PDGFBB was used to promote migration and invasion. f The ratio of TNFSF11 (RANKL) and TNFRSF11B (OPG) in mRNA level in expanded fibroblast subsets. g Protein expression of TNFRSF11B in the supernatant of expanded fibroblast subsets. h Number of osteoclasts developed in the co-culture of osteoclast progenitors with fibroblast subsets in the presence of 20 ng/mL M-CSF and 5 ng/mL RANKL. i Protein expression of cytokines and MMPs in the supernatant of expanded fibroblast subsets in the presence of 1 ng/mL TNF-α. j The number of migrated monocytes in monocyte recruitment assay using supernatant of cultured fibroblast subsets. In bar graphs, blue, red and green color coding corresponds to the fibroblast subsets shown in a. White bar in j indicates culture media without fibroblasts. Bars show means and error bars show standard deviations. One-way ANOVA, Tukey's post hoc: *P < 0.05; ** P < 0.005