UEF//eRepository
DSpace https://erepo.uef.fi
Rinnakkaistallenteet Luonnontieteiden ja metsätieteiden tiedekunta
2017
New Algorithm for Simulation of Proteoglycan Loss and Collagen
Degeneration in the Knee Joint: Data from the Osteoarthritis Initiative
Mononen Mika
Wiley-Blackwell
info:eu-repo/semantics/article
info:eu-repo/semantics/acceptedVersion
© Orthopaedic Research Society All rights reserved
http://dx.doi.org/10.1002/jor.23811
https://erepo.uef.fi/handle/123456789/5147
Downloaded from University of Eastern Finland's eRepository
1
New Algorithm for Simulation of Proteoglycan Loss and Collagen Degeneration in the Knee Joint: Data from the Osteoarthritis
Initiative
Mika E Mononen1, Petri Tanska1, Hanna Isaksson2 and Rami K Korhonen1,3
1Department of Applied Physics, University of Eastern Finland, Kuopio, Finland
2Department of Biomedical Engineering, Lund University, Lund, Sweden
3Diagnostic Imaging Centre, Kuopio University Hospital, Kuopio, Finland
Keywords: Knee joint, Articular cartilage, Osteoarthritis, Finite element analysis, Degeneration, Collagen, Proteoglycan
Corresponding author:
Mika E. Mononen, Ph.D.
Department of Applied Physics University of Eastern Finland POB 1627, FI-70211 Kuopio Finland
Tel. +358 40 705 0416 Fax. +358 17 162131
E-mail: mika.mononen@uef.fi
AUTHOR CONTRIBUTIONS
Mika E. Mononen; design of the study, design and coding the degeneration algorithm, conducting simulations and analyses, preparing all figures, main writer (responsible for the integrity of the work as a whole).
Petri Tanska; design of the study, helping in designing the degeneration algorithm, interpretation of the results, critical revision of the article for important intellectual content;
Hanna Isaksson; design of the study, helping in the design of the degeneration algorithm, interpretation of the results, critical revision of the article for important intellectual content;
Rami K. Korhonen; original design of the concept, helping in the design of the degeneration algorithm, interpretation of the results, critical revision of the article for important intellectual content;
2 ABSTRACT
Osteoarthritis is a harmful joint disease but prediction of disease progression is problematic.
Currently, there is only one modeling framework which can be applied to predict the progression of knee osteoarthritis but it only considers degenerative changes in the collagen fibril network. Here, we have developed the framework further by considering all of the major tissue changes (proteoglycan content, fluid flow and collagen fibril network) occurring in osteoarthritis. While excessive levels of tissue stresses controlled degeneration of the collagen fibril network, excessive levels of tissue strains controlled the decrease in proteoglycan content and the increase in permeability. We created four knee joint models with increasing degrees of complexity based on the depth-wise composition. Models were tested for normal and abnormal, physiologically relevant, loading conditions in the knee. Finally, the predicted depth-wise compositional changes from each model were compared against experimentally observed compositional changes in vitro. The model incorporating the typical depth-wise composition of cartilage produced the best match with experimental observations. Consistent with earlier in vitro experiments, this model simulated the greatest proteoglycan depletion in the superficial and middle zones, while the collagen fibril degeneration was located mostly in the superficial zone. The presented algorithm can be used for predicting simultaneous collagen degeneration and proteoglycan loss during the development of knee osteoarthritis.
3 INTRODUCTION
Osteoarthritis (OA) is a degenerative, complex, multi-faceted joint disease, causing eventually substantial pain in the joint.The first signs of cartilage degeneration are proteoglycan (PG) loss and collagen degeneration, leading to higher permeability of the tissue1-6. The mechanisms driving these changes during the progression of OA are still not fully understood, but it has been suggested that excessive stresses, strains and cumulative loadings in cartilage play a major role7-10.
When OA is initiated and progressed due to external loading (e.g. continuous overloading, injurious loading leading to anterior cruciate ligament rupture or some other form of joint injury) or internal factors (e.g. aging, inflammation due to joint injury), homeostasis inside cartilage is disturbed; the cartilage tissue starts to fail to perform its normal functions primarily, due to the superficial collagen fibrillation and depletion of PGs1-3 (Figure 1). General understanding of the initiation and progression of cartilage degeneration is still vague. Previous studies have indicated that cartilage overloading can result in cell death and release of inflammatory cytokines, inhibiting PG production and leading to PG loss11-14. This was suggested to subject the collagen network to damage and degradation by changing the loading of the network12, 13. PG loss to synovial fluid in early OA has also been suggested to occur secondarily due to collagen network damage and altered integrity15, 16, which can be initiated by tissue overloading15, 17.
Computational models of the knee joint are able to simulate knee joint reaction forces, stresses and strains within the knee joint tissues during physiological joint loading18-22. Especially in finite element (FE) based approaches20-22, constitutive material models of soft tissues in the knee have varied from a simple isotropic linear elastic material22 to very complex anisotropic materials, which are able to take into account the effects of composition and structure (proteoglycan, collagen, interstitial fluid) inside the knee joint cartilage20, 21. Nonetheless, although the above models can
4 provide an estimation of stresses and strains inside cartilage, they do not offer any solution to predict changes in tissue composition and structure occurring during the progression of OA.
Currently, there are only a few computational approaches available which are able to simulate cartilage degeneration15, 17, 23. In these approaches, the degeneration is triggered either by excessive levels of strains or stresses or a combination of both. In a recent study17, a cartilage degeneration algorithm was implemented for the first time in the geometry of the knee joint. In this algorithm, cartilage degeneration was controlled by the excessive and cumulatively accumulated maximum principal stresses within the knee joint cartilage during physiologically relevant loading. This parameter is typically considered to reflect tissue failure. Collagenous cartilage tissue is often modeled as a fibril reinforced material24-26, in which local maximum principal stresses are primarily oriented along the tensile direction of the collagen fibrils. Even though the model predictions (degree of collagen fibril network degeneration) were in good agreement with experimentally determined follow-up radiographic grades of OA, the PG depletion and increase in tissue permeability were omitted in this approach6, 27, 28.
It has been shown that cell death and PG depletion of cartilage may occur as a result of excessive levels of compressive and/or shear strains29, 30. On the other hand, deviatoric strain takes into account several strain components and was earlier implemented in a cartilage model to control the nonfibrillar matrix damage in indentation31. It has also been shown that the PG loss is related to increased cartilage permeability4-6.
Cartilage is also a highly inhomogeneous tissue with an arcade-like collagen architecture; the collagen and PG contents increase as a function of tissue depth in human cartilage32. This inhomogeneity is responsible for the depth-wise tissue tension and compression, which should be
5 directly reflected in altered maximum principal stresses and deviatoric strains. Experiments have also shown cartilage failure to be dependent on collagen content33, 34. Some of these aforementioned mechanisms have already been implemented and tested in computational damage models of cartilage15, 23. However, they have not been tested simultaneously in joint level computational models.
The aim of this study was to further develop the previously presented degeneration algorithm17 by applying a previously validated fibril reinforced poroviscoelastic material model of cartilage24-26 in the knee in order to investigate the effect of physiologically relevant loading on the PG depletion and altered fluid flow in conjunction with collagen fibril degeneration during the progression of knee OA. Results from the developed algorithm were compared with the experimentally determined depth-wise histological data (Figure 1). We hypothesized that experimentally observed depth- dependent PG depletion and collagen fibril degeneration2, 3, 32, 35-37
can be simultaneously simulated at the joint level based on excessive tissue deviatoric strains and maximum principal stresses, adopted from earlier experimental in vitro studies10, 29, 30, 34
. We also hypothesized that this result can only be obtained if the material model for cartilage includes depth-dependent properties of the tissue. The presented results may reveal potential mechanisms and interactions behind PG depletion and collagen degeneration caused by knee joint overloading.
METHODS Model geometry
The image data used for generating the knee joint model was obtained from the Osteoarthritis Initiative database (OAI - https://oai.epi-ucsf.org). See more details in the supplementary material.
6 Material models
Similarly as in our previous study17, cartilage tissue was considered as a fibril reinforced poroviscoelastic (FRPVE) material. The total stress tensor of the FRPVE material is:
𝛔t = 𝛔nf+ 𝛔f− 𝑝𝐈 , (1)
where 𝛔t is the total stress tensor, 𝛔nf is the stress tensor of the non-fibrillar matrix, 𝛔f is the stress tensor of the fibril network (sum of primary and secondary fibril stresses), 𝑝 is the fluid pressure and 𝐈 is the unit tensor. The non-fibrillar matrix can be expressed by the Neo-Hookean model as follows:
𝛔nf = 𝐾mln(𝐽)𝐽 I + 𝐺𝐽m(F ∙ 𝐅𝐓− 𝐽23𝐈), (2)
where 𝐽 is the determinant of the deformation gradient tensor F. Bulk modulus (𝐾m) and shear modulus (𝐺m) are characterized with the non-fibrillar matrix Young’s modulus (Em) and Poisson’s ratio (υm). Fluid flow inside the cartilage tissue is controlled by Darcy’s law with permeability (k).
The stress for each individual collagen fibril in the fibril direction is calculated as follows:
𝜎f = −2√(𝜎 η
f−𝐸0𝜀f)𝐸ε𝜎̇f+ 𝐸0𝜀f+ ( 𝜂 + 2√(𝜎𝜂𝐸0
f−𝐸0𝜀f)𝐸ε )𝜀̇f, if 𝜀f> 0, (3)
𝜎f = 0, if 𝜀f≤ 0, (4)
where the fibril network properties are characterized by the initial (E0) and the strain-dependent fibril network modulus (Eε) as well as the damping coefficient (η), 𝜀fis the fibril strain, and 𝜎̇f and
7 𝜀̇f are stress- and strain rates, respectively. Each integration point of the continuum element has 17 separate fibrils, which are divided into 4 primary and 13 secondary fibrils. The secondary fibrils represent collagen cross-links and are organized systematically randomly, while the primary fibrils represent the arcade-like organization. The organization of the primary fibrils changes according to tissue depth similarly as in the previous studies24, 25, 38. Fibril stresses in the individual primary 𝜎f,p and secondary 𝜎f,s fibrils are calculated as follows:
𝜎f,p = 𝜌𝐶𝜎f, (5)
𝜎f,s = 𝜌𝜎f, (6)
where 𝜎f is the fibril stress, ρ is the normalized fibril density and C (= 12.16) is the fraction between the primary and secondary fibrils. The total fibril stress tensor is then the sum of individual fibril stresses
𝛔f= ∑𝑡𝑜𝑡𝑓𝑖=1 𝜎f,i 𝑒̅f,𝑖⨂ 𝑒̅f,𝑖, (7)
where ⨂ denotes the dyadic product, totf is the total number of fibrils and 𝑒̅f,𝑖 is the fibril direction vector (orientation).
Since only short-term loading was applied, and because our main focus was on analyzing cartilage degeneration, only forces through the meniscal tissues (meniscal support forces) during loading were needed (not stresses/strains within the menisci). Therefore, meniscal tissues were considered as a transverse isotropic linear elastic material, similarly as in several previous knee joint modeling studies17, 39-41. In order to mimic the compression-tension behavior in the meniscus, the axial and radial Young’s moduli were designated as 20 MPa while the circumferential Young’s modulus was
8 140 MPa. The in-plane and out-of-plane Poisson’s ratios were 0.3 and 0.2, respectively, and the shear modulus was set to 57.7 MPa17, 39-41.
Degeneration algorithm
Collagen fibril network degeneration Based on previous studies7, 9, 42-44
, the collagen fibril network was assumed to be degenerated due to excessive and cumulatively accumulated tensile tissue stresses (positive maximum principal stresses from the Cauchy stress tensor) within the knee joint cartilage during physiologically relevant loading by local changes in the primary and secondary fibril stresses (see equations 5 and 6)17:
𝜎f,p,𝑒𝑙,𝑖 = 𝜌𝐷𝑒𝑙,𝑖𝐶𝜎f, (8)
𝜎f,s,𝑒𝑙,𝑖 = 𝜌𝐷𝑒𝑙,𝑖𝜎f, (9)
where el is the location (element number) and i is the current iteration number. The current degree of fibril degeneration 𝐷𝑒𝑙,𝑖 (1 = no degeneration; 0 = total degeneration) was written as follows17:
𝐷𝑒𝑙,𝑖 = 1, if 𝑖 = 1, (10)
𝐷𝑒𝑙,𝑖 = 𝐷𝑒𝑙,𝑖−1− (𝐷𝑒𝑙,𝑖−11.5√∑𝑇𝑂𝑇𝑡=1𝐷𝑒𝑙,𝑡× 𝐼𝑁𝐶𝑡) , if 𝑖 > , (11)
where TOT is the total number of required time points (t) during each iterative loading step and INCt is the duration of each time increment. The calculated fibril degeneration factor 𝐷𝑒𝑙,𝑡 was assumed to be dependent on the time increment (phase of the gait load), i.e., 𝐷𝑒𝑙,𝑡 was calculated separately for each element and time point during the stance phase as follows17:
9 𝐷𝑒𝑙,𝑡 = ((𝑆𝑒𝑙,𝑡−𝑇100f)/𝑇f) , if 𝑆𝑒𝑙,𝑡 > 𝑇f, (12)
𝐷𝑒𝑙,𝑡 = 0, if 𝑆𝑒𝑙,𝑡 ≤ 𝑇f, (13)
where 𝑇f is the threshold limit for the collagen fibril degeneration and 𝑆𝑒𝑙,𝑡 is the local tensile stress value (positive maximum principal stress) at each time point.
Proteoglycan depletion
A similar principle as used for collagen degeneration was applied for PG depletion. Instead of maximum principal stress, in the present study, excessive deviatoric strain was assumed to trigger PG loss, similarly as described before in the computational damage model of cartilage31.Since the PG content has been associated with the non-fibrillar matrix modulus (Em)28, the PG depletion was related to a decrease in the non-fibrillar matrix modulusafter each iteration (i) as follows:
𝐸m𝑒𝑙,𝑖 = 𝐸m𝑒𝑙,𝑖−1∙ 𝑃𝐺𝑒𝑙,𝑖, (14)
where 𝑃𝐺𝑒𝑙,𝑖 is the current degree of PG degeneration. Similarly as in our previous study and collagen fibril degeneration (above)17, and consistent with other earlier suggestions7, 9, 43, 44
, both excessive (i.e. when the threshold is exceeded at a certain time point) and cumulatively accumulated mechanical signals were considered to cause PG depletion. Thus, the contribution of excessive and accumulated tissue strains to PG depletion was considered with a similar principle as in the case of collagen degeneration (see equation 11):
𝑃𝐺𝑒𝑙,𝑖 = 𝑃𝐺𝑒𝑙,𝑖−1− (𝑃𝐺𝑒𝑙,𝑖−11.5√∑𝑇𝑂𝑇𝑡=1𝑃𝐺𝑒𝑙,𝑡× 𝐼𝑁𝐶𝑡) , (15)
10 where the calculated PG degeneration factor 𝑃𝐺𝑒𝑙,𝑡 was related to deviatoric strain 𝜀dev𝑒𝑙,𝑡 and assumed to be dependent on the time point (phase of the gait load) and was calculated as follows:
𝑃𝐺𝑒𝑙,𝑡 = 13√𝜀dev𝑒𝑙,𝑡 − 𝑇εdev , if (𝜀dev𝑒𝑙,𝑡 > 𝑇𝜀dev), (16)
𝑃𝐺𝑒𝑙,𝑡 = 0, if (𝜀dev𝑒𝑙,𝑡 ≤ 𝑇𝜀dev), (17)
where 𝑇𝜀dev is the threshold level for the initiation and progression of the PG depletion. The equation for deviatoric strain is defined as follows:
𝜀dev𝑒𝑙,𝑡 = 13 √(𝜀p,1− 𝜀p,2)2+ (𝜀p,1− 𝜀p,3)2+ (𝜀p,2− 𝜀p,3)2 , (18)
where 𝜀p,𝑖 are the principal strains of the “integrated” total strain tensor (output variable E in Abaqus, see more details in supplementary materials) for an element el at time point t (phase of the gait load). Finally, in order to avoid convergence problems with zero stiffness values in elements, we assumed that the non-fibrillar matrix modulus would never reach zero. Thus, the maximum amount of PG depletion was set to 90% (see the initial values for Em from Table 1).
Cartilage permeability
According to previous experimental data4-6, the increase in cartilage permeability was related to the decrease in PG content. Otherwise permeability 𝑘𝑒𝑙,𝑖 was constant during each loading cycle and was implemented as follows:
𝑘𝑒𝑙,𝑖 = 𝑘0(1 + 𝑤(1 − 𝑃𝐺𝑒𝑙,𝑖)), (19)
11 where 𝑘0 is the initial permeability, 𝑤 is the factor defining the increase in cartilage permeability as a result of a decrease in the PG content. Based on selected experimental data of the relationship between the non-fibrillar matrix modulus and permeability due to cartilage degeneration1, 6, the factor was set to 𝑤 = 3. On the other hand, the maximum increase in permeability was set to 3-fold compared to the initial values (see the initial values from Table 1 and more details in the supplementary material).
Cartilage inhomogeneities and degeneration thresholds
To understand the importance of cartilage inhomogeneity and different mechanisms during cartilage degeneration, four different models with increasing complexities were constructed (See material properties and differences in Table 1). In all models, in terms of the material parameters in the FRPVE model, the fraction between individual primary and secondary fibrils (C), fluid fraction (nf), viscoelastic damping coefficient (η), Poisson’s ratio (νm) and initial (E0) and strain-dependent fibril network moduli (Eε) remained unchanged, whereas permeability (k), collagen degeneration level (D) and PG depletion level (PG) changed due to cartilage degeneration caused by excessive levels of maximum principal stresses and deviatoric strains (see equations 11, 15, 19).
Inhomogeneities (see Table 1)
Model 1: This was the reference model with the same level of homogeneity as in the previous study17; arcade-like collagen orientation, depth-dependent fluid fraction, homogeneous distribution of proteoglycans (considered with the distribution of Em28), homogeneous distribution of collagen fibril density (ρ), constant threshold initiating collagen network degeneration and PG loss.
12 Model 2: In order to investigate the importance of the depth-wise PG distribution on cartilage degeneration, a depth-wise distribution of Em was applied to Model 1 according to an earlier experimental study27. Other properties were the same as in Model 1.
Model 3: In order to investigate the importance of depth-wise collagen distribution on cartilage degeneration, depth-dependent distribution of ρ32 was applied to Model 2. Other properties were the same as in Model 2.
Model 4: Since cartilage failure stress in tension has been shown to be dependent on the collagen content, the depth-dependent threshold for initiating collagen fibril network degeneration was added to Model 3 based on the depth-dependent fibril density33, 45. Other properties were the same as in Model 3. Different depth-dependent threshold levels for PG depletion were not tested similarly as with the collagen fibril degeneration, since currently there are limited experimental data available to reveal how different strain levels in different depth-wise zones affect PG damage and loss.
Thresholds for collagen and PG degeneration (see Table 1)
The threshold limit for the initiation and progression of collagen degeneration was based on earlier experiments and our theoretical predictions, and it was assumed to be 7 MPa10, 17 in all models, except in Model 4, in which the depth-wise threshold value was related to the normalized depth- wise collagen fibril volume density 𝜌𝑧𝑜𝑛𝑒 as follows:
𝑇𝑧𝑜𝑛𝑒 = (𝜌𝜌𝑧𝑜𝑛𝑒
f,s ) 𝑇f, (𝑧𝑜𝑛𝑒 = superficial (f, s), middle (f, m) or deep (f, d)) (20)
where 𝑇𝑧𝑜𝑛𝑒is the threshold stress level for the fibril network failure (maximum principal stress) at different depth-wise zones. Linear relationship between the collagen density and limit for the fibril
13 network failure was justified by experimental measurements for bovine and human cartilage33, 45. The normalized fibril volume densities in different zones (𝜌f,s, 𝜌f,m, 𝜌f,d) were based on a previous study32. Now, the zone-dependent fibril degeneration factor 𝐷𝑒𝑙,𝑡 in equations (12) and (13) was calculated as follows:
𝐷𝑒𝑙,𝑡 = ((𝑆𝑒𝑙,𝑡−𝑇𝑧𝑜𝑛𝑒100 )/𝑇𝑧𝑜𝑛𝑒) , if 𝑆𝑒𝑙,𝑡 > 𝑇𝑧𝑜𝑛𝑒, (21)
𝐷𝑒𝑙,𝑡 = 0, if 𝑆𝑒𝑙,𝑡 ≤ 𝑇𝑧𝑜𝑛𝑒. (22)
In the previous theoretical study simulating indentation loading31, the deviatoric strain limit of 30%
was used as a threshold for initiation of PG depletion. This strain limit was based on a previous experimental study46, where strain-rates between 1.7% s-1 and 17% s-1 were applied for PVA hydrogel samples. During physiologically relevant loading of the knee, the peak strain-rate can be as high as ~100% s-1, which may influence the strain threshold for the initiation and progression of PG depletion. However, this kind of data is currently not available. Furthermore, a much lower shear strain (~7-12%) was proposed recently for cell death, and presumably for subsequent PG loss, in another experimental study applied for cartilage30. Therefore, different threshold levels for the initiation and progression of the PG depletion (𝑇𝜀m) were first tested (from 10% to 30%) and the results are shown for 15% and 17% strain thresholds, because they produced experimentally observed results.
Simulations
Similarly as in our previous study17, we simulated 100 consecutive degeneration iterations (arbitrary time) to evaluate collagen degeneration and proteoglycan depletion during the progression of knee OA within the medial compartment during simplified gait (Figure 2). Detailed boundary conditions and their justifications can be found from our previous study17. After the simulations, degenerated
14 volumes were calculated in each model from those elements in which the degrees of PG depletion (𝑃𝐺𝑒𝑙,𝑖, 𝑖 = 100) and fibril degeneration (𝐷𝑒𝑙,𝑖, 𝑖 = 100) were less than one (1 = no degeneration, 0
= total degeneration). Average degenerative changes from those elements (PG depletion, fibril degeneration and increase in permeability) were also calculated. Degenerated volumes and average degenerative changes were calculated in a depth-dependent manner (superficial, middle and deep) by dividing the tissue into three equally thick layers (33% of the total thickness, for more details see in supplementary material). These volumes reflect those elements that experienced degenerative changes in tissue properties, while the true volume of each element remained unchanged (see more details from the supplementary material). Finally, the results were compared against experimental observations (Figure 1, 1, 37).
RESULTS
Under normal joint loading, collagen fibril degeneration was negligible in all models (data not shown). A slight PG depletion occurred at the superficial zone (surface), but it was less than 2% in all models.
During overloading of the knee joint, the volumes of the PG depletion were substantial with both thresholds (i.e., 15% and 17%) and they became altered with increasing complexity of the model (Figure 3, top left and right). The most substantial change occurred in the deep zone, where the volumes of PG depletion were decreased with increasing complexity of the model (Model 1 vs Models 2-3). However, implementation of the depth-dependent threshold level for initiating collagen fibril damage did not alter PG depletion (Model 3 vs. Model 4). In the analysis of the fibril network degeneration (Figure 3, bottom left and right), only slight changes in the volume of the degenerated fibril network were observed due to implementation of the depth-wise PG distribution (Model 1 vs. Model 2). Implementation of the depth-dependent collagen fibril density decreased slightly the volume of the degenerated fibril network in the superficial zone (Model 1 vs Model 3),
15 whereas the implementation of the depth-dependent threshold level for initiating collagen fibril damage decreased the volume of the degenerated fibril network down to zero in the deep zone (Model 3 vs Model 4).
In the analysis of the average PG depletion, increase of permeability and average degeneration in collagen network, the degenerative changes were altered when increasing complexity was introduced into the model (Figure 4). By percent, the simulated average PG depletions and average increases of permeability were more severe than the simulated average collagen degeneration. The threshold of 15% for PG depletion produced severe PG depletion (>40%) and a major increase in permeability (>90%) in each depth-wise zone (Figure 4, left column). The threshold of 17% for PG depletion produced more depth-dependent changes in those parameters, and these were consistent with the experimental observations (Figure 4, middle column vs right column). A more detailed evaluation revealed that inclusion of the depth-dependent PG content changed the maximum PG depletion and maximum increase of permeability from the superficial zone to the middle zone (Figure 4, middle-column, Model 1 vs. Model 2). The additional increase in model complexity (Models 3 and Model 4) had only a minor influence on the depth-dependent degeneration.
The trends between the volume of degenerated collagen fibril network and average degeneration of collagen fibril network were similar in each model (Figure 3 vs. Figure 4). In Models 1-3, the simulated average degeneration in collagen fibril network was ~21-24% in the superficial zone and
~5% in the deep zone, respectively, but in contrast, no fibril degeneration occurred in the middle zone. When the depth-dependent thresholds were included (Model 4), the fibril degeneration was observed only in the superficial zone (Figure 4, bottom-row). Once again, the best match with the experimental observations was found from Model 4.
16 DISCUSSION
In the current study, a novel model for simulating cartilage degeneration was developed and implemented into a knee joint geometry. The amount of degeneration of the nonfibrillar matrix (PGs) was controlled by the deviatoric strain and that of the collagen network by the maximum principal stress. Fluid flow was assumed to be directly related to the PG matrix properties. While the maximum principal stress threshold for collagen damage and degeneration was adopted from earlier studies10, 17, 34, the sensitivity of different strain values for the PG depletion were tested here.
Our results show that the most complex, inhomogeneous model with a 17% strain threshold for initiation and progression of PG depletion captured well the experimentally observed depth-wise PG depletion and collagen degeneration. One particularly interesting finding was that the middle zone displayed most of the PG depletion both in the experiments and model predictions. These results were in accordance with our main hypothesis and suggest that simultaneous PG depletion and collagen fibril degeneration in the knee joint cartilage can be predicted with our novel model.
In the subject during conditions of normal loading, there were no signs of OA during the follow-up period. In an obese subject, the radiographic KL grade increased to 3 during the 4-year follow-up.
Our model also showed that no PG depletion or collagen fibril degeneration occurred during normal loading, whereas a substantial PG depletion and collagen degeneration was simulated due to repetitious overloading (due to overweight). However, it must be noted that we did not have specific information of the depth-wise tissue changes from these subjects. Therefore, here our experimental comparison is conducted against general depth-wise cartilage changes during the progression of OA (Figures 1 and 4). Furthermore, our aim was not to validate the algorithm against clinical data, but to present a new degeneration model, in which parameters can be adjusted when more patient-specific data becomes available.
17 While the threshold for initiation and progression of the collagen fibril network degeneration (tensile stress of 7MPa) was based on previous experimental and computational studies10, 17, different thresholds for initiating the PG depletion (deviatoric strain thresholds from 15 to 30%) were tested based on the literature29-31, 46 and our preliminary tests. The thresholds over 20% did not cause any PG depletion in our model for the chosen physiologically relevant loading. On the other hand, the threshold of 15% did not cause a typical depth-wise PG depletion as reported in many experimental measurements (Figure 1 vs Figure 4, top left). Based on the current knowledge from previous studies2, 3, 36, the PG depletion as well as the collagen fibril network degeneration is initiated typically in the superficial zone prior to reaching down to the middle and deep zones. At later stages of OA, the magnitude of PG loss may be more similar in each zone35, which was actually predicted by our simulation using the threshold of 15% for PG depletion. When the deviatoric strain threshold for initiating PG depletion was increased to 17%, the characteristic depth-wise degeneration was apparent even in the homogeneous model (Model 1); high PG depletion at the superficial zone compared to other zones. After implementing the depth-wise PG and collagen distributions (Models 3 and 4), the characteristic depth-wise PG depletion became clearer. In these models, the well reported progressive PG depletion from the surface till the cartilage-bone interface was observed. This was in good agreement with the hypothesis and the characteristic tissue changes in early or moderate osteoarthritis concluded from experimental studies2, 3, 36.
In contrast to our secondary hypothesis, depth-wise collagen density (in the range that is typical in intact tissue) had only a minor influence on the simulated collagen fibril degeneration; only small changes in the superficial tissue could be observed (Models 2 and 3). However, by adding the depth-dependent threshold for the initiation of the collagen fibril network degeneration (Model 4), as expected, changes in the degeneration of the collagen fibril network were observed.
18 Degeneration of the superficial and middle zone collagen fibril network remained the same, while no collagen fibril network degeneration was observed in the deep zone. Several experimental measurements1-3, 47 show a similar trend (Figures 3 and 4). On the other hand, models without the depth-dependent threshold were able to capture minor collagen network degeneration in the deep zone which is also consistent with the altered collagen orientation in selected experiments (Figure 4). This could point to some role of the subchondral bone in the initiation of OA48. Therefore, based on these results, we cannot say conclusively whether it is important to take into account the depth- dependent collagen network failure limit. In order to resolve this question, well-controlled in vitro studies with specific information of the depth-wise cartilage properties would be needed.
The relationship between the volume of the elements with collagen network degeneration (volume of the elements where 𝐷𝑒𝑙,𝑖 < 1) and degree of fibril degeneration (average decrease in 𝐷𝑒𝑙,𝑖 in the elements where 𝐷𝑒𝑙,𝑖 < 1) (Figure 3 vs. Figure 4) were in accordance with each other and obeyed the trends found from the literature1-3, 47. However, the degenerated PG volume (volume of the elements where 𝑃𝐺𝑒𝑙,𝑖 < 1) and degree of PG loss (average decrease in 𝑃𝐺𝑒𝑙,𝑖 in the elements where 𝑃𝐺𝑒𝑙,𝑖 < 1) (Figure 3 vs. Figure 4) displayed some differences. In the inhomogeneous models 3-4, the volume of the elements experiencing a PG depletion was concentrated within the superficial zone, but the amount of average PG depletion was highest in the middle zone (Figure 3 vs. Figure 4). This latter observation was in good agreement with the previous experimental observation37. These mechanisms can be attributed primarily by the collagen fibril architecture. Obviously, the largest areas of high tensile stresses (maximum principal stresses) are observed in the superficial zone with collagen fibrils oriented in parallel to the cartilage surface. Therefore, both the degree of fibril degeneration and the volume of the elements with fibril degeneration showed the same depth- wise behavior. However, in the middle zone, collagen fibrils do not resist tension or shear forces as much as in the superficial zone (due to more random fibril orientation). This led to higher deviatoric
19 strains of the tissue, causing larger amounts of PG depletion in the middle zone compared to the superficial zone.
Since the change in permeability was linked directly to the amount of local PG content, the increase in permeability was in accordance with the depleted PG amount (Figure 4 and the supplementary material). This assumption was based on earlier experiments1, 6 and the predicted relationship between PG depletion and increased permeability is consistent with other experiments4-6. The algorithm did not consider any increase in the fluid content due to the progression of OA. However, changes in fluid fraction are typically reflected in the permeability and alterations in this parameter were considered in the algorithm.
In the current model for cartilage degeneration, the PG loss was controlled iteratively by the excessive levels of tissue strains (deviatoric strains) under physiologically relevant loading. Due to short-term (dynamic) loading, the collagen fibril network and fluid inside cartilage generate the primary mechanism to resist overloading (strains) within the tissue. For this reason, in the current approach, severe degeneration in the collagen fibril network (softening) will eventually lead to PG loss even if higher thresholds (𝑇𝜀m ≥ 20%) for initiation of the PG depletion would be used. As suggested in earlier studies15, 16, this kind of behavior might also be the result of the PG release into the synovial fluid due to damage in the collagen fibril network.
On the other hand, some researchers have hypothesized that injurious loading or overloading causes initially a PG loss which softens cartilage, leading to increases in collagen fiber strain/stress and later to collagen fibrillation12. Our results were not fully consistent with this hypothesis and showed that a higher amount of PG depletion and softer nonfibrillar matrix did not change the collagen fibril network degeneration (Figure 3). This was mainly due to high fluid pressure during short-term
20 loading resisting tissue deformation, while the nonfibrillar matrix supported only a small fraction from the generated forces1, 26, 49
. This behavior might change by substantially altered threshold for PG depletion and subsequent PG loss, or by different mechanism for collagen degeneration (such as fibril strain50, see the supplementary material). Nonetheless, fibril strains in our model were rather small compared to the possible failure limits proposed in other in vitro studies. On the other hand, other loading conditions, such as creep loading, might be needed (see more details from the supplementary material). In creep loading, fibril degeneration could be rather controlled by fibril strains instead of maximum principal stresses50.
Our model is also purely biomechanical, for instance, it does not take into account potential cytokine release due to injury/overloading and subsequent PG loss12, or possible repair mechanisms of PGs51. These mechanisms would be potential improvements which could be incorporated into the current method if the parameters (biomechanical, biological) were known.
The presented algorithm was able to show some characteristic variations in the depth-wise compositional degeneration of cartilage tissue in the knee joint. The results with the deviatoric strain threshold for PG depletion and maximum principal stress for collagen degeneration produced similar depth-wise degenerative changes in the collagen fibril network and PG content as reported in previous experimental studies2, 3, 32, 35-37
. In the future, the computational simulation presented here could be used for predicting PG loss and collagen degradation during OA and for planning suitable treatment operations to prevent OA (e.g. weight loss). However, in order to validate the model thoroughly, specific data of the changes in PG and collagen distribution during the progression of OA will be needed from diverse subject groups.
21 ACKNOWLEDGEMENTS
European Research Council under the European Union’s Seventh Framework Programme (FP/2007–2013), ERC Grant Agreement no. 281180, the Academy of Finland (projects 286526, 305138), Sigrid Juselius Foundation. Simo Ojanen, M.Sc., is acknowledged for providing experimental digital densitometry data of PG depletion.
CONFLICT OF INTEREST
The authors have no conflicts of interest.
22 References
1. Makela JT, Huttu MR, Korhonen RK. 2012. Structure-function relationships in osteoarthritic human hip joint articular cartilage. Osteoarthritis Cartilage 20: 1268-1277.
2. Turunen SM, Han SK, Herzog W, Korhonen RK. 2013. Cell deformation behavior in mechanically loaded rabbit articular cartilage 4 weeks after anterior cruciate ligament transection. Osteoarthritis Cartilage 21: 505-513.
3. Makela JT, Rezaeian ZS, Mikkonen S, et al. 2014. Site-dependent changes in structure and function of lapine articular cartilage 4 weeks after anterior cruciate ligament transection.
Osteoarthritis Cartilage 22: 869-878.
4. Armstrong CG, Mow VC. 1982. Variations in the intrinsic mechanical properties of human articular cartilage with age, degeneration, and water content. J Bone Joint Surg Am 64: 88-94.
5. Setton LA, Elliott DM, Mow VC. 1999. Altered mechanics of cartilage with osteoarthritis:
Human osteoarthritis and an experimental model of joint degeneration. Osteoarthritis Cartilage 7: 2-14.
6. Makela JT, Han SK, Herzog W, Korhonen RK. 2015. Very early osteoarthritis changes sensitively fluid flow properties of articular cartilage. J Biomech 48: 3369-3376.
7. Sadeghi H, Shepherd DE, Espino DM. 2015. Effect of the variation of loading frequency on surface failure of bovine articular cartilage. Osteoarthritis Cartilage 23: 2252-2258.
8. Patwari P, Cook MN, DiMicco MA, et al. 2003. Proteoglycan degradation after injurious compression of bovine and human articular cartilage in vitro: Interaction with exogenous cytokines. Arthritis Rheum 48: 1292-1301.
9. Seedhom BB. 2006. Conditioning of cartilage during normal activities is an important factor in the development of osteoarthritis. Rheumatology (Oxford) 45: 146-149.
10. Danso EK, Honkanen JT, Saarakkala S, Korhonen RK. 2014. Comparison of nonlinear mechanical properties of bovine articular cartilage and meniscus. J Biomech 47: 200-206.
11. Hauselmann HJ, Flechtenmacher J, Michal L, et al. 1996. The superficial layer of human articular cartilage is more susceptible to interleukin-1-induced damage than the deeper layers.
Arthritis Rheum 39: 478-488.
12. Li Y, Wang Y, Chubinskaya S, et al. 2015. Effects of insulin-like growth factor-1 and dexamethasone on cytokine-challenged cartilage: Relevance to post-traumatic osteoarthritis.
Osteoarthritis Cartilage 23: 266-274.
13. Kar S, Smith DW, Gardiner BS, et al. 2016. Modeling IL-1 induced degradation of articular cartilage. Arch Biochem Biophys 594: 37-53.
23 14. Fick JM, Ronkainen AP, Madden R, et al. 2016. Early in situ changes in chondrocyte biomechanical responses due to a partial meniscectomy in the lateral compartment of the mature rabbit knee joint. J Biomech 49: 4057-4064.
15. Gardiner BS, Woodhouse FG, Besier TF, et al. 2016. Predicting knee osteoarthritis. Ann Biomed Eng 44: 222-233.
16. Rieppo J, Toyras J, Nieminen MT, et al. 2003. Structure-function relationships in enzymatically modified articular cartilage. Cells Tissues Organs 175: 121-132.
17. Mononen ME, Tanska P, Isaksson H, Korhonen RK. 2016. A novel method to simulate the progression of collagen degeneration of cartilage in the knee: Data from the osteoarthritis initiative. Sci Rep 6: 21415.
18. Alexander N, Schwameder H. 2016. Lower limb joint forces during walking on the level and slopes at different inclinations. Gait Posture 45: 137-142.
19. Guess TM, Stylianou AP, Kia M. 2014. Concurrent prediction of muscle and tibiofemoral contact forces during treadmill gait. J Biomech Eng 136: 021032.
20. Adouni M, Shirazi-Adl A, Shirazi R. 2012. Computational biodynamics of human knee joint in gait: From muscle forces to cartilage stresses. J Biomech 45: 2149-2156.
21. Korhonen RK, Tanska P, Kaartinen SM, Fick JM, Mononen ME. 2015. New concept to restore normal cell responses in osteoarthritic knee joint cartilage. Exerc Sport Sci Rev 43: 143- 152.
22. Yang NH, Nayeb-Hashemi H, Canavan PK, Vaziri A. 2010. Effect of frontal plane tibiofemoral angle on the stress and strain at the knee cartilage during the stance phase of gait. J Orthop Res 28: 1539-1547.
23. Hosseini SM, Wilson W, Ito K, van Donkelaar CC. 2014. A numerical model to study mechanically induced initiation and progression of damage in articular cartilage. Osteoarthritis Cartilage 22: 95-103.
24. Wilson W, van Donkelaar CC, van Rietbergen B, Ito K, Huiskes R. 2004. Stresses in the local collagen network of articular cartilage: A poroviscoelastic fibril-reinforced finite element study. J Biomech 37: 357-366.
25. Wilson W, van Donkelaar CC, van Rietbergen B, Ito K, Huiskes R. 2005. Erratum to
“Stresses in the local collagen network of articular cartilage: A poroviscoelastic fibril-reinforced finite element study”[J biomech 37 (2004) 357–366] and “A fibril-reinforced poroviscoelastic swelling model for articular cartilage”[J biomech 38 (2005) 1195–1204]. J Biomech 38: 2138- 2140.
26. Julkunen P, Kiviranta P, Wilson W, Jurvelin JS, Korhonen RK. 2007. Characterization of articular cartilage by combining microscopic analysis with a fibril-reinforced finite-element model. J Biomech 40: 1862-1870.
24 27. Schinagl RM, Gurskis D, Chen AC, Sah RL. 1997. Depth-dependent confined compression modulus of full-thickness bovine articular cartilage. J Orthop Res 15: 499-506.
28. Korhonen RK, Laasanen MS, Toyras J, et al. 2003. Fibril reinforced poroelastic model predicts specifically mechanical behavior of normal, proteoglycan depleted and collagen degraded articular cartilage. J Biomech 36: 1373-1379.
29. D'Lima DD, Hashimoto S, Chen PC, Lotz MK, Colwell CW,Jr. 2001. Cartilage injury induces chondrocyte apoptosis. J Bone Joint Surg Am 83-A Suppl 2: 19-21.
30. Bonnevie ED, Delco ML, Bartell LR, et al. 2017. Chondrocyte death, mitochondrial dysfunction, and apoptosis are mediated by friction and local shear strain. Trans Orth Res Soc 42: 76.
31. Hosseini SM, Wilson W, Ito K, van Donkelaar CC. 2014. A numerical model to study mechanically induced initiation and progression of damage in articular cartilage. Osteoarthritis Cartilage 22: 95-103.
32. Saarakkala S, Julkunen P, Kiviranta P, et al. 2010. Depth-wise progression of osteoarthritis in human articular cartilage: Investigation of composition, structure and biomechanics.
Osteoarthritis Cartilage 18: 73-81.
33. Kempson GE, Muir H, Pollard C, Tuke M. 1973. The tensile properties of the cartilage of human femoral condyles related to the content of collagen and glycosaminoglycans. Biochim Biophys Acta 297: 456-472.
34. Kempson GE. 1982. Relationship between the tensile properties of articular cartilage from the human knee and age. Ann Rheum Dis 41: 508-511.
35. Yin J, Xia Y. 2014. Proteoglycan concentrations in healthy and diseased articular cartilage by fourier transform infrared imaging and principal component regression. Spectrochim Acta A Mol Biomol Spectrosc 133: 825-830.
36. Hayami T, Pickarski M, Zhuo Y, et al. 2006. Characterization of articular cartilage and subchondral bone changes in the rat anterior cruciate ligament transection and meniscectomized models of osteoarthritis. Bone 38: 234-243.
37. Saarakkala S, Julkunen P. 2010. Specificity of fourier transform infrared (FTIR) microspectroscopy to estimate depth-wise proteoglycan content in normal and osteoarthritic human articular cartilage. Cartilage 1: 262-269.
38. Benninghoff A. 1925. Form und bau der gelenkknorpel in ihren beziehungen zur function. II.
der aufbau des gelenkknorpels in seinen beziehungen zur funktion. Zeitschrift fue Zellforschvly 2: 783-862.
39. Vaziri A, Nayeb-Hashemi H, Singh A, Tafti BA. 2008. Influence of meniscectomy and meniscus replacement on the stress distribution in human knee joint. Ann Biomed Eng 36: 1335- 1344.
25 40. Zielinska B, Donahue TL. 2006. 3D finite element model of meniscectomy: Changes in joint contact behavior. J Biomech Eng 128: 115-123.
41. Klodowski A, Mononen ME, Kulmala JP, et al. 2015. Merge of motion analysis, multibody dynamics and finite element method for subject-specific analysis of cartilage loading patterns during gait - differences between rotation and moment driven models of human knee joint.
Multibody Syst Dyn 37: 271-290.
42. Andriacchi TP, Mundermann A, Smith RL, et al. 2004. A framework for the in vivo pathomechanics of osteoarthritis at the knee. Ann Biomed Eng 32: 447-457.
43. Horisberger M, Fortuna R, Valderrabano V, Herzog W. 2013. Long-term repetitive mechanical loading of the knee joint by in vivo muscle stimulation accelerates cartilage degeneration and increases chondrocyte death in a rabbit model. Clin Biomech (Bristol, Avon) 28: 536-543.
44. Miller RH, Edwards WB, Brandon SC, Morton AM, Deluzio KJ. 2014. Why don't most runners get knee osteoarthritis? A case for per-unit-distance loads. Med Sci Sports Exerc 46:
572-579.
45. Williamson AK, Chen AC, Masuda K, Thonar EJ, Sah RL. 2003. Tensile mechanical properties of bovine articular cartilage: Variations with growth and relationships to collagen network components. J Orthop Res 21: 872-880.
46. Stammen JA, Williams S, Ku DN, Guldberg RE. 2001. Mechanical properties of a novel PVA hydrogel in shear and unconfined compression. Biomaterials 22: 799-806.
47. Poole AR, Kobayashi M, Yasuda T, et al. 2002. Type II collagen degradation and its regulation in articular cartilage in osteoarthritis. Ann Rheum Dis 61 Suppl 2: ii78-81.
48. Shirazi R, Shirazi-Adl A. 2008. Deep vertical collagen fibrils play a significant role in mechanics of articular cartilage. J Orthop Res 26: 608-615.
49. Armstrong CG, Lai WM, Mow VC. 1984. An analysis of the unconfined compression of articular cartilage. J Biomech Eng 106: 165-173.
50. Heck TA, Wilson W, Foolen J, et al. 2015. A tissue adaptation model based on strain- dependent collagen degradation and contact-guided cell traction. J Biomech 48: 823-831.
51. van Haaften EE, Ito K, van Donkelaar CC. 2016. The initial repair response of articular cartilage after mechanically induced damage. J Orthop Res 35: 1265-1273.
26 Table 1. Implemented material parameters for femoral and tibial cartilage as the model evolved (Model 1-4). Gray color in the box indicates changes with respect to Model 1. All parameters are based on previous studies17, 27, 32.
FEMORAL CARTILAGE
Material
properties Model 1 Model 2 Model 3 Model 4
Collagen fibril network architecture
Depth-wise arcade-like
Depth-wise arcade-like
Depth-wise arcade-like
Depth-wise arcade-like Em (MPa) 0.215 / 0.215 / 0.215 0.0538 / 0.215 / 0.860 0.0538 / 0.215 / 0.860 0.0538 / 0.215 / 0.860
ρf,s / ρf,m / ρf,d 1 / 1 /1 1 / 1 /1 0.89 / 0.975 /1.16 0.89 / 0.975 /1.16
𝑻𝐟,𝐬/ 𝑻𝐟,𝐦/ 𝑻𝐟,𝐝
(MPa) 7 / 7 / 7 7 / 7 / 7 7 / 7 / 7 7 / 7.67 / 9.12
𝑻𝛆𝐝𝐞𝐯 15% or 17% 15% or 17% 15% or 17% 15% or 17%
𝑬𝟎 (MPa) 0.92 0.92 0.92 0.92
𝑬𝜺 (MPa) 150 150 150 150
𝜼 (MPa s) 1062 1062 1062 1062
C 12.16 12.16 12.16 12.16
νm 0.15 0.15 0.15 0.15
nf 0.8-0.15*hz 0.8-0.15*hz 0.8-0.15*hz 0.8-0.15*hz
k0 (10-15m4/Ns) 6 6 6 6
TIBIAL CARTILAGE
Material
properties Model 1 Model 2 Model 3 Model 4
Collagen fibril Depth-wise arcade-like
Depth-wise arcade-like
Depth-wise arcade-like
Depth-wise arcade-like network architecture
Em (MPa) 0.106 / 0.106 / 0.106 0.0215 / 0.106 / 0.424 0.0215 / 0.106 / 0.424 0.0215 / 0.106 / 0.424
ρf,s / ρf,m / ρf,d 1 / 1 /1 1 / 1 /1 0.89 / 0.975 /1.16 0.89 / 0.975 /1.16
𝑻𝐟,𝐬/ 𝑻𝐟,𝐦/ 𝑻𝐟,𝐝
(MPa) 7 / 7 / 7 7 / 7 / 7 7 / 7 / 7 7 / 7.67 / 9.12
𝑻𝛆𝐝𝐞𝐯 15% or 17% 15% or 17% 15% or 17% 15% or 17%
𝑬𝟎 (MPa) 0.18 0.18 0.18 0.18
𝑬𝜺 (MPa) 23.6 23.6 23.6 23.6
𝜼 (MPa s) 1062 1062 1062 1062
C 12.16 12.16 12.16 12.16
νm 0.15 0.15 0.15 0.15
nf 0.8-0.15*hz 0.8-0.15*hz 0.8-0.15*hz 0.8-0.15*hz
k0 (10-15m4/Ns) 18 18 18 18
hz = normalized cartilage thickness (0 = cartilage surface; 1 = bone-cartilage interface)
27 Figure 1. Experimentally observed PG depletion (left) and collagen fibrillation as estimated by altered fibril orientation (right) in osteoarthritis, as characterized in earlier studies2, 3, 32, 35-37
.
28 Figure 2. Workflow from geometry creation to iterative degeneration algorithm. The top row describes interphases which are considered when calculating tissue stresses and strains, whereas the bottom row presents how predicted stresses and strains (1) affect collagen fibril degeneration (2), proteoglycan depletion (3) and increase of permeability (4) in the algorithm.
29 Figure 3. Depth-dependent PG depletion (top-column) and collagen fibril network degeneration as simulated with our algorithm. The top row figures indicate the volume elements with PG depletion in each cartilage zone (superficial, middle and deep), whereas the bottom row figures indicate the volumes of elements with collagen fibril degeneration. See changes in composition and material properties between different model evolutions (Model 1-4) from Table 1.
30 Figure 4. Average PG depletion (top-row), change in cartilage permeability (middle-row), and collagen fibril network degeneration (bottom row) in each model evolution. Average values were calculated from those elements which experienced PG depletion or collagen fibril degeneration (see volumes from Figure 3). The right column indicates typical experimentally observed changes in the amount of PGs and collagen fibril network orientation in osteoarthritis1, 37. In the experimental observations, different depth-dependent zone thicknesses are considered to be the same as in the models (~33%). See changes in composition and material properties between different model evolutions (Model 1-4) from Table 1.
1 Supplementary material for:
New Algorithm for Simulation of Proteoglycan Loss and Collagen Degeneration in the Knee Joint: Data from the Osteoarthritis Initiative
Mika E Mononen1, Petri Tanska1, Hanna Isaksson2 and Rami K Korhonen1,3
1Department of Applied Physics, University of Eastern Finland, Kuopio, Finland
2Department of Biomedical Engineering, Lund University, Lund, Sweden
3Diagnostic Imaging Centre, Kuopio University Hospital, Kuopio, Finland
METHODS Model geometry
First, the Osteoarthritis Initiative dataset 0.C.2 (OAI - https://oai.epi-ucsf.org) was used to obtain image data from a sagittal dual echo steady-state (SAG 3D DESS) magnetic resonance imaging sequence (TR = 16.32 ms, TE = 4.71 ms, in-plane resolution = 0.36 mm, slice thickness = 0.7 mm) from right knee joints of two test subjects (Figure 2, top-left); normal weight (BMI 24.3) and obese (BMI 35.4), simulating normal loading and overloading. The former subject did not develop OA during a 4-year follow-up period (Kellgren-Lawrence grade 0), while the latter subject had initially healthy cartilage but it became osteoarthritic during the same follow-up period (Kellgren-Lawrence grade from 0 to 3). Then, soft tissues needed for computational simulations (cartilage and meniscus tissues) were manually segmented using Mimics v12.3 (Materialise, Leuven, Belgium) (Figure 2, top-left). The segmented geometries were saved as a surface mesh (.stl format) and converted into solid geometries (.sat format) with Matlab (R2012b, The Mathworks, Inc., Natick, MA, USA).
Finally, the solid geometries were opened in Abaqus finite element package (v6.13-3, Dassault Systèmes, Providence, RI, USA) and cartilages and menisci were meshed identically to our previous study17.
2 DISCUSSION
Limitations
The aim of the current study was to extend a previously developed knee joint model and algorithm based on the current knowledge of the assumed parameters (such as parameters initiating collagen degeneration and PG depletion). Several mechanisms were tested to examine the progression of knee OA in terms of depth-wise collagen fibril network degeneration, PG depletion and increase in permeability. Since direct in vivo data about tissue degeneration was not available (such as could be obtained from MRI), the focus was to compare simulation results to commonly accepted changes in OA1-6. A variety of very specific, highly controlled in vitro experiments could be a better starting point for model validation and to find different threshold values for initiating the progression of OA. Such data might even prove different results. However, our assumptions incorporated into the degeneration algorithm were based on earlier, primarily experimental and computational in vitro studies7-23.
Instead of using zonal thicknesses obtained from experimental studies (e.g., ~15% for superficial,
~35% for middle and ~50% for deep zone11, 24, 25
), the division of cartilage into 3 equally thick layers was mainly based on the model convergence and computational time. The use of thinner superficial zone thickness requires a higher spatial mesh density in order to enable adequate model convergence due to increased strains as degenerative changes progress. As one simulation takes currently about 1 week, adding spatially more elements would increase substantially the computational time. On the other hand, the range of the zone thicknesses is rather large in the literature11, 24-27, and the fundamental behavior of the collagen degeneration and PG loss would not change with any other subdivision of the layers.
