Osteoarthritis (OA) is a chronic disease characterized by persistent joint pain and limited mobility, especially for those who experience from OA of the knee [1]. In the United States alone, the overall economic burden of OA is estimated at almost US $140 billion annually, highlighting the substantial societal and personal costs associated with the condition [2]. Pain is the most important outcome for patients with OA, with moderate to severe pain associated with greater functional limitations, treatment dissatisfaction, and reduced quality of life [3]. In managing OA, patient education and exercise therapy, which are the recommended first-line treatments according to clinical guidelines, play a critical role [4-6]. Treatment programs that combine therapeutic exercise with other treatments such as weight loss [7], manual therapy, and pain education have shown promising results in managing mild to moderate OA symptoms and provide complementary options to improve patient’s pain, function, and quality of life [8,9]. Patient education further empowers individuals by enhancing their understanding of the disease and emphasizing the importance of lifestyle adjustments such as a balanced diet and regular exercise [10].

However, the implementation of recommended treatments could be improved [11]. The discrepancy between recommended and provided treatments has been attributed to several factors, including health care providers offering no evidence-based care and patients’ lack of motivation and awareness about the long-term benefits of the recommended therapies [12,13]. Furthermore, some patients are hesitant to participate in exercise therapy programs, either because they are uninterested or face logistical difficulties [14,15]. It is relevant to develop strategies that assist clinicians and patients in an engaging way in the treatment decision process such as offering personalized outcome predictions.

Digital technology presents a promising solution to address such challenges in health care. Prediction models, particularly, have emerged as powerful tools to support patients and health care professionals in the treatment decisions and management of chronic diseases like OA [16-18].

The GLA:D (Good Life with osteoArthritis in Denmark) program has garnered significant attention as a pioneering initiative in OA management [19]. Comprising three components—training for physiotherapists, a patient education and neuromuscular exercise therapy program delivered to patients, and clinical data collection—the program stands out as an example of evidence-based care [19]. To support shared decision-making, the acceptance, and participation rate of programs like GLA:D, we aimed to validate an existing model and introduce an updated concise personalized prediction model that estimates changes in knee pain intensity for patients with OA considering participation in the GLA:D program [20].

This paper follows the guidelines of TRIPOD (Transparent Reporting of a Multivariable Prediction Model for Individual Prognosis or Diagnosis) [21].

In this paper, we reproduced and validated the results of previous work [20] with some methodological changes as described below.

We used data from the Danish GLA:D initiative for patients with knee and hip OA. The initiative consists of a 2-day course for physiotherapists, a patient treatment program delivered in clinical practice, and a registry of data reported by patients and clinicians. The patient program combines 2 patient education and 12 supervised exercise therapy sessions aiming at improving symptoms, physical function, and quality of life. Patients can join the GLA:D program through referrals from health care professionals (eg, general practitioners or orthopedic surgeons) or via self-referral, followed by an assessment by a certified GLA:D clinician to confirm the OA diagnosis [22-25]. By 2023, in total, 257 clinics were offering the GLA:D program for knee or hip OA [26]. Of these, most clinics were private practices, but 20 clinics were public municipalities, and 2 municipalities offered the treatment for employees.

The GLA:D registry collects data from participants, which includes demographic details, medical history, pain intensity, and physical function measures [27]. The data are collected at three time points: at baseline, immediately following the program, and at 12 months follow-up. GLA:D was established in 2013 and is regularly updated to include new evidence. Comprehensive details about GLA:D’s education, neuromuscular exercise program, and general information are previously available [19]. For this paper, we followed the inclusion criteria of a previous publication [20] but with an extended inclusion period from October 9, 2014, to November 12, 2022, instead of ending on August 31, 2017. We included patients who indicated their knee as the primary joint of complaint and provided complete data. We excluded participants for specific time periods due to technical problems in the registry, meaning one of our variables of interest was not collected. This was the case between May 13, 2016, and November 12, 2016, and April 10, 2018, and August 05, 2019.

To develop and evaluate our predictive model, we began by splitting our dataset into two subsets: a training set for model development and a test set for validation. The training process was consistently applied on all models, regardless of training on all variables, on the top k predictor variables, or on a concise set of variables. To enhance the robustness and generalizability of our models, we integrated 10-fold cross-validation with random forest regression. Using the test sets, we validated our findings, ensuring that our results were not biased or limited to specific data instances by averaging results across the cross-validation test sets.

To evaluate the performance of the models, we used the root-mean-squared error (RMSE) and R² metrics. RMSE measures the average deviation between the predicted and actual values, providing an assessment of the model’s accuracy. The R² represents the proportion of variance the model explains and indicates how well the model fits the data. To compare our personalized outcome predictions with predictions of the average model, we evaluated the number of correct predictions, in contrast to the previous publication, which used the absolute mean differences in a test dataset [20].

The predicted changes in pain intensity by our random forest regression constitute the personalized predictions. The average model prediction was defined as a mean function of the VAS pain change score of the training samples, μ=∑n=1NVASChangeScoreN. The average model always provides μ, regardless of the features of the unseen samples, as the predicted change in pain intensity. In this paper, we are always focused on the VAS pain change score, rather than the pain score at a specific time point. Recognizing that predictions inherently have some margin of error, we incorporated in our evaluation a margin of error that corresponds to a clinically relevant threshold.

In is the indicator value that will be one for cases where their predicted change in pain value is inside the interval. For all other samples, it is zero. Ln=VASChangeScore− R and Un=VASChangeScore + R show the lower and upper bound based on the clinical relevance threshold R. Finally, we calculate the percentage of predictions inside the interval ρ to evaluate and compare our model against the average model:

Since the value of a clinical relevance R is debatable and ranges up to 20 mm [33-35] points, we evaluated the number of correct predictions with different tolerable margins of error, ranging from 5 to 20 points. Thus, we calculated the number of correct predictions for several cases, the first allowing a deviation of 5 points and the last allowing a deviation of 20 points from the true change in pain.

All data analyses were performed in Python (Python Software Foundation) using Scikit-learn, Pandas, and NumPy libraries [36]. We used Pandas and NumPy for data preprocessing tasks, such as cleaning and variable scaling. Scikit-learn was used for validation techniques, including train-test split, cross-validation, and stratified sampling. We used standard hyperparameters for the random forest regression: n_estimators set to 100, max_depth limited to 10, and random_state fixed at 42. We chose to use them to provide a baseline performance that can be compared with other standard implementations.

To evaluate if the inclusion of recently added variables to the GLA:D registry improved our outcome predictions, we additionally developed and evaluated models following the same process as described above, but including 12 additional variables (load-related pain, reduced functional capacity, morning stiffness, crepitus, reduced knee movement, bony enlargement, previous joint injury, occupational or recreational overuse, family members with OA, Knee Injury and Osteoarthritis Outcome Score-12 pain, functioning, and Knee Injury and Osteoarthritis Outcome Score-12 summary score, which were added to the registry in May 2018). A flowchart for this additional analysis is provided in Figure S2 in Multimedia Appendix 1 [37].

This study did not require ethics approval or clinical trial registration, as determined by the local ethics committee of the North Denmark Region, since it was a register-based study rather than a clinical trial. The GLA:D registry has been approved by the Danish Data Protection Agency, ensuring compliance with data protection regulations. In accordance with the Danish Data Protection Act, patient consent was not required, as personal data were processed solely for research and statistical purposes. To protect participant confidentiality, all data were anonymized before analysis, with no access to personally identifiable information. Data handling followed the regulations of the Danish Data Protection Agency. No compensation was provided to participants, as this study relied on secondary data analysis of the GLA:D registry.

A total of 38,547 patients fulfilled our inclusion criteria for participating in the GLA:D program between October 9, 2014, and November 12, 2022. We considered patients with knee pain as the primary complaint and our analytical dataset was reduced to 10,216 patients. Details on the process of patient inclusion and exclusion are available in Figure 2.

The Gini impurity scores for all variables are presented in Figure S2 in Multimedia Appendix 1 [37]. The characteristics of all included and excluded samples are provided in Table 1. Based on these findings, we applied the elbow method and included the top 11 variables (Figure 3).

In total, we built 3 models each including a different set of variables (Table 2). The first includes all 34 variables (full model). The second one includes the 11 most important variables identified with the Gini impurity and elbow method (continuous model), highlighting the importance of these baseline variables in making predictions. The third includes the 6 predictive variables, which were chosen based on their importance and clinical reasoning—it is assumed that answering the underlying 11 questions would be applicable in routine clinical practice. These 6 variables are baseline pain, duration of symptoms, the EQ-5D score, time to complete a 40-meter walking test, age, and BMI. Each of those 6 variables contains 1 item but BMI and EQ-5D include 2 and 5 items, respectively.

The model performances are presented in Table 3 based on the average performance on the test sets in a 10-fold cross-validation process. Accordingly, the full model’s average RMSE was 18.56 (SD 0.59), while the coefficient of determination (R²) was 0.32. The continuous model displayed an RMSE of 18.64 (SD 0.58) and an R² of 0.32. Moreover, when considering variables in the concise model, the RMSE and R² were 18.76 (SD 0.64) and 0.31, respectively.

The comparison between true predictions based on our concise model including the top six most important variables and predictions using the average model is illustrated in Figure 4. It shows a small but distinct advantage of our personalized predictions over average-based predictions. Our personalized prediction method, based on random forest regression, performed better than the average prediction model by a margin of about 7% when allowing a difference of 15 points in each direction between the predicted and the true change on the test data [35]. Specifically, our results showed that the average model predicts 51.47% and our personalized model predicts about 57.82% of the cases correctly. We found that personalized predictions were slightly and consistently better than those based on the average improvement, regardless of the deviation allowed (5 points or 20 points). Table 3 shows correct personalized predictions and correct average predictions within the interval of ±15 points for all models. An illustration using these numerical values is as follows: if a patient, for instance, has a pain score of 45 before GLA:D and our personalized model estimates a change of 20 points, this means that it is estimated that the patient will have pain between 10 and 40 after the program with a 58% certainty.

Including the 12 additional variables reduced our sample with complete cases to n=1458, since these variables were first added later. However, the evaluation of the variable importance followed a similar pattern highlighting that the now 14 continuous variables were the most important (Figure S2 in Multimedia Appendix 1 [37]). Similarly, the model performances and evaluations remained comparable among all three models (Table 4).

In this study, we developed personalized prediction models for changes in knee OA pain after supervised patient education and exercise therapy (GLA:D) and compared them to existing prediction and average-based models. Our findings validate the two previously developed models from Baumbach et al [20] and introduce a new concise personalized prediction model for estimating changes in knee pain intensity among patients with knee OA considering participating in the GLA:D program. Neither the increase in sample size nor the inclusion of additional variables improved the previous prediction model. Nonetheless, our concise model correctly predicted 58% of the cases and demonstrated a 7% improvement in the rate of correct predictions over the use of currently used average values for informing patients about their expected changes in knee pain.

In clinical decision-making, predictive models serve as pivotal tools across various disciplines. These models, akin to those used in cancer prognosis or the optimization of exercise regimens, underline the significance of personalized predictions in enhancing patient outcomes [38-42]. This study delves into this paradigm by comparing prediction models and validating the findings of an earlier study. As in the previous study [20], we used the GLA:D data and presented a full and continuous model to predict personalized outcomes. The latter was chosen based on the variable importance. While we used the Gini impurity and elbow method, the previous study determined the most important variables by the reduction in RMSE for out-of-bag cases and applying the elbow method [20]. Despite the difference in the used method, our findings that the continuous variables were most important align with the previous study. Our outcomes match those of the previous study, despite our model’s sample size being doubled and the variables being reduced by 17 and 4 in the full and continuous model, respectively. For the full and continuous models, there was no difference, and a difference of 0.01 in the R² values, respectively. In addition, the RMSE aligned closely. In our supplementary analysis, we incorporated 12 new variables, included in the GLA:D registry between 2017 and 2018, which reduced our sample size to 1458 but maintained consistent variable importance and model performance. Thus, we found that neither an increase in sample size nor the added variables could improve the personalized outcome predictions compared to the earlier study [20].

Recent research [43,44] highlights that patients using analgesics, experiencing constant pain, or preferring the GLA:D program over saline injection demonstrated greater benefits from the program. Unfortunately, our data did not include “constant pain” or “preference for the GLA:D program.” However, we considered “intake of analgesics,” although it was not among the top predictive variables. The random forest method inherently models nonlinear relationships and interactions, enabling relevant subgroup distinctions within the data. This approach avoids reducing sample sizes and maintains robustness by incorporating all relevant variables in the model [45]. Including the variables “constant pain” and “preference for the GLA:D program” in future analyses could provide deeper insights into personalized outcomes and patient stratification. Additionally, it is important to note that our sample probably differs from the randomized controlled trial–based group in the study by Henriksen et al [43,44], as not all patients in our real-world cohort would have agreed to participate in a controlled trial.

To compare the personalized predictions with predictions based on average values, we again used a different approach than the previous prediction study. The previous study took the absolute mean difference and concluded that the difference between their model and the average improvements was not clinically relevant. In this study, we focused on the number of correct predictions, assessing a range of clinical relevance thresholds to determine whether our personalized model surpasses the average model’s performance, despite changing the clinical relevance threshold. Here, we emphasize this incremental progress of increasing the clinical relevance threshold. As a result, allowing a higher deviation from the true change resulted in a higher number of correct predictions, for both our personalized predictions and the average prediction. Overall, we found that our personalized prediction model was able to predict about 7% more cases correctly compared to using average predictions, independent of the threshold allowed from the true change.

However, the question of whether these findings are sufficient for implementation into clinical practice remains open. Existing shared decision-making tools often concentrate on comparing different treatments and assisting patients in making informed choices [18]. Specifically, in the context of OA, where education and exercises are universally recommended and surgery is considered a last resort after the failure of previous treatments, the potential of improved prediction models to alter this recommendation landscape is significant [46]. If future models could reliably predict the likelihood of initial treatment failures for specific patient groups, clinical guidelines to treat OA could evolve. Future models could benefit from including more comprehensive patient characteristics, such as constant pain and treatment preferences [43,44], to further explore patient-specific responses to the GLA:D program. Yet, it is crucial to acknowledge that the current performance of our models is too preliminary for such decisive changes in clinical practice. This underscores the necessity for further research to enhance model accuracy and reliability, as well as to understand the needs and preferences of patients and clinicians, potentially reshaping treatment protocols for OA based on predictive insights.

This study also has some limitations. Selection bias is possible in any registry-based study due to the cases lost to follow-up. However, we only excluded 28.44% due to missing information on the outcome at follow-up, suggesting no serious threats to its external validity [47]. Furthermore, from previous studies using the GLA:D data, we know that these patients do not differ to a clinically relevant level from the included patients [20]. Two periods of technical flaws in data collection were identified, and these missing values, unrelated to patient characteristics or variables, meet the criteria for missing completely at random. Participation in the GLA:D program may also be influenced by factors contributing to selection bias, such as limited accessibility, financial constraints, awareness and referral issues, perceived severity of OA, and language barriers [22,25,48-52]. Our sample’s greater representation of female participants reflects known trends of greater knee OA prevalence and greater participation in exercise programs like GLA:D [53-55]. While this enhances relevance to real-world treatment populations, it should be considered when generalizing findings. This study suggests that the increased sample size does not significantly affect the predictive performance outcome [20]. However, it should be noted that we included fewer variables in comparison to the previous study [20]. These variables were deemed of less relevance and therefore excluded from the registry, nonetheless, they contributed to the predictions of the previous models. Therefore, to be precise, we can only conclude that the increased sample size was able to compensate for the reduced number of variables. Hence, to improve our models’ predictability, new variables would probably need to be integrated. The newly integrated variables (12 additional variables) in our supplementary analyses, led to a reduction in the sample size (n=1458), however, the model performance remained similar. Consequently, we can (only) conclude that the newly added variables could compensate for the reduced sample size.

Based on these observations, we can point toward the limitation of available predictor variables. Implementing more detailed questionnaire-based measures for anxiety and depression might enhance pain predictions [56]. Weight loss is another recommended first-line treatment for individuals with overweight or obesity with knee OA [57] that is known to moderate systemic inflammation [58], extending the questionnaire to enquire about patients’ dietary habits might thus be another possibility to increase the predictive value [59]. However, due to the clinical applicability, the extent of the GLA:D registry questionnaire needs to balance the collection of information against patient burden, adding additional variables with minor importance is not desired. In conclusion, we currently do not know which variables could have an additional major positive influence on the personalized outcome predictions over those already included [60].

The precision of our model aligns with existing clinical studies, suggesting a need for future investigations into the practical value of such predictive accuracy within clinical settings [61]. Current initiatives like the A.S.K. report and Movement is Life have begun to explore the integration of predictive models and patient-reported outcomes to guide shared decision-making, demonstrating the utility of these tools in clinical practice [18,62]. However, there is a noticeable gap in the literature concerning the specific performance benchmarks required for predictive models in physiotherapy to be deemed clinically relevant, particularly in the management of OA. Therefore, before a decision on implementing our model as a web tool into clinical practice can be drawn, guidance on which metrics and predictive values need to be achieved in a model to be feasible for clinical practice needs to be investigated. Furthermore, it remains uncertain if the implementation of a model, like ours, which provides more correct predictions than what is currently used, but still several wrong predictions, is valuable enough for patients and clinicians to be implemented in clinical practice. Further, to the best of our knowledge, it is yet unknown if personalized outcome predictions regarding a treatment lead to different expectations than average outcome predictions in patients. Moreover, on ethical aspects, the consequences of wrong predictions in general on a patient’s treatment, and mental and emotional well-being should be investigated. For patients with estimated worsening pain predicting other health outcomes additionally could be essential to motivate them to participate in the program. For example, prediction tools within GLA:D for secondary outcomes such as physical activity, which changes independently of changes in pain, could support participation rates, given the overarching health benefits associated with increased physical activity [27,60]. Similarly, patients might be more motivated if they knew about the expected changes in physical function and quality of life. Finally, any prediction tool needs to be comprehensively evaluated based on its role in treatment decisions, patient satisfaction, and overall treatment outcomes in clinical settings, before implementation.

The methodology used in this study, which leverages machine learning techniques to analyze questionnaire data and functional test data for predicting treatment outcomes has been used before and is a promising advancement in the field of personalized medicine [63,64]. This approach is characterized by the selection of predictive variables from data and can be used in similar questioning. Thus, our method extends beyond the specific context of the GLA:D program and could be applied to other exercise therapy programs targeting knee OA, and potentially, for a broader range of conditions [62,65-67].

Generalizing our findings, the model predictions of the GLA:D program for knee OA to other exercise and physiotherapy interventions are challenging. The unique combination of group-based exercises and education in GLA:D may be crucial and not present in all physiotherapeutic settings. Therefore, our outcome predictions specific to GLA:D cannot be expected universally. However, if the core principle of structured exercise therapy and education are incorporated into the treatment of a patient with OA, our model might be applicable. Nonetheless, the personalized predictions would likely deviate even further from the true changes than in patients participating at GLA:D, therefore, they should not be generalized.

We developed and validated a personalized prediction model for pain intensity changes in patients with knee OA participating in the GLA:D program. Our model’s 58% accuracy is 7% better than using average pain improvement as a prediction. This likely represents the best possible prediction given current knowledge and available data. Although performance is improved, further guidance is needed to determine when such models are suitable for clinical implementation. This study highlights both the potential and current limitations of personalized prediction in knee OA management.