Circulating microRNAs differentiate fast-progressing from slow-progressing and non-progressing knee osteoarthritis in the Osteoarthritis Initiative cohort.

Introduction: The objective of this study is to identify circulating microRNAs that distinguish fast-progressing radiographic knee osteoarthritis (OA) in the Osteoarthritis Initiative cohort by applying microRNA-sequencing.
Methods: Participants with Kellgren-Lawrence (KL) grade 0/1 at baseline were included (N = 106). Fast-progressors were defined by an increase to KL 3/4 by 4-year follow-up (N = 20), whereas slow-progressors showed an increase to KL 2/3/4 only at 8-year follow-up (N = 35). Non-progressors remained at KL 0/1 by 8-year follow-up (N = 51). MicroRNA-sequencing was performed on plasma collected at baseline and 4-year follow-up from the same participants. Negative binomial models were fitted to identify differentially expressed (DE) microRNAs. Penalized logistic regression (PLR) analyses were performed to select combinations of DE microRNAs that distinguished fast-progressors. Area under the receiver operating characteristic curves (AUC) were constructed to evaluate predictive ability.
Results: DE analyses revealed 48 microRNAs at baseline and 2 microRNAs at 4-year follow-up [false discovery rate (FDR) < 0.05] comparing fast-progressors with both slow-progressors and non-progressors. Among these were hsa-miR-320b, hsa-miR-320c, hsa-miR-320d, and hsa-miR-320e, which were predicted to target gene families, including members of the 14-3-3 gene family, involved in signal transduction. PLR models included miR-320 members as top predictors of fast-progressors and yielded AUC ranging from 82.6 to 91.9, representing good accuracy.
Conclusion: The miR-320 family is associated with fast-progressing radiographic knee OA and merits further investigation as potential biomarkers and mechanistic drivers of knee OA.

Introduction

Characterizing disease phenotypes remains an outstanding challenge for osteoarthritis (OA), particularly for early stages of disease and for rapid disease progression. Although biochemical (e.g. collagen breakdown products) and imaging (e.g. magnetic resonance imaging) biomarkers have been explored for OA, these markers typically become robust in established disease. The ability to identify patients with early-stage and fast-progressing OA would support investigation of disease mechanisms, delivery of preventive interventions, and recruitment to clinical trials, and therefore is essential to advance the field. Unfortunately, there are currently no validated methods for detecting early-stage or fast-progressing knee OA, but there is promise in exploring circulating microRNAs.

MicroRNAs have emerged as powerful candidate biomarkers for musculoskeletal diseases (e.g. osteoporosis)[2,3] and are known to play important mechanistic roles in OA.[4,5] These small, stable, noncoding RNA molecules can be comprehensively profiled in liquid biopsies using microRNA-sequencing. Sequencing offers advantages with respect to sensitive and specific detection of microRNAs, especially compared with previous techniques such as microarrays which are limited to semi-quantification of preselected microRNAs. Recent sequencing studies in blood identified differentially expressed microRNAs in OA versus controls and in early-stage versus late-stage radiographic knee OA.[7-10] Among these, a well-powered study revealed a signature of seven circulating microRNAs to be associated with early-stage radiographic knee OA. However, sequencing remains to be applied to large-scale longitudinal cohorts to identify microRNAs that are associated with disease progression.

The purpose of this study is to apply our method for microRNA-sequencing to identify circulating microRNAs that are associated with radiographic knee OA progression. We hypothesize that unique microRNAs in plasma are associated with fast-progressing radiographic knee OA as compared with slow-progressing and non-progressing radiographic knee OA. To test this hypothesis, we leverage the Osteoarthritis Initiative (OAI) which is a well-characterized knee OA cohort with imaging, biospecimens, and clinical data available from baseline to 10-year follow-up, with ongoing collection. We apply a data-driven strategy to define groups that advance from early to late stages of radiographic knee OA, then perform microRNA-sequencing at both baseline and 4-year follow-up in the same participants to identify circulating microRNAs that are differentially expressed over time (Figure 1). These microRNAs may be clinically relevant as potential biomarkers and/or mechanistic drivers of disease.

Figure 1.

Schematic overview of study design. Participants selected from the Osteoarthritis Initiative (OAI) cohort were characterized by the rate of radiographic knee osteoarthritis progression and plasma samples collected at baseline and 4-year follow-up were subjected to microRNA (miRNA)-sequencing. Created with BioRender.com.

Methods

OAI cohort

The OAI includes participants with varying disease trajectories that can be defined based on established features such as Kellgren–Lawrence (KL) radiographic grades for knee OA severity. From among the 1926 participants in the OAI with KL grades 0/1 at baseline, we used literature recommendations to guide our selection of a sample size of at least 20 participants per group which can reliably detect a 1.5-fold change in microRNA expression with a power of 80%. Based on individual KL-grade trajectories over 8 years, fast-progressors were defined as participants who increased to KL 3/4 by 4-year follow-up (N = 20). Slow-progressors were defined by KL 0/1 at baseline and 4-year follow-up, with KL 2/3/4 by 8-year follow-up (N = 35). Non-progressors remained at KL 0/1 throughout the 8-year follow-up and N = 51 were selected using propensity score methods. The OAI cohort is described in detail at https://nda.nih.gov/oai and has been approved by the Institutional Review Board for the University of California, San Francisco (FWA approval # 00000068; IRB approval # 10-00532), and its affiliates. The following OAI datasets were used: XRay[00-10], AllClinical[00-10], Kxr_sq_bu[00-10], Enrollees, Biomarkers[00-10], and Outcomes99. All participants provided written informed consent to the OAI.

MicroRNA-sequencing

Frozen aliquots of 200 µL plasma collected at baseline (N = 106) and 4-year follow-up (N = 102; 4 samples were not available from the OAI) were obtained from the OAI and stored at −80°C until use. After thawing on ice, RNA was isolated using the miRNeasy Serum/Plasma Advanced kit (QIAGEN, Hilden, Germany) and microRNA libraries were created with the QIAseq miRNA Library kit (QIAGEN, Hilden, Germany). This protocol is described in detail in our previous publication. Single-end 76-base sequencing on the Illumina NextSeq550 platform was performed to an average depth of 11.6 million (± 2.6 million standard deviation) reads per sample. Sequencing data were deposited in the Gene Expression Omnibus database under accession number GSE183188.

Bioinformatics

Quality control of sequencing data (e.g. Q score ⩾ 30), reference-based alignment, counts generation, and novel microRNA discovery were performed as previously described.[6,8] No samples were excluded from the present analysis due to insufficient data quality. For microRNAs shortlisted by differential expression (DE) analyses, predicted target genes were retrieved from mirDIP version 4.1 (http://ophid.utoronto.ca/mirDIP). Pathway enrichment analysis was performed for each set of targets using pathDIP 4 API in R 4.0.3. From the pathways common among at least three of four lists, the gene-pathway links were further investigated. Protein interactions for the gene targets were retrieved from Integrative Interactions Database (IID), version 2021-05 (http://ophid.utoronto.ca/iid). The network was visualized using NAViGaTOR, version 3.0.16, and exported to Adobe Illustrator, version 26.0.2, to be finalized with legends.

Statistical analysis

DE analyses were performed at baseline and 4-year follow-up using negative binomial models cross-sectionally and a negative binomial mixed model longitudinally.[19,20] A propensity score constructed using age, sex, body mass index (BMI), Western Ontario and McMaster Universities Osteoarthritis Index (WOMAC) and arthroscopy status (whether the participant underwent arthroscopy prior or during the 8-year follow-up) was introduced in the models as a covariate; additionally, the cross-sectional models were adjusted for sequencing batch. MicroRNAs with false discovery rate (FDR) < 0.05 were considered in penalized logistic regression (PLR) models for each of the studied pairs (fast- versus slow-progressors, fast- versus non-progressors, and slow- versus non-progressors) to determine their predictive performance. PLR hyperparameters were tuned using threefold 50-repeated cross-validation with Cohen’s kappa as a performance measure under the ‘one-SE’ rule for model selection. Resulting microRNA-based models were benchmarked against clinicodemographics-only models using area under the receiver operating characteristic curve (AUC).

Results

Characterizing participants by rate of radiographic knee OA progression

Using a data-driven phenotyping approach, we characterized the rate of radiographic knee OA progression in 106 participants and identified three distinct groups (Figure 2(a)). By 4-year follow-up, fast-progressors exhibited KL 3/4 (N = 20), whereas slow-progressors remained at KL 0/1 and only exhibited KL 2/3/4 by 8-year follow-up (N = 35). Non-progressors remained at KL 0/1 through to 8-year follow-up (N = 51). Exploring clinicodemographic factors across groups, differences were found in worsening total WOMAC scores from baseline to 4-year follow-up and arthroscopy status (Table 1). More slow-progressors (42.9%) showed worsening WOMAC score compared with non-progressors (17.6%), despite both groups being KL 0/1 at 4-year follow-up. Furthermore, more fast-progressors and slow-progressors underwent arthroscopy as compared with non-progressors, suggesting that our definition of non-progressors was internally consistent.

Figure 2.

Circulating microRNAs are differentially expressed among fast-progressors in the OAI cohort. (a) Kellgren–Lawrence (KL) grade (XRKL) is plotted on the y-axis and time (years) is plotted on the x-axis. Dashed lines connect individual trajectories; solid lines represent fitted second-order polynomials describing XRKL as a function of time. Fast-progressors increase from KL grade 0/1 at baseline to KL 3/4 by 4-year follow-up. Slow-progressors remain at KL 0/1 from baseline to 4-year follow-up then increase to KL 2/3/4 by 8-year follow-up. Non-progressors remain at KL 0/1 through to 8-year follow-up. (b, c) Volcano plots showing changes in microRNA expression resulting from cross-sectional analysis of fast-progressors versus slow-progressors, fast-progressors versus non-progressors, and slow-progressors versus non-progressors at baseline (b) and 4-year follow-up (c). Red points indicate microRNAs showing statistically significant changes with a log2 fold change (FC) of ±1 or more. (d, e) Venn diagrams showing unique and overlapping differentially expressed microRNAs among each of the three comparisons performed at baseline (d) and 4-year follow-up (e).

NS, not significant.

Table 1.

Clinicodemographic characteristics of progression groups.

	Non	Slow	Fast	p-value
N	51	35	20
Age, mean years (SD)	59.2 (9.2)	57.5 (8.4)	61.0 (8.1)	0.374
Age, N (% > 60 years)	20 (39.2)	10 (28.6)	11 (55.0)	0.153
Sex, N (% males)	13 (25.5)	9 (25.7)	7 (35.0)	0.696
Race, N (% non-white)	8 (15.7)	2 (5.7)	4 (20.0)	0.247
BMI, mean kg/m2 (SD)	28.8 (4.1)	29.2 (5.1)	28.5 (5.1)	0.855
BMI categories:				0.513
Healthy weight = 18.5–24.9 kg/m2, N (%)	10 (19.6)	7 (20.0)	5 (25.0)
Overweight = 25–29.9 kg/m2, N (%)	26 (51.0)	12 (34.3)	9 (45.0)
Obese = 30 + kg/m2, N (%)	15 (29.4)	16 (45.7)	6 (30.0)
WOMAC worsening a , N (% worsen)	9 (17.6)	15 (42.9)	6 (30.0)	0.038*
Arthroscopy, N (% yes)	3 (5.9)	12 (34.3)	6 (30.0)	0.002*
Propensity score, mean (SD)	0.43 (0.15)	0.64 (0.22)	0.55 (0.25)	<0.001*

Each factor was summarized across progression groups by their means (standard deviations) or counts (proportions) for continuous or categorical factors, respectively. Differences among groups were ascertained using Welch (continuous) and chi-square (categorical) tests.

BMI, body mass index; SD, standard deviation; WOMAC, Western Ontario and McMaster Universities Osteoarthritis Index.

Total WOMAC score with worsening defined by an increase of 7 or more points from baseline to 4-year follow-up.

Statistically significant.

Identifying microRNAs associated with fast-progressing knee OA

Following microRNA-sequencing, we identified 532 and 411 microRNAs at baseline and 4-year follow-up with greater than 10 counts-per-million in at least two samples among 106 and 102 samples, respectively. Cross-sectional DE analyses comparing the three groups pairwise within each time-point identified 78 microRNAs at baseline (Figure 2(b)) and 13 microRNAs at 4-year follow-up (Figure 2(c)) at FDR < 0.05 (Supplementary File 1). Among these, 48 microRNAs at baseline (Figure 2(d)) and 2 microRNAs at 4-year follow-up (Figure 2(e)) were upregulated when comparing fast-progressors versus slow-progressors and non-progressors, but not slow-progressors versus non-progressors. The two microRNAs at 4-year follow-up were hsa-miR-320c and hsa-miR-320d. Notably, hsa-miR-320d was also among the 48 microRNAs identified at baseline, along with hsa-miR-320b and hsa-miR-320e (Figure 3(a)). Therefore, hsa-miR-320d was significantly increased among fast-progressors as compared with both slow-progressors and non-progressors at both baseline and 4-year follow-up. Longitudinal DE analyses identified seven microRNAs between fast-progressors versus non-progressors (FDR < 0.05) and no significant differences between fast-progressors versus slow-progressors or slow-progressors versus non-progressors (Supplementary File 1). Our analysis for novel microRNAs identified four sequences at baseline and two at 4-year follow-up that were in  > 50% of fast-progressors, though the sequences were unique to each time-point (Supplementary Figure 1).

Figure 3.

MiR-320 family members are distinguishing factors of fast-progressing knee OA. (a) Log2 fold changes for members of the miR-320 family identified as differentially expressed at FDR < 0.05 in at least one comparison at baseline and 4-year follow-up. (b) ROC curves comparing models using baseline clinicodemographic variables alone (orange) including age (continuous), sex, BMI (continuous), and WOMAC (continuous) against microRNA-based models (green) from each comparison.

NS, not significant.

Exploring the miR-320 family in fast-progressing knee OA

To determine the extent to which the DE microRNAs could distinguish fast-progressors from slow-progressors and non-progressors, we built and tested a PLR model per group comparison. Using the comparison-specific DE microRNAs as input, we found that AUCs for the microRNA-based models were greater than those based on baseline age, sex, BMI, and WOMAC (Figure 3(b)). Of note, hsa-miR-320b/c/d/e were all among the top 20 predictors for the fast-progressor versus non-progressor model (Supplementary File 1). We then sought to predict the potential function of the miR-320 family. Bioinformatic analyses revealed that these miR-320 members had both unique and overlapping putative target genes and pathways, with 137 genes (Figure 4(a); Supplementary File 2) and 14 pathways (Figure 4(b); Supplementary File 3) common to miR-320b/c/d. Among these, four of the seven members of the 14-3-3 gene family comprising YWHAZ, YWHAQ, YWHAH, YWAHE, YWHAB, YWHAG, and SFN were identified to be annotated with eight common pathways each (Figure 4(c); Supplementary File 4). We further investigated the connection between the miR-320 family and the 14-3-3 gene family in mirDIP and found that all four members of the miR-320 family were predicted to target all seven members of the 14-3-3 gene family, albeit some predictions had lower confidence (Figure 4(d)). In addition, we found multiple protein–protein interactions connecting all seven members of the 14-3-3 gene family (Figure 4(d)).

Figure 4.

MiR-320 family members are potential mechanistic players in fast-progressing knee OA. (a) Unique and overlapping mirDIP gene targets for four members of the miR-320 family, filtered for top 1%. (b) Unique and overlapping pathways for four members of the miR-320 family, filtered for p < 0.01. (c) Number of common enriched pathways (x-axis) with which the gene targets (y-axis) are annotated. (d) Predicted targeting of the 14-3-3 gene family by the miR-320 family using all mirDIP results. Edges with arrows denote the score class with V = very high, H = high, M = medium, and L = low according to the legend. Blue edges denote protein–protein interactions (PPIs) from IID.

Discussion

Leveraging the OAI cohort, we used a data-driven approach to characterize three distinct disease progression groups that provide greater phenotypic resolution than dividing participants into radiographic progressors versus non-progressors alone. Applying microRNA-sequencing to samples from 106 participants at 2 time-points, we identified circulating microRNAs that are associated with fast-progressing radiographic knee OA over time, including members of the miR-320 family, which are predicted to target members of the 14-3-3 gene family among other genes. While future studies are required to externally validate these findings in independent longitudinal knee OA cohorts, our data suggest microRNAs may be key molecular markers and mechanistic drivers of fast-progressing knee OA.

Imaging and biochemical markers have been explored as prognostic biomarkers in the OAI cohort, yet are not currently in clinical use.[23,24] As these approaches require structural changes to have occurred before markers can be detected, opportunities for prevention are reduced. As epigenetic regulators of gene expression, microRNAs may appear earlier in disease pathogenesis and therefore could precede irreversible structural damage in knee OA. This is a strength of microRNAs as biomarkers, and our data suggest unique combinations of microRNAs offer strong predictive ability. Through use of the OAI cohort, our microRNA-sequencing data can be integrated with other imaging and biochemical data to support discovery of composite biomarkers for OA.

To note limitations of this study, we characterized OA based on radiographic severity (i.e. KL grade) alone. Future studies may explore other structural (e.g. joint space narrowing) or symptomatic (e.g. pain) features that typify knee OA. We also selected a subset of the OAI cohort enriched for participants who progressed from early-stage radiographic knee OA at baseline (KL 0/1) to late-stage radiographic knee OA by 8-year follow-up (KL 3/4). Future studies with larger sample sizes could investigate broader and alternative definitions of fast-progressing knee OA (e.g. KL 2 to KL 4 in 4 years) while accounting for confounders such as comorbidities.

We report four of the five members of the miR-320 family to be associated with fast-progressing radiographic knee OA, a family of microRNAs our predictions show to have both unique and overlapping gene targets. Previous studies suggest miR-320a and miR-320c protect against OA cartilage degeneration,[26-28] potentially through downregulation of the Wnt signaling pathway.[27,29] In addition, plasma miR-320a along with miR-98-5p was found to distinguish between neuropathic and nociceptive musculoskeletal pain. Our data indicate miR-320d is upregulated in fast-progressors versus both slow-progressors and non-progressors at both baseline and 4-year follow-up, suggesting it is associated with fast-progressing knee OA over time and may have mechanistic functions through regulation of its targets, including the 14-3-3 gene family.

The 14-3-3 gene family comprises seven conserved protein isoforms which function to regulate signal transduction. Members of this family have been explored as biomarkers and mechanistic drivers of rheumatoid arthritis (RA) and OA. Of the seven isoforms, YWHAH and YWHAG were detected in synovial fluid from inflammatory joint disease patients, and YWHAH was elevated in the serum compared with healthy controls.[34,35] Furthermore, YWHAH was found to correlate with other biomarkers of RA, and like YWHAE could induce inflammatory factors and signaling pathways linked to joint damage.[36,37] Given the observed upregulation of miR-320 family members in fast-progressing radiographic knee OA, members of the 14-3-3 gene family would likely be downregulated in a direct regulatory relationship. Taken together, our data suggest nuanced interactions among individual miR-320 and 14-3-3 family members may be a mechanism for fine-tuning signal transduction, a hypothesis that merits further investigation for its potential role in fast-progressing knee OA.

Click here for additional data file.

Supplemental material, sj-docx-1-tab-10.1177_1759720X221082917 for Circulating microRNAs differentiate fast-progressing from slow-progressing and non-progressing knee osteoarthritis in the Osteoarthritis Initiative cohort by Shabana Amanda Ali, Osvaldo Espin-Garcia, Andy K. Wong, Pratibha Potla, Chiara Pastrello, Madison McIntyre, Starlee Lively, Igor Jurisica, Rajiv Gandhi and Mohit Kapoor in Therapeutic Advances in Musculoskeletal Disease

Click here for additional data file.

Supplemental material, sj-png-1-tab-10.1177_1759720X221082917 for Circulating microRNAs differentiate fast-progressing from slow-progressing and non-progressing knee osteoarthritis in the Osteoarthritis Initiative cohort by Shabana Amanda Ali, Osvaldo Espin-Garcia, Andy K. Wong, Pratibha Potla, Chiara Pastrello, Madison McIntyre, Starlee Lively, Igor Jurisica, Rajiv Gandhi and Mohit Kapoor in Therapeutic Advances in Musculoskeletal Disease

Click here for additional data file.

Supplemental material, sj-txt-1-tab-10.1177_1759720X221082917 for Circulating microRNAs differentiate fast-progressing from slow-progressing and non-progressing knee osteoarthritis in the Osteoarthritis Initiative cohort by Shabana Amanda Ali, Osvaldo Espin-Garcia, Andy K. Wong, Pratibha Potla, Chiara Pastrello, Madison McIntyre, Starlee Lively, Igor Jurisica, Rajiv Gandhi and Mohit Kapoor in Therapeutic Advances in Musculoskeletal Disease

Click here for additional data file.

Supplemental material, sj-txt-2-tab-10.1177_1759720X221082917 for Circulating microRNAs differentiate fast-progressing from slow-progressing and non-progressing knee osteoarthritis in the Osteoarthritis Initiative cohort by Shabana Amanda Ali, Osvaldo Espin-Garcia, Andy K. Wong, Pratibha Potla, Chiara Pastrello, Madison McIntyre, Starlee Lively, Igor Jurisica, Rajiv Gandhi and Mohit Kapoor in Therapeutic Advances in Musculoskeletal Disease

Click here for additional data file.

Supplemental material, sj-txt-3-tab-10.1177_1759720X221082917 for Circulating microRNAs differentiate fast-progressing from slow-progressing and non-progressing knee osteoarthritis in the Osteoarthritis Initiative cohort by Shabana Amanda Ali, Osvaldo Espin-Garcia, Andy K. Wong, Pratibha Potla, Chiara Pastrello, Madison McIntyre, Starlee Lively, Igor Jurisica, Rajiv Gandhi and Mohit Kapoor in Therapeutic Advances in Musculoskeletal Disease

Click here for additional data file.

Supplemental material, sj-xlsx-1-tab-10.1177_1759720X221082917 for Circulating microRNAs differentiate fast-progressing from slow-progressing and non-progressing knee osteoarthritis in the Osteoarthritis Initiative cohort by Shabana Amanda Ali, Osvaldo Espin-Garcia, Andy K. Wong, Pratibha Potla, Chiara Pastrello, Madison McIntyre, Starlee Lively, Igor Jurisica, Rajiv Gandhi and Mohit Kapoor in Therapeutic Advances in Musculoskeletal Disease