Cellular composition of the IFP

IFPs from 21 study participants (6 healthy control and 15 KOA) were analysed using snRNA-seq and bioinformatic analyses to identify cell types11 (figure 1A; online supplemental figure 1,2,3A, online supplemental table 1). A total of 82 924 nuclei were sequenced and, after filtering, 73 808 nuclei were analysed. Eight cell types were identified based on canonical markers: fibroblasts, adipocytes, macrophages, endothelial cells, dendritic cells, smooth muscle cells, lymphocytes and mast cells (figure 1B–D; online supplemental figure 3A). The major cell populations contributing to the IFP (>10% of nuclei analysed) were fibroblasts (44.35%), macrophages (19.44%), adipocytes (16.41%) and endothelial cells (12.12%) (figure 1C,D).

• Download figure
• Open in new tab
• Download powerpoint

Figure 1

Infrapatellar fat pad (IFP) cell populations and fibroblast subclusters identified by single-nucleus RNA sequencing (n=21). (A) Schematic workflow to identify cell types and subclusters within human IFPs (created with BioRender.com). (B) Uniform manifold approximation and projection (UMAP) of clusters identified in IFPs. (C) Dot plot with the average and percent population expression of canonical markers for each cell population. (D) Heatmap displaying expression of the top two differentially expressed genes (DEGs) within each cell population. (E) UMAP of fibroblast subclusters and bar plot depicting fibroblast proportion of each subcluster. Top five unique DEGs defining each subcluster are indicated. (F) Violin and density plots displaying the expression of the top two unique DEGs within each fibroblast subcluster. Density plots display expression density ranging from dark purple (low) to yellow (high). (G) Two-dimensional plot showing inferred trajectory across pseudotime of integrated fibroblast subclusters using SCORPIUS. Heatmap demonstrating the change in average scaled expression of the top 94 important genes defining fibroblast subcluster trajectory across pseudotime (q<0.05). Genes in the heatmap are clustered into four modules, each representing a set of genes expressed at different time points of the trajectory. Jitter plots depict the expression of the top marker gene for distinct subclusters across pseudotime (see online supplemental figures 2–4, online supplemental figure 6 and online supplemental table 2A).

Cell subclusters and transcriptomic profiles of the IFP

We next aimed to determine cell subclusters of each major cell type identified within the IFP (n=21).11 Subclustering fibroblasts determined five subclusters, with proportions ranging from 30.32% (subcluster 0) to 9.63% (subcluster 4) (total nuclei: 32 862, resolution 0.2) (figure 1E; online supplemental figure 3B). Unique transcriptomic profiles for each fibroblast subcluster were generated by identifying differentially expressed genes (DEGs) in one subcluster as compared to others. These unique profiles ranged from 225 genes (subcluster 2) to 43 genes (subcluster 3) [q<0.05, log2 fold change (FC) ≥0.5, min.pct≥0.25] (online supplemental figure 4; online supplemental table 2A). The top five genes for each fibroblast subcluster, based on decreasing Log2FC, included subcluster 0: FBN1, PXDNL, PIEZO2, GUCY1A2 and PTGIS; subcluster 1: ABCA10, KCND2, ABCA6, ABCA8 and ABCA9; subcluster 2: CRTAC1, CLIC5, PRG4, FN1 and SEMA3A; subcluster 3: GPAM, FABP4, F13A1, FTL and FTH1 and subcluster 4: MGAT4C, GRIP1, BTBD11, IGF1 and C1GALT1 (figure 1F, online supplemental figure 3C; online supplemental table 2A). We then used flow cytometry to identify the major fibroblast subclusters at the protein level using an additional n=12 IFP samples (online supplemental figure 5A). Our flow cytometry panel included CD26/DPP4, CD10/MME/CALLA and HLA-C, cell surface proteins that we identified as key markers for the three major fibroblast subclusters, subcluster 0, subcluster 1 and subcluster 2, respectively (online supplemental figure 5B). These markers were assessed on fibroblasts (PDPN+CD45−), revealing subcluster 0 and subcluster 1 populations to be the most prevalent (online supplemental figure 5C,D), as we observed in our snRNA-seq dataset.

Next, pseudotime trajectory analysis using SCORPIUS V.1.0.912 showed a clear transition of fibroblast subclusters from subcluster 0 to subcluster 2, through subclusters 3, 1 and 4 (figure 1G). We identified genes that predicted the ordering of fibroblast subclusters along pseudotime, ranked by importance scores. There were distinct changes in module expression profiles, corresponding to transcriptomic profiles of fibroblast subclusters (figure 1G). Additionally, top markers of subcluster 0 (FBN1), subcluster 1 (ABCA10), subcluster 2 (CRTAC1) and subcluster 4 (MGAT4C) showed peak expression in respective subclusters along pseudotime (figure 1G). Fibroblasts expressing DPP4, PI16 and CD34 have previously been characterised as universal fibroblasts able to differentiate into more specialised fibroblast subtypes.13 In our dataset, these markers were highly expressed in subcluster 0, suggesting that it may act as precursors for other fibroblast subclusters we identified (figure 1G; online supplemental figure 6).

Using a similar bioinformatic strategy, we uncovered four subclusters of macrophages within the IFP, with proportions ranging from 49.85% (subcluster 0) to 8.06% (subcluster 3) (total nuclei: 14 351, resolution 0.1) (figure 2A, online supplemental figure 7A). Each macrophage subcluster had unique transcriptomic profiles spanning from 232 genes (subcluster 1) to 74 genes (subcluster 3) (online supplemental figure 8; online supplemental table 2B). The top five genes for each macrophage subcluster based on decreasing log2FC (q<0.05, Log2FC≥0.5, min.pct≥0.25), included subcluster 0: P2RY14, MAN1A1, FGF13, DSCAML1 and PID1; subcluster 1: FN1, TPRG1, KCNQ3, PDE3A and FMNL2; subcluster 2: TTN, STAB1, AHNAK, RBM25 and PRRC2C; subcluster 3: NOVA1, ADH1B, DLC1, FBN1 and NEGR1 (figure 2B; online supplemental figure 7B; online supplemental table 2B).

• Download figure
• Open in new tab
• Download powerpoint

Figure 2

Macrophage, adipocyte and endothelial cell subclusters in infrapatellar fat pad identified by single-nucleus RNA sequencing (n=21). (A, C, E) UMAPs and bar plots depicting cell type subclusters and their respective proportions. Top five DEGs within each cluster are indicated. (B, D, F) Violin plots and density plots depicting the expression of the top two unique transcriptomic markers within each macrophage, adipocyte and endothelial cell subcluster. Density plots present the expression density from dark purple (low) to yellow (high) (see online supplemental figures 7–13, online supplemental table 2B–D). DEGs, differentially expressed genes; UMAP, uniform manifold approximation and projection.

Subclustering of adipocytes identified four subclusters, with proportions ranging from 72.97% (subcluster 0) to 2.25% (subcluster 3) (total nuclei: 12 111, resolution 0.08) (figure 2C; online supplemental figure 9A). All adipocyte subclusters had unique transcriptomic markers, ranging from 160 genes (subcluster 2) to 89 genes (subcluster 0) (online supplemental figure 10; online supplemental table 2C). The top five genes for each subcluster based on decreasing log2FC (q<0.05, log2FC≥0.5, min.pct≥0.25), included subcluster 0: CLSTN2, WDPCP, PDE11A, ESR2 and LEPR; subcluster 1: CCDC200, SCD, TEX14, PDE4D and F13A1; subcluster 2: ITGBL1, ANO5, GPC6, STK32A and SEMA3C; subcluster 3: ATP1B3, CMSS1, PGAP1, PHLDB2 and CTH (figure 2D; online supplemental figure 9B; online supplemental table 2C). Pseudotime trajectory analysis of both adipocytes and fibroblasts together demonstrated a bidirectional differentiation pattern. Fibroblast subcluster 0 showed the greatest expression of universal fibroblast markers DPP4, PI16 and CD34, as described within fibroblast trajectory alone (figure 1G). This fibroblast population showed a trajectory towards adipocyte subcluster 1, through subclusters 3 and 0, terminating in subcluster 2, a more differentiated adipocyte population expressing PPARG and ADIPOQ (online supplemental figure 11). Fibroblast subcluster 3 was located between fibroblast subcluster 0 and all adipocyte subclusters within the trajectory suggesting that fibroblast subcluster 3 may contain a population of preadipocyte cells which further differentiates into specialised adipocyte populations (online supplemental figure 11). Additionally, fibroblasts expressing universal markers could also differentiate into more specialised fibroblasts, as demonstrated through CD55 and PDGFRA expression (online supplemental figure 11). This suggests that more differentiated adipocytes and fibroblasts can be derived from fibroblasts that express universal fibroblast markers.

Finally, we identified five endothelial cell subclusters within the IFP, with proportions ranging from 40% (subcluster 0) to 4.22% (subcluster 4) (total nuclei: 8944, resolution 0.05) (figure 2E; online supplemental figure 12A). Each endothelial cell subcluster had unique transcriptomic markers, ranging from 225 genes (subcluster 3) to 153 genes (subcluster 0) (online supplemental figure 13; online supplemental table 2D). The top five genes for each endothelial cell subcluster based on decreasing log2FC (q<0.05, log2FC≥0.5, min.pct≥0.25), included subcluster 0: CADM2, BTNL9, CD36, KIAA1217 and PPARG; subcluster 1: PLCXD3, TLL1, GPM6A, ZNF521 and MYRIP; subcluster 2: NR4A1, UBC, ZFP36, DEPP1 and ACTA2; subcluster 3: NKAIN2, SULF1, ARL15, PCSK5 and TOX; subcluster 4: NOVA1, ABI3BP, DCN, DCLK1 and BICC1 (figure 2F; online supplemental figure 12B; online supplemental table 2D).

Using previously published datasets, we determined that fibroblast subclusters within the IFP mapped to fibroblast populations identified from other adipose tissues (online supplemental figure 14A).7 13–16 In combined IFP and synovial tissue, four fibroblast subsets were uncovered,7 consistent with that determined using IFP only in our study (online supplemental figure 14A). Fibroblasts identified can also be referred to as preadipocytes, adipocyte stem and progenitor cells and mesenchymal stromal cells as they mapped to populations using this nomenclature.7 13–16 In addition, we mapped macrophage, adipocyte and endothelial cell subclusters we identified using reported white adipose tissue markers (online supplemental figure 14 B–D).13

Spatial profiling of major cell types, fibroblast subclusters and transcriptomic profiles of the KOA-IFP

After revealing major cell types, subtypes and their transcriptomic profiles within the IFP, we sought to resolve their spatial distribution using spatial transcriptomics (figure 3A). We analysed n=12 KOA-IFPs patient-matched to tissues analysed by snRNA-seq using Visium CytAssist spatial transcriptomics17 (online supplemental table 1). Spatial data were clustered and annotated independently from snRNA-seq data. We confirmed that fibroblasts, macrophages, adipocytes and endothelial cells were the major cell populations within IFPs, with fibroblasts distributed across the IFP alongside all other major cell types (Figure 3B; online supplemental figure 15A).

• Download figure
• Open in new tab
• Download powerpoint

Figure 3

Spatial profiling of major cell types identified by snRNA-seq data. (A) Schematic workflow to spatially resolve cell types within IFPs using spatial sequencing (created with BioRender.com). (B) UMAP (left) of spatial sequencing data of n=12 knee osteoarthritis (KOA)-IFP. Spatial resolution of cell populations identified using spatial sequencing across multiple IFPs (right). (C) UMAP (left) of spatial deconvolution of fibroblast subclusters across n=12 IFPs using snRNA-seq. Spatial resolution of deconvolved fibroblast subclusters within a representative IFP sample (middle). Density plots (right) displaying the module scores of the average gene expression of the top 10 transcriptomic markers of each deconvolved fibroblast subcluster. Module scores within density plots range from yellow (high) to dark purple (low). (D) Dot plots show the average expression of the top 10 transcriptomic markers within each fibroblast subcluster, from blue (low average expression) to red (high average expression). Spatial plots show the resolution of the module scores of coexpression of the top 10 transcriptomic markers of each fibroblast subcluster within a representative IFP, from blue (low module score) to yellow (high module score). (E) Spatial resolution of individual transcriptomic markers of each fibroblast subcluster include subcluster 0—FBN1 (DEG 1), subcluster 1—COL15A1 (DEG 10), subcluster 2—CRTAC1 (DEG 1), subcluster 3—FABP4 (DEG 2), subcluster 4—ENAH (DEG 5). (F) Bar plot showing the probability of each fibroblast subcluster and major cell types residing within a neighbourhood (see online supplemental figures 15,16). DEG, differentialy expressed gene; IFP, infrapatellar fat pad; snRNA-seq, single-nucleus RNA sequencing; UMAP, uniform manifold approximation and projection.

Fibroblasts from spatial transcriptomic data were deconvoluted using fibroblast subcluster annotations derived from snRNA-seq data. Each fibroblast subcluster was visualised in areas independently identified as fibroblasts by spatial transcriptomics (figure 3C). Subsequently, we spatially resolved the transcriptomic profiles of identified fibroblast subclusters through coexpression of the top 10 genes of each subcluster (figure 3D). We also visualised one of the top DEGs of each fibroblast subcluster including: subcluster 0: FBN1; subcluster 1: COL15A1; subcluster 2: CRTAC1; subcluster 3: FABP4; subcluster 4: ENAH (figure 1F; figure 3E; online supplemental figure 3C,15B).

Neighbourhood analysis using hoodscanR18 identified the probability of colocalisation of all cell types within the KOA-IFP. In five potential neighbourhoods, we found unique neighbourhoods containing all cell populations (figure 3F; online supplemental figure 16). Based on probabilities, adipocytes and endothelial cells were most likely to be found across all neighbourhoods (figure 3F). All fibroblast subclusters had a probability of being found across all neighbourhoods, with subcluster 1 likely being the dominant fibroblast subcluster within most neighbourhoods (figure 3F). Fibroblast subclusters 0 and 1 were likely the most dominant fibroblast subcluster within neighbourhood 4 (figure 3F). While we have identified five distinct neighbourhoods where all fibroblast subclusters reside, at this point, we are unable to specifically identify which anatomical features individual neighbourhoods are associated with. Overall, these data suggest that fibroblast subclusters are distributed across the IFP alongside adipocytes, endothelial cells and macrophages, with the potential to communicate.

Cell−cell communications within the KOA-IFP

Since fibroblasts had the largest proportion of nuclei, and all major cell types were spatially organised proximal to fibroblasts across IFPs, we investigated putative communications that occur between fibroblasts and all other major cell types within our snRNA-seq data. Using CellChat V.2.0,19 we identified major ligand-receptor cell signalling pathways, with macrophages, adipocytes and endothelial cells as ‘sender’ cells expressing ligands and fibroblasts as ‘receiver’ cells expressing receptors. Multiple cell communications to fibroblasts were identified, with the highest proportion of ligand-receptor pairs originating from adipocytes and the highest proportion of interactions occurring with fibroblast CD44 and ITGAV/ITGB8 receptors (figure 4A; online supplemental table 3).

• Download figure
• Open in new tab
• Download powerpoint

Figure 4

Putative cell–cell communication patterns between fibroblasts and major cell types. (A) Chord diagram (left) showing putative cell–cell communication between fibroblasts and other identified major cell types using snRNA-seq data from n=21 IFP samples. Fibroblasts are the ‘receiver’ cell type, expressing receptors and all other cell types are the ‘sender’ cells, expressing ligands. Individual chord diagrams (right) showing cell–cell communication between fibroblasts and other individual cell types. (B) River plot displays putative outgoing signalling patterns from each identified cell type. (C) Spatial plot showing putative collagen signalling between four major cell types along with another spatial plot displaying the expression of COL4A2-CD44 ligand-receptor pair in a representative sample. (D) Bubble plots displaying putative upregulated (left) and downregulated (right) signalling interactions between fibroblasts expressing receptors and other major cell types expressing ligands within KOA-IFPs compared to healthy control donor IFPs (see online supplemental table 3, online supplemental figures 17,18). IFP, infrapatellar fat pad; KOA, knee osteoarthritis; snRNA-seq, single nucleus RNA sequencing.

Additional investigations uncovered outgoing signalling patterns from all major cell types, without defining ‘receivers’ or ‘senders’. Regardless of directionality, fibroblasts, macrophages, adipocytes and endothelial cells had unique outgoing signalling patterns composed of multiple signalling pathways (figure 4B) suggesting intercommunication through multiple outgoing signals including collagen, PECAM1, PDGF and ANGPTL, which are known to stimulate regeneration and repair, extracellular matrix (ECM) deposition and angiogenesis.20–23

Since the collagen signalling pathway had the highest proportion of signalling going towards fibroblasts across all cell types, we spatially plotted the connectivity of this signalling network (figure 4C; online supplemental figure 17A). Furthermore, we were also able to spatially colocalise the COL4A2-CD44 ligand-receptor pair, which had one of the highest proportions of signalling within the collagen signalling pathway (figure 4A,C; online supplemental figure 17B).

We next analysed cell–cell communications within KOA-IFPs versus healthy control donor IFPs. Macrophages had a higher proportion of signalling towards fibroblasts within KOA-IFPs compared to healthy control donor IFPs (online supplemental figure 18A). Furthermore, KOA-IFPs had a higher proportion of upregulated communication interactions across all cell type pairings, compared to healthy control donor IFPs (figure 4D). The CD44 and ITGAV/ITGB8 receptors were involved in most upregulated and downregulated signals within KOA-IFPs and healthy control donor IFPs (Figure 4D). Regardless of directionality, all major cell types had some unique outgoing signalling patterns based on OA status, with a larger number of signalling interactions identified within KOA-IFPs (online supplemental figure 18B).

Transcriptomic differences in fibroblasts within the IFP based on KOA status

Since signalling to fibroblasts appeared modified in KOA-IFPs compared to healthy control donor IFPs, we investigated differences in cell composition and transcriptomic signatures based on KOA status. Fibroblasts, having the highest proportion of cells contributing to the IFP, were selected for subsequent analysis. Bioinformatic investigations of snRNA-seq data determined alterations in fibroblast subclusters within KOA-IFPs (n=15) compared to healthy control donor IFPs (n=6) (online supplemental table 1). We found all fibroblast subclusters within individual KOA-IFPs and healthy control donor IFPs (figure 5A). However, subcluster 3 had a higher proportion of nuclei contributed from healthy control donor IFPs compared to KOA-IFPs (q=0.032, multiple unpaired t-tests with false discovery rate (FDR) correction; figure 5B). Additionally, 38 DEGs (20 upregulated, 18 downregulated) were identified within KOA-IFPs versus healthy control donor IFPs, with the top five upregulated genes: FN1, ZNF385B, PRG4, ANKH and KAZN, and the top five downregulated genes: JUN, ABCA1, VIM, GSN and SAMD4A (figure 5C; online supplemental table 4A). The top five upregulated genes were highly expressed within subclusters 2 and 4 while the top five downregulated genes were highly expressed within subclusters 0, 1 and 3 (online supplemental figure 19). As downregulated and upregulated genes were localised to fibroblast clusters consistent with universal fibroblasts and the terminal fibroblast subclusters, respectively (figure 1G; online supplemental file 6), these results suggest that KOA-IFP fibroblasts have the expression profile of a more terminal differentiation state.

• Download figure
• Open in new tab
• Download powerpoint

Figure 5

Differences in IFP fibroblast subclusters based on OA status. (A) UMAPs of fibroblast subclusters within n=21 IFPs split by disease state (n=15 KOA-IFPs, n=6 healthy control donor IFPs). UMAPs of fibroblast subclusters within each of the n=21 IFPs. 15 IFPs (top) have KOA while 6 (bottom) are healthy control donor IFPs. (B) Stacked bar plot displaying the proportion of fibroblast subclusters within KOA-IFPs compared to healthy control donor IFPs. Fibroblast subcluster 3 had significantly more nuclei contributed from healthy control donor IFPs compared to KOA-IFPs (q=0.032). Fibroblast subcluster proportions were analysed by arcsin-transformation of the absolute proportions across all samples and performing multiple unpaired t-tests with false discovery rate (FDR) correction using the Benjamini, Krieger and Yekutieli two stage step up method. Top five DEGs within each cluster are indicated. (C) Volcano plot showing the log2 fold change (FC) of DEGs in KOA-IFPs compared to healthy control donor IFPs. DEGs were defined by a minimum of 50% of nuclei expressing a marker, with a log2FC>0.5, q<0.05 are upregulated while genes with a log2FC<−0.5, q<0.05 are downregulated. (D) Enriched PathDIP pathways associated with the upregulated DEGs from KOA-IFPs compared to healthy control donor IFPs. (E) Enriched PathDIP pathways associated with the downregulated DEGs from KOA-IFPs compared to healthy control donor IFPs. (F) Interaction network showing protein–protein interactions and transcription factor-gene interactions with enriched pathDIP pathways linked to DEGs (see online supplemental table 4A-C and online supplemental figure 19). DEGs, differentially expressed gene; IFP, infrapatellar fat pad; KOA, knee osteoarthritis; UMAP, uniform manifold approximation and projection.

We next used pathDIP analysis (pathDIP V.5, https://ophid.utoronto.ca/pathDIP)24 to identify enriched pathways associated with DEGs in KOA-IFPs compared to healthy control donor IFPs (figure 5D,E; online supplemental table 4B,C). Network analysis revealed connections between DEGs and enriched pathways from KOA-IFPs versus healthy control donor-IFPs, and putative interactions between transcription factors (Figure 5F). Of note, transcription factors SOX5 and CUX1 were upregulated in KOA-IFPs compared to healthy control donor IFPs while JUN and RBPJ were downregulated, with putative interactions between JUN and CUX1, SOX5 and JUN, and RBPJ and CUX1 (figure 5C,F). This suggests that transcription factor interactions may regulate fibroblast gene expression in KOA compared to healthy IFP.

Sex-based transcriptomic differences in fibroblasts within the KOA-IFP

To understand how sex influences fibroblasts within KOA-IFPs, we performed bioinformatic analyses of snRNA-seq data of female (n=8) compared to male (n=7) KOA-IFPs (online supplemental table 1). We observed no significant differences in the presence or proportions of fibroblast subclusters (q<0.05, multiple unpaired t-tests with FDR correction; figure 6A,B; online supplemental figure 20). However, we identified differences in transcriptomic profiles of fibroblasts within female versus male KOA-IFPs (figure 6C). We determined 105 DEGs: 35 upregulated and 70 downregulated. The top five upregulated genes include CLU, HTRA1, PRG4, TMEM196 and FN1 and the top five downregulated genes include UTY, USP9Y, LAMA2, AFF3 and NAV3 (figure 6C; online supplemental table 5A).

• Download figure
• Open in new tab
• Download powerpoint

Figure 6

Differences in KOA-IFP fibroblast subclusters based on sex. (A) UMAPs of fibroblast subclusters within n=15 KOA-IFPs split by sex (n=8 female KOA-IFPs, n=7 male KOA-IFPs). (B) Stacked bar plot displaying the proportion of fibroblast subclusters within female KOA-IFPs compared to male KOA-IFPs. No significant differences were found in proportions of nuclei contributed to subclusters from KOA-IFPs from female versus male study participants (q<0.05). Fibroblast subcluster proportions were analysed by arcsin-transformation of the absolute proportions across all samples and performing multiple unpaired t-tests with FDR correction using the Benjamini, Krieger and Yekutieli two stage step up method. Top five DEGs within each cluster are indicated. (C) Volcano plot showing the log2 fold change (FC) of DEGs in female KOA-IFP samples compared to male KOA-IFP samples. DEGs were defined by a minimum of 50% of nuclei expressing a marker, with a log2FC>0.5, q<0.05 are upregulated while genes with a log2FC<−0.5, q<0.05 are downregulated. (D) Gene Ontology (GO) biological processes and pathDIP pathways enriched in the upregulated genes of female compared to male KOA-IFPs. (E) GO biological processes and pathDIP pathways enriched for the downregulated genes within female compared to male KOA-IFPs. (F) Interaction network showing protein–protein interactions and transcription factor–gene interactions with enriched biological processes or pathDIP pathways linked to DEGs (see online supplemental figure 20 and online supplemental file 5A-E). DEGs, differentially expressed genes; FDR, false discovery rate; KOA, knee osteoarthritis; IFP, infrapatellar fat pad; UMAP, uniform manifold approximation and projection.

We used clusterProfiler25 and pathDIP24 to perform Gene Ontology biological process and pathway enrichment analysis, identifying processes and pathways associated with DEGs of female versus male KOA-IFPs (figure 6D,E; online supplemental table 5B–E). Network analysis determined connections between DEGs and associated biological processes and pathways, with putative interactions between transcription factors (figure 6F). Notably, transcription factors SOX5, FOXP2 and CREB5 were upregulated within female versus male KOA-IFPs, while ZEB1, BCL6, ZBTB16, EBF2 and RBPJ were downregulated, with putative interactions between transcription factor pairs ZBTB16 and BCL6, and ZEB1 and BCL6 (figure 6C,F), suggesting transcriptional gene regulation by transcription factor interactions in males.

Transcriptomic differences of fibroblasts within the KOA-IFP based on obesity status

To discern how obesity status influences fibroblasts of KOA-IFPs, we performed bioinformatic analysis of snRNA-seq data comparing fibroblast subclusters within obese (n=8; BMI 30–40 kg/m²) versus normal BMI KOA-IFPs (n=7; BMI 18.5–24.9 kg/m²) (online supplemental table 1). We found no significant differences in the presence or proportions of fibroblast subclusters between normal and obese BMI KOA-IFPs (q<0.05, multiple unpaired t-tests with FDR correction) (figure 7A,B; online supplemental figure 21). However, 21 DEGs were identified; 10 upregulated, 11 downregulated (figure 7C; online supplemental table 6A). Upregulated genes from obese versus normal BMI KOA-IFPs included: HMCN1, LHFPL6, ABCA6, PRRX1, ABCA8, ITGA11, DST, STARD9, SRRM2 and FBLN1, while downregulated genes included PRG4, CLU, FN1, ITGB8, FKBP5, UGP2, TIMP3, GPC6, ZBTB16, ITM2B and CREB5 (figure 7C; online supplemental table 6A).

• Download figure
• Open in new tab
• Download powerpoint

Figure 7

Differences in KOA-IFP fibroblast subclusters based on obesity status. (A) UMAPs of fibroblast subclusters within n=15 KOA-IFPs, split by obesity status (n=7 normal BMI KOA-IFPs (18.5–25 kg/m²), n=8 obese BMI KOA-IFPs (30–40 kg/m²)). (B) Stacked bar plot displaying the proportion of fibroblast subclusters within obese BMI KOA-IFPs compared to normal BMI KOA-IFPs. No significant differences were found in proportions of nuclei contributed to subclusters from KOA-IFPs from obese versus normal BMI study participants (q<0.05). Fibroblast subcluster proportions were analysed by arcsin-transformation of the absolute proportions across all samples and performing multiple unpaired t-tests with FDR correction using the Benjamini, Krieger and Yekutieli two stage step up method. Top five DEGs within each cluster are indicated. (C) Volcano plot showing the log2 fold change (FC) of genes that are uniquely differentially expressed in obese BMI compared to normal BMI KOA-IFPs. DEGs were defined by a minimum of 50% of nuclei expressing a marker, genes with a log2FC>0.5, q<0.05 are upregulated while genes with a log2FC<−0.5, q<0.05 are downregulated. (D) Dot plot (left) with average expression and percent of population expressing each DEG based on obesity status identified in snRNA-seq data (top) and within spatial sequencing data (bottom). Spatial resolution (right) of DEGs significantly upregulated within obese BMI compared to normal BMI KOA-IFPs in both snRNA-seq and spatial sequencing, visualised within spatial sequencing data. (E) GO biological processes and pathDIP pathways enriched for the upregulated genes from obese BMI compared to normal BMI KOA-IFP samples. (F) GO biological processes and pathDIP pathways enriched for downregulated genes from obese compared to normal BMI KOA-IFP samples. (G) Interaction network showing protein–protein interactions and transcription factor-gene interactions with enriched biological processes or pathDIP pathways linked to DEGs (see online supplemental figure 21, online supplemental table 6A-F). BMI, body mass index; DEGs, differentially expressed genes; ECM, extracellular matrix; FDR, false discovery rate; GO, gene ontology; KOA, knee osteoarthritis; IFP, infrapatellar fat pad; snRNA-seq, single-nucleus RNA sequencing; UMAP, uniform manifold approximation and projection.

We spatially resolved the 10 upregulated DEGs in obese compared to normal BMI KOA-IFPs using spatial transcriptomics (n=12; n=6 normal and n=6 obese BMI KOA-IFPs). All 10 upregulated genes identified in snRNA-seq data were also upregulated within spatial sequencing data based on log2FC (figure 7D; online supplemental table 6B). Interestingly, DST has a higher percent expression within snRNA-seq while FBLN1 has higher percentage expression within spatial data (figure 7D). Of the 10 upregulated genes, ITGA11, STARD9 and ABCA6 were significantly upregulated in spatial transcriptomic analyses (q<0.05; figure 7D; online supplemental table 6B).

We used clusterProfiler25 and pathDIP24 to perform Gene Ontology biological process and pathway enrichment analysis and identified enriched biological processes and pathways associated with DEGs in obese compared to normal BMI KOA-IFPs (figure 7E,F; online supplemental table 6C-F). Network analysis determined connections between DEGs and enriched biological processes and pathways, with putative interactions between transcription factors (figure 7G). Of note, the transcription factor PRRX1 was upregulated in obese versus normal BMI KOA-IFPs while ZBTB16 and CREB5 were downregulated, with a putative link between PRRX1 and ZBTB16 (figure 7C,G), suggesting transcription factor interactions may mediate transcriptional changes related to obesity status. We also identified enriched pathways related to metabolic functions, including lipid transport and carbohydrate metabolism (figure 7G).

Differences in metabolite levels of fibroblasts isolated from obese and normal BMI KOA-IFPs

Since metabolic dysregulation is a key element of obesity, and transcriptional differences in fibroblasts by obesity status had enriched pathways linked to metabolism, we sought to elucidate if alterations in supernatant metabolite levels of cultured fibroblasts from obese (n=5) compared to normal BMI KOA-IFPs (n=5) existed using targeted metabolomics (figure 8A). Overall, 445 metabolites were detected with 5 metabolites found at significantly different levels (p<0.05) from obese versus normal BMI KOA-IFPs, including upregulated: triglycerides (TG) (54:2) and TG (50:2), C18:2 (linoleic acid) and choline; downregulated: homoarginine (figure 8B; online supplemental table 7A). Using pathDIP,24 we found that differentially secreted metabolites from KOA-IFP fibroblasts by obesity status were annotated to metabolism-related pathways, including lipid metabolism, with enriched pathways associated with molecular transport (figure 8E; online supplemental figure 22).

• Download figure
• Open in new tab
• Download powerpoint

Figure 8

Differences in metabolite secretome in fibroblasts from KOA-IFPs based on obesity status. (A) Schematic workflow to identify differences in the metabolite levels within supernatants of fibroblast cultures of obese compared to normal BMI KOA-IFPs, stimulated with or without profibrotic or proinflammatory stimulus (created with BioRender.com). (B) Violin plot (left) and bar plot (right) showing the log2 fold change (FC) of metabolites differentially found in supernatants of fibroblast cultures from obese compared to normal BMI KOA-IFPs. Positive log2FC indicates an increase in the metabolite level within supernatants of cultures from obese compared to normal BMI KOA-IFP fibroblasts. (C) Volcano plot (left) showing the change in log2FC of metabolites differentially found in supernatants of obese compared to normal BMI KOA-IFP fibroblast cultures when treated with TGF-β compared to vehicle control. Positive values indicate a greater upregulation of metabolite level in obese versus normal BMI KOA-IFP fibroblast culture supernatants when treated with TGF-β versus vehicle control. Bar plot (right) showing the log2FC of metabolites differentially found in supernatants of obese compared to normal BMI KOA-IFP fibroblast cultures when treated with TGF-β compared to vehicle control, with positive values indicating an upregulation in levels found in supernatant. (D) Volcano plot (left) showing the change in Log2FC of metabolites differentially found in supernatants of obese compared to normal BMI KOA-IFP fibroblast cultures when treated with TNF-α compared to vehicle control. Positive values indicate a greater upregulation of a metabolite level in obese versus normal BMI KOA-IFP fibroblast culture supernatants when treated with TNF-α versus vehicle control. Bar plot (right) showing the log2FC of metabolites differentially found in supernatants of obese compared to normal BMI KOA-IFP fibroblast cultures treated with TNF-α compared to vehicle control, with positive values indicating an upregulation in secretion. Please refer to C for the legend. (E) Dot plot displaying pathways enriched for the differentially secreted metabolites by fibroblasts from obese compared to normal BMI KOA-IFPs and after treated with TGF-β versus vehicle control or TNF-α versus vehicle control (see online supplemental figure 22, online supplemental table 7A-C). BMI, body mass index; KOA, knee osteoarthritis; IFP, infrapatellar fat pad; TG, triglycerides; TGF-β, transforming growth factor β; TNF-α, tumour necrosis factor α.

Transforming growth factor β (TGF-β) can be found within KOA synovial fluid and is a key mediator involved in fibroblast activation and joint fibrosis.26 27 Furthermore, tumour necrosis factor α (TNF-α) is also found in KOA synovial fluid and is involved in promoting inflammatory responses in synovial fibroblasts.28 29 To mimic the fibrotic and inflammatory response in IFP fibroblasts, we stimulated cells with TGF-β or TNF-α, respectively. We then investigated how metabolite levels in culture supernatants of obese versus normal BMI KOA-IFP fibroblasts were modified when stimulated by TGF-β or TNF-α, compared to stimulation by vehicle control (phosphate-buffered saline). When comparing fibroblast cultures from obese versus normal BMI KOA-IFPs, nine metabolites were significantly differentially detected (p<0.05) in response to TGF-β treatment, compared to vehicle control (figure 8C; online supplemental table 7B). Mean levels of homocitrulline, ethanolamine, C16:1OH (2-hydroxypalmitoleylcarnitine), TG (52:5) and Hex2Cer (34:1) (di-hexosyl ceramide) were increased by TGF-β-treated fibroblasts from obese BMI KOA-IFPs versus vehicle control and decreased in normal BMI KOA-IFPs. In contrast, mean levels of cholesteryl esters (17:0) and (20:1), alpha-ketoglutaric acid and phosphatidylcholine diacyl (PCaa) C26:0 were increased in TGF-β-treated fibroblast cultures from normal BMI KOA-IFPs and decreased in obese BMI KOA-IFPs (figure 8C; online supplemental table 7B). Pathway enrichment analysis identified enriched metabolism and neuron-related pathways associated with the differentially secreted metabolite signature by BMI after profibrotic TGF-β stimulation (figure 8E).

Additionally, when comparing fibroblast cultures from obese versus normal BMI KOA-IFPs treated with TNF-α, levels of 10 metabolites were significantly differentially detected (p<0.05) compared to vehicle control (figure 8D; online supplemental table 7C). Mean levels of serotonin, 3-hydroxybutyric acid, 3-hydroxyisobutyric acid, glutaric acid, TGs (54:7) and (52:5) and N-acetyl-glycine were increased in TNF-α-treated fibroblast cultures from obese BMI KOA-IFPs and decreased in normal BMI KOA-IFPs (figure 8D; online supplemental table 7C). In contrast, mean levels of sphingomyelins C24:0, C16:0 and C16:1 were increased in TNF-α-treated normal BMI KOA-IFPs and decreased in obese BMI KOA-IFPs (figure 8D; online supplemental table 7C). Pathway enrichment analysis identified a single enriched pathway, cAMP signalling, associated with the differentially secreted metabolite signature from fibroblasts by obesity status with proinflammatory TNF-α stimulation (figure 8E). Overall, these data show differential metabolite levels in supernatants of fibroblasts from obese compared to normal BMI KOA-IFPs under normal conditions and in response to profibrotic and proinflammatory stimuli, with enrichment for metabolism and/or neuronal-related pathways.