Changes in the Serum Metabolome in an Inflammatory Model of Osteoarthritis in Rats

1. Introduction

Osteoarthritis (OA) is the most common form of arthritis, and one of the most prevalent diseases in middle-aged and older people. Knee osteoarthritis is one of the leading causes of physical disability in adults [1,2,3,4]. Initially, it was seen as a disease in which only a mechanical degradation of the cartilage occurred, but nowadays it is considered a very complex disease involving different tissues. Thus, alterations in the joint happen at different levels, namely, in the metabolism and architecture of the subchondral bone, in the morphology and metabolism of the articular cartilage, presenting periarticular osteophytosis, inflammation, and fibrosis of the synovial membrane in different degrees. In addition, it is related with changes in other tissues, such as ligaments, tendons and surrounding musculature [3,5,6,7,8,9].

The epidemiology of this disorder is complex and multifactorial with genetic, biological and biomechanical components. The main risk factors are age, obesity, abnormal mechanical joint loading and altered joint morphology [5,7,9,10,11,12].

Despite its worldwide importance, there is no standard treatment that cures, reverses, or slows down the development of the disease. This could be explained due to the poor understanding on its pathogenesis [1,10]. For this reason, treatment of OA has traditionally consisted on management of the primary cause, followed by treatment of pain, control of clinical signs, and surgical intervention in some late stages [1,5,12].

For OA diagnosis , a new tool called metabolomics has emerged, and it might provide more valuable information, given that current imaging techniques offer a late diagnosis lacking information on the functional adaptation of cartilage. Metabolomics consists in the study of small biological molecules in a system and holds great potential for early diagnosis, monitoring therapies and for the understanding of the pathogenesis of many diseases [1,3,4,5,13].

For this reason, it has become an ideal method for the identification of OA biomarkers in a variety of biological samples. Different studies have reported several metabolites and metabolic pathways that can be altered in OA, such as amino acid metabolism, fatty acid and lipids metabolism, phospholipids, arginine, phosphatidylcholine, L-tryptophan, tyrosine, carnitine, and arachidonic acid [1,2,4,9,14,15,16,17,18,19].

This current study has been performed through non-targeted metabolome and gene expression detection in samples obtained from rats using monoiodoacetate (MIA)-induced OA model. The hypothesis of the present study was that metabolomics could allow detection of changes in the serum metabolome in a patient with early osteoarthritis, which have not yet been described. The main objective was to evaluate serum metabolome and gene expression and to identify the main altered metabolites in a patient with osteoarthritis of inflammatory origin.

2. Results

2.1. Metabolomic Study

Blood samples were taken from all animals and analyzed by mass spectrometry at time before OA induction (T0), 28 days (T28), 56 days (T56) and 84 days (T84) after OA induction.

2.1.1. T0 vs T28 paired analysis

A paired analysis of serum metabolite results at T0 and T28, T56 and T84 was performed.

For the selection of significant variables (significant differences between T0 and T28) a Volcano Plot combining a Fold Change (FC) method with a t-test was performed With this data analysis , the aim was to obtain an overview and to select those potentially significant variables capable of discriminating between the conditions or categories of the study.

The paired analysis was performed using a script developed in the Analytical Unit with R software. Figure 1 shows the Volcano plot obtained. Note that both the FC and the p-value are in logarithmic scale (log10). Variables were selected (in red) with a threshold for FC=2 and a t-test p-value <0.05.

With the significant variables selected in the previous analysis, to detect the variables with the greatest discriminating power carried out Orthogonal Projections to Latent Structures Discriminant Analysis (OPLS-DA) multivariate analysis.

Figure 2 shows how the model discriminates between the two times T0 and T28 in the score plot obtained after this analysis.

The model diagnostics were adequate with R2Y=0.994 and Q2Y=0.967 and the model is further validated with a p-value coefficient of variation (CV)-ANOVA <0.001 and a permutation test of 1000 iterations (Figure 3).

The Variable importance in projection plot (VIP Plot) was performed from the model to select the most important discriminant variables.

From the VIP Plot, variables with low interaction coefficient (CI) (not including 0) and which were among the first 30 variables ranked by discriminant order (VIP score) were selected. Extracting each ion (m/z) in one of the quality control (QC) raw data and checking the peak shape and retention time also verified each variable

Table 1 shows the selected variables for identification using the Human Metabolome Database (HMDB) [20,21] and Metlin databases [22], as well as the mass spectrometry (MS/MS) analyzed in the equipment.

Eight organic acids and derivatives, benzenoids, organoheterocyclic compounds and lipid molecules variables were selected. Including lactacystin from the carboxylic acid class and derivatives; Taurine from the organic sulfonic acids class and derivatives; Styrene oxide and Tyramine from the benzene class and substituited derivatives; Setanaxib from the pyridines class and derivatives; Norsalsolinol from the tetrahydroisoquinolines class; Ganoderic acid V from the prenol lipids class and Ganglioside GM1 from the class sphingolipids (Table 2).

2.1.2. T0 vs T56 paired analysis

For the selection of significant variables (significant differences between T0 and T56) a Volcano Plot combining a Fold Change method with a t-test was constructed. With this analysis of the data, the aim was to obtain an overview and to select those potentially significant variables capable of discriminating between the conditions or categories of the study.

The paired analysis was carried out using a script developed in the Analytical Unit with R software. Figure 4 shows the Volcano plot obtained. Note that both the FC and the p-value are in logarithmic scale (log10). Variables were selected (in red) with a threshold for FC=2 and a t-test p-value <0.05.

In this case, the significant variables were only 5 variables, which were selected to build the OPLS-DA model, whose Score plot is presented in Figure 5.

As can be seen in the graph, samples 1, 4, 8 and 11 seem to form a cluster between them, which indicates that there are other variables not considered, apart from time, that may be influencing the results (intra-group variability).

Despite this, the VIP plot is constructed to arrange the variables by discriminant power selecting and identifying the most important ones. Out of 6 variables, one of them was discarded due to a high confidence interval of the VIP score, and another one because the chromatographic peak was not verified with the retention time.

The Table 3 shows an abridgementof the results obtained for the 4 variables finally obtained.

The identification results of the selected variables were three organic acids and derivatives, including gamma-Glutamylcysteine, Tyrosyl-Serine and Acetyl citrate from the carboxylic acid class and derivatives (Table 4).

2.1.3. T0 vs T84 paired analysis

For selection of significant variables (significant differences between T0 and T84) a Volcano Plot combining a Fold Change method with a t-test was performed. With this analysis of the data, the aim was to obtain an overview selecting those variables potentially significant capable of discriminating between conditions or categories of the study.

The paired analysis was carried out using a script developed in the Analytical Unit with R software. Figure 6 shows the Volcano plot obtained. Note that both the FC and the p-value are in logarithmic scale (log10). Variables were selected (in red) with a threshold for FC=2 and a t-test p-value <0.05.

With the significant variables selected in the previous analysis, to detect the variables with the greatest discriminating power, a OPLS-DA multivariate analysis was carried out.

Figure 7 shows how the model discriminates between the two times T0 and T84 in the score plot obtained after this analysis. The sample of animal 8 (T84) must be considered as it differs from the rest of the samples in that same group.

The model diagnostics were adequate with R2Y=0.946 and Q2Y=0.925 and the model is further validated with a p-value CV-ANOVA <0.001 and a permutation test of 1000 iterations (Figure 8).

The Variable importance in projection Plot (VIP Plot) was constructed from the model to select the most important discriminant variables.

From the VIP Plot, those variables with low CI (not including 0) and which were among the first 30 variables ranked by discriminant order (VIP score) were selected. Extracting each ion (m/z) in one of the QC raw data and checking the peak shape and retention time also verified each variable.

Table 5 presents the variables finally selected for identification using the The Human Metabolome Database (HMDB) [20,21] and Metlin databases [22], as well as the mass spectrometry (MS/MS) analyzed in the equipment.

The results following identification of the selected variables were seven lipid molecules, including: Ginsenoside Rh1 and Theasapogenol A from the prenol lipids class; Phosphatidic acid from the glycerophospholipids class; Polyporusterone F, Brassinolides, Ursodeoxycholic acid from the steroids class and steroid derivatives; and 10-Hydroperoxy-H4-neuroprostane related to prostaglandins and from the fatty acyls class (Table 6).

From the results obtained, can be concluded that there are differences between the times T0 and T84.

In the discriminant analysis observed intra-group variability, which indicates that there are other variables not taken in account in the study, in addition to time that are influencing the distribution of metabolites.

2.2. Histologic study

As expected, all rats developed degenerative and inflammatory changes associated with OA following MIA injection. This occurred at every time point, after 28, 56 and 84 days (Figure 9).

Histologic sings of osteoarthritis were observed microscopically in all OA induced samples. The main alterations presented were the reduction of the cartilage stain intensity, and the presence of an irregular cell density along the samples (Figure 10).

3. Discussion

OA induction in rats by intraarticular MIA produces changes in serum metabolome, as observed in the results of this study. Furthermore, these differences are associated with the time elapsed since the induction of OA.

Being the pathophysiology of OA complex remaining unknown now a days, it was decided to investigate metabolomics as an advanced diagnostic method to see whether it could provide further insights of the disease. This technique allows us to observe the progress at a biochemical level and monitor the evolution of treatments [4,9,10,13]. In addition, as the study period was relatively short, we needed a sampling technique that would allow us to make an early diagnosis compared to other techniques such as radiography or MRI where diagnosis is made at more advanced stages of the disease[1,5,9,23,24]. The effectivenessof the OA model was verified by histological examinations.

Metabolomics is supported by several studies in which significant differences in serum metabolites were observed between healthy and OA rats, as observed in this study [1,2,9,10,13,14]. In the study by Chen et al. in rats with OA, plasma samples were analyzed by metabolomics and it was shown that density labelling mass spectrometry, the same technique used in this study, had a high sensitivity for detecting metabolites in rat plasma [1].

On the other hand, one of the factors to consider was the use of serum to analyze the metabolites present in the animal. The study carried out by Zhang et al. compared blood and synovial fluid (SF) observing how the range of metabolites varied considerably and out of 168 only 8 were consistently related. This study suggests that metabolic changes are joint specific and other inflammatory processes may influence the concentration of metabolites in serum [4], as observed atthe study performed by Guma et al. which states that there is a fairly modest correlation between plasma and SF [13].

Thus, as discussed in several articles, it is necessary to know which metabolites are altered for a better understanding of the joint status [3,13,27,28]. Regarding the metabolites found in this study, they differ depending on if the sample has been collected at T28, T56 and T84. Different benzenoids, organoheterocyclic compounds, organic acids, and lipid molecules were detected at T28, whereas only organic acids were observed at T56. At T84 mainly lipid molecules where detected.

The differential benzenoids classes found in this study were benzene and substituted derivatives. The organoheterocyclic compounds classes found in this study were tetrahydroisoquinolines and pyridines and derivatives. Previous studies had not observed any association of these metabolites with OA.

On the other hand, differences in organic acids and derivatives were detected at T0 vs T28 and T0 vs T56. The main organic acids classes found in this study were carboxylic acid and derivatives (T28 and T56), and organic sulfonic acids and derivatives (T28). The organic sulfonic acid detected in the analysis at T0 vs T28 was taurine. This result correlates with other studies where taurine metabolism was found to be one of the main metabolic pathways involved in OA, as taurine is related to the pathophysiology of OA, correlating with subchondral bone sclerosis and playing a vital regulatory role [15,23,24,29,30,31]. Anderson et al. observed that elevated taurine in OA processes could indicate an increased subchondral bone sclerosis [31]. Yang et al. showed that taurine levels in sclerotic subchondral bone was positively regulated [24]. Taken together, these studies in synovial fluid revealed that altered taurine metabolism in subchondral bone has a direct correlation with subchondral bone sclerosis in osteoarthritis. In this study, taurine was analyzed from blood samples, so this metabolite could be used as an early biomarker of subchondral bone sclerosis in OA, but further studies are required to confirm this hypothesis.

Within the organic acid group, we have also detected carboxylic acids at T28 and T56, just as Swank et al. detected hippurate carboxylic acid in urine after 18 months of OA progression [15]. Within this group, we found that at T56 the metabolite acetyl citrate was observed. Regarding to this metabolite, different studies have detected the presence of citrate in urine and synovial fluid, related to OA. Citrate is an intermediate in the tricarboxylic acid cycle and its increase indicates an OA-related alteration in the cycle [29,31,32,33,34]. Therefore, it would be interesting to investigate further the presence of acetyl citrate in blood as a biomarker of OA.

Lipid molecules were found at T0 vs T28 and T0 vs T84. These possible biomarkers were related to the findings of several studies, where there is an alteration of lipid metabolism associated with OA, due to its pro-inflammatory properties [16,17,18,35,36]. On the other hand, although in synovial fluid, Kosinska et al. in several studies detect alterations at the level of phospholipids and sphingolipids at different stages of the disease [37,38,39], thus, understanding the relationship between OA and lipid molecules analysis may be helpfull in future treatments [8].

The main lipid classes in which differences were found in this study are sphingolipids at T28, prenol lipids at T28 and T84, and glycerophospholipids, steroids and fatty acyls at T84. These results correlate with prior reports in which the same lipid biomarkers have been detected, such as a study by Pousinis et al. which describes the presence of glycerophospholipids, sphingolipids and fatty acyls in plasma from OA rat model at 112 days [36].

Regarding prenol lipids, these are a group of natural molecules, and their most biologically relevant classes are fat-soluble vitamins (vitamins E, A and K) [40]. Neogi et al. observed an association between low plasma levels of vitamin K and increased prevalence of OA manifestations in hand and knee, because vitamin K supports calcium homeostasis, facilitating bone mineralization among other effects [41,42,43]. Regarding vitamin E, different studies have shown its potent anti-inflammatory properties, as well as in the prevention and regulation of the progression of age-related diseases [44,45]. Therefore, further research on the relationship between OA and prenol lipids would be needed.

On the other hand, sphingolipids detected at T28 had been previously detected in other studies and were found to be related to subchondral bone sclerosis in OA playing an important regulatory role in the pathological process of sclerotic subchondral bone[18,24,36,38]. Tootsi et al. found changes in serum sphingolipid levels in humans with OA, confirming their involvement in the pathogenesis of OA [18]. Kosinska et al. found that sphingolipids could alter synovial inflammation and repair responses in damaged joints [38]. Thus, it would be interesting to use sphingolipids as blood biomarkers of OA.

Phospholipids are molecules associated with inflammation and increased cartilage damage at synovial fluid level and may be associated with the pathogenesis of OA [37]. The glycerophospholipids form the essential lipid bilayer of all biological membranes, and changes in glycerophospholipid concentrations and composition are associated with OA development , as shown in multitude studies [4,16,18,36,46]. Therefore, changes in the concentration of lipid molecules, more specifically glycerophospholipids, may indicate risk of OA.

The biomarker Ursodeoxycholic acid detected in this study belongs to the steroids class and super class lipids and lipid-like molecules (HMDB). Ursodeoxycholic acid is a naturally occurring dihydroxy hydrophilic bile acid, where Moon et al. demonstrated that this bile acid has a preventive potential as treatment in a model of induced OA by reducing pain and ameliorating cartilage destruction [47]. On the other hand, Carlson et al. detected metabolites from steroid hormone biosynthesis in synovial fluid in humans with rheumatoid arthritis [48]. No studies have been found on detection of alterations in ursodeoxycholic acid at the metabolomic level in animals with OA, so this metabolite should be considered in future investigatigations.

4-Hydroperoxy-H4-neuroprostane, also known as 14-H4-NeuroP, is a member of the class of compounds known as prostaglandins, and related compounds, of the class fatty acyls and super class lipids and lipid-like molecules (HMDB). Similarly, in the study by Zhao et al., it was observed that serum levels prostaglandin estradiol2 were significantly increased in the OA group. In addition, several metabolites of the class fatty acyls and super class lipids such as aminobutyric acid, stearic acid or L-carnitine were increased [9]. Attur et al. examined plasma lipids prostaglandins E2 (PGE2) and found PGE2 elevated in symptomatic knee OA patients [49]. Similarly, Gierman et al. associated changes in PGE2 levels with the development of OA [50]. On the other hand, the study by Shi et al., and Pausinis et al., also shows changes in arachidonic acid or linoleic acid, metabolites within the same classification [2,36]. Rregarding acylcarnitines, of the fatty acyls class, several studies have shown changes of their concentration in serum of animal models with OA [46,51]. Thus, changes in prostaglandin concentrations or fatty acyls could indicate presence of OA.

There are some limitations to this study. Firstly, the sample size was small and did not allow strong validation of these potential biomarkers. Therefore, a larger sample size would be necessary in future research. Secondly, only blood samples were used, and in the future, it would be interesting to correlate serum with synovial fluid measurements for a better understanding of these observed metabolic changes. Thirdly, the possible relationship between OA and arthritic diseases and whether these biomarkers are useful or not for identifying other forms of arthritis. Fourthly, variables other than time were not considered and there is intra-group variability in these results, for this reason other variables will need to be considered in the future. Finally, it should be considered that the identified metabolites in this study should be directly evaluated in subsequent targeted studies, aiming to confirm or rule out their role in the modification of serum metabolome in an inflammatory model of osteoarthritis.

4. Materials and Methods

4.1. Experimental model:

4.1.1. Experimental design:

A prospective, experimental, randomized, and double-blinded study was designed. The study was conducted at the Hospital Universitari i Politècnic La Fe, within the animal facility of the Instituto de Investigación Sanitaria La Fe (IISLaFe), Valencia, Spain. This experimental study was approved by the Ethics and Animal Welfare Committee of the Hospital Universitari i Politècnic La Fe and authorized by Valencian Government with 2017/VSC/PEA/00177 type 2 code, in accordance with the provisions of article 31 of Royal Decree 51/2013.

4.1.2. Experimental trial:

To carry out the study, thirty ten-weeks-old Wistar rats of similar age and body condition entered the study. Animals were distributed in three groups depending on their survival time (28, 56 and 84 days). Out of 10 animals in each group , metabolomic study was performed in 6 animals and histological analysis in 4 animals.

Osteoarthritis was induced in all 30 subjects by intra-articular infiltration of 0.4 mg of monoiodoacetate (Sodium iodoacetate^®, Sigma Aldrich) into the right knee. All animals were sedated with buprenorphine (0.03mg/kg) (Buprex^®, Indivior, city, country), ketamine (65 mg/kg) (Ketolar^®, Pfizer, city, country) and medetomidine (0.01 mg/kg) (Sedator^®, Dechra, city, country) intraperitoneally.

Once each subgroup reached their survival times, rats were euthanised. Thereafter, all right and left stifles were photographed for macroscopic analysis. Subsequently, the knees of 4 aleatory subjects from each subgroup were assigned for histological study, and the other 6 for metabolome serum study. Blood samples were obtained of each group at day 0 (before MIA infiltration) and just before euthanasia (28, 56 and 84 days after MIA infiltration).

4.2. Obtaining the metabolomic results:

4.2.1. Sample preparation:

Serum sample preparation was performed following the protocols established in the Analytical Unit, as detailed below. First, 50 μL of serum, plus 150 μL of acetonitrile (ACN) 0.1% formic acid (FA) (cold) vortexed for 30 min at -20ºC, were collected. The sample was centrifuged for 10min at 4ºC and 13000g; the supernatant extract was collected in an eppendorf tube and stored at -80ºC. Subsequently, 20μL of the extract was collected and placed in a 96-well plate for Liquid chromatography - Quadrupole time-of-flight-6550 (LC-QTOF-6550); 100 μL of FM (Agμa 0.1% FA) + 10 μL of MIX internal standard (ISTD) (20μM) was added. Once the plate was prepared, 5 μL of each sample was taken and the QC was prepared. The reagent blank was prepared using the same blood collection tube used in the water study and following the same preparation procedure as the serum samples (looking for artefacts in the tube, reagents and other material).

4.2.2. Sample analysis:

Figure 11 depicts the protocoled that was followed to perform the analysis of the samples.

To avoid intra-baTch variability, as well as to ensure the quality and reproducibility of the analysis, we proceeded as follows. First, a random injection order. Subsequently, an analysis of at least 5 quality control condition (QC_cond) at the beginning of the sequence, to condition column and equipment (these data was not used in the multivariate analysis of the data). Finally, an analysis of a quality control pool (QC_pool) every 5-7 samples.

4.2.3. Ultra-performance liquid chromatography, time-of-flight, mass spectrometry (UPLC-ToF-MS) method

For analysis, liquid chromatography equipment was used coupled to a time-of-flight mass spectrometer. Standard procedures of the Analytical Unit laid down the chromatographic and mass spectrometric conditions, summarized below: mode positive and negative electrospray ionisation (ESI); range m/z: 100-1700 Da; UPLC column: Acquity UPLC BEH C18 (100 x 2.1 mm, 1.7 μm); injection volume (V_inj): 5 μL; column temperature: 45°C; autosampler temperature 4°C; flow rate: 500 μL/min; mobile phase A=H2O (0.1%v/v HCOOH); mobile phase B=CH3CN (0.1% v/v HCOOH)

4.3. Data analysis of metabolomic results:

4.3.1. Pre-processing of the metabolome data:

Before performing the multivariate analysis of the data, pre-processing of the acquired data was required. This pre-processing consists on a series of processes such as filtering, molecular feature detection, peak alignment and clustering, and data normalization. Based on the results obtained in this study, several parameters were selected for each treatment. R statistical software and the XCMS library were used to carry out the data pre-processing.

At the end of the processing, we obtained two tables of "molecular features", with all the variables extracted and normalized, Table in negative mode (9712 molecular features: m/z and Rt), and Table in positive mode (4045 molecular features: m/z and Rt). Data analysis performed from these tables.

4.3.2. Analysis of the quality of the metabolome results:

4.3.2.1. Evaluation of the response of internal standards

To detect possible problems in the injection or in the preparation of the samples, the variability of the internal standards as well as intra-batch variability was assessed for each sequence. The control chart of reserpine and leucine-enkephalin (leuEnk) used as internal standard (STDI) for the positive ESI mode is attached in Figure 12 as an example.

No tendency in the response of these patterns is observed throughout the sequence, thus assuming low intra-batch variability and correct analysis.

4.3.2.2. QC evaluation

Before proceeding to the QC evaluation, the data were adjusted for the area of the internal standard LeuEnk in positive mode and reserpine in negative mode.

For the positive and negative mode sequences, the coefficients of variation (CV %) of the QCpool are calculated and those variables where the CV≥ 30% are removed. In this way, analytical variability is eliminated, and the remaining variables are considered to come from possible biological/metabolic variability.

4.4. Histological evaluation

Following sacrifice, left and right stifles were dissected. A craniolateral approach of the skin and dissection of the soft tissues around the femur and tibia were performed. A 1.5 cm osteotomy proximal to the femoral trochlea and distal to the tibial plateau was done to retrieve the stifles. Once the stifle was isolated, periarticular soft tissues were meticulously dissected to retrieve the biological samples. After removing all the stifle soft tissue, direct visualization of the joint was performed to score the samples following the macroscopic scale described by Laverty et al. [52].

Once anatomical samples were retrieved and after macroscopical evaluation had been performed, femoral condyles samples were extracted and fixed in formaldehyde 4%. Following fixation and decalcification with ethylenediaminetetraacetic acid (EDTA) (Osteodec^®, LABOLAN, Navarra, Spain), paraffin inclusion and 4µm longitudinal section cuts using a microtome was performed. Sections were obtained from three different anatomical zones: lateral condyle, femoral trochlea, and medial condyle. Samples were stained using hematoxylin & eosin and Masson’s Trichrome stain. Following staining, slides were digitalized for evaluation using a specific slide viewer software (CaseViewer 2.2^®, 3DHISTECH Ltd, Budapest, Hungary).

Lastly, microscopic evaluation was done using the Osteoarthritis Research Society International (OARSI) semi-quantitative scale described by Laverty et al. to evaluate matrix stain, cartilage structure, chondrocyte density and cluster formation [52]. On top of this, researchers added an additional parameter resulting from the addition of the results of all measured variable as a total score altogether. Structure of subchondral bone was evaluated using semi-quantitative scales from OARSI described by Gerwin et al. [53].

4.4.1. Analysis of the Histological results.

A descriptive analysis was done for the histologic variables. The data available were obtained from semiquantitative scales, and are presented numerically as median, minimum to maximum, and they are graphically presented as median, inter-quartile range, and minimum and maximum.

5. Conclusions

Based on the results obtained in this study, it can be concluded that metabolomics is a promising tool for a better understanding the pathogenesis of OA, as MIA-induced osteoarthritis resulted in changes in the serum metabolome of rats. Furthermore, it can be observed how the metabolites that may be involved in these changes differ according to the time elapsed since the induction of osteoarthritis. Eighteen potential biomarkers were identified, predominantly from the lipid molecules, organic acids, benzenoids and organoheterocyclic compounds classes. The lack of data on the functions of most of these metabolites helps to focus future studies and highlights the need for further clinical studies in knee OA.