1. Introduction

Articular cartilage covers the bone ends inside a joint to provide lubrication and help withstand loading in locomotion. Its superior mechanical properties are attributed to both its composition and structure. In particular, the collagen network in cartilage forms a unique depth-dependent zonal organization, which can be characterized based on collagen fiber orientation as revealed in polarized light microscopy (PLM) [1–4]. In the superficial zone (SZ), the collagen fibers are aligned in parallel to the surface where the proteoglycan concentration is low. The transitional zone (TZ) contains proteoglycans and thicker collagen fibers with oblique or randomized orientations. The deeper radial zone (RZ) has the highest concentration of proteoglycans, where the collagen fibers are usually perpendicular to the cartilage-bone interface.

Cartilage diseases such as osteoarthritis are often associated with depth progressive structural and functional changes. In early stages of osteoarthritis, fibrillation and disruption occur in the SZ [5] that is a key region to disperse the load applied at the cartilage surface. With the disease progression, the SZ becomes further eroded and fiber structural changes spread into deeper zones, which eventually results in functional failure. Being avascular and aneural, the cartilage cannot heal itself once damaged. The emerging tissue engineering approach is promising in making engineered tissue constructs to repair and regenerate damaged cartilages. To mimic the native cartilage tissues, it is crucial to understand fully the structure-function relationship in cartilage, in particular in the context of depth-dependent zonal structure [6,7].

The depth-dependent mechanical response in cartilage has been investigated extensively [8–13]. Most of these studies used elastography based methods with conventional microscopic imaging [10,11] or magnetic resonance imaging [12,13] to calculate mechanical loading induced strains. These studies have established that the compressive modulus is significantly greater in the deep radial zone than the SZ. The collagen fiber orientation may affect the local cartilage mechanical response [10,14], and may change in response to mechanical loading [11,15]. Fiber reorientation and the associated tension release were proposed as the cause of the “cartilage softening” phenomenon under unconfined (non-laterally bounded) compression [16]. However, most previous elastography studies did not obtain precise zonal structural changes while the sample was under mechanical loading due to the difficulty in quantifying mechanical responses and zonal organization simultaneously.

It has been reported that the zonal collagen organization in cartilage may change in response to external mechanical loading. While reorientation of collagen fibers in the RZ was often observed [8,9,17], the results on superficial cartilage deformation have been inconsistent. For example, compression induced SZ thickening was reported in studies using scanning electron microscopy (SEM) [18], PLM [19], and MRI [9]; whereas another SEM study suggested a constant thickness of the “upper” zone (or SZ + TZ) during the compression [20]. On the contrary, a small angle X-ray scattering study [17] and a MRI study [21] reported a decrease in the upper zone thickness under compression. The conflicting results can be attributed to the lack of a practical technology that can visualize precisely the detailed zonal cartilage structure during the mechanical tests.

Due to the highly organized collagen matrix, cartilage samples show significant optical birefringence that can be quantified using optical polarization based methods such as PLM or polarization-sensitive optical coherence tomography (PSOCT) [22–24]. Recently, a PSOCT based tractography imaging method has shown potential for high-resolution visualization of detailed fiber organization in tissue [25–28]. This method, referred to as optical polarization tractography (OPT), is similar to quantitative PLM in its capability of measuring the optic axis of the cartilage sample to determine its zonal organization [29,30]. However, OPT works in the reflection mode, does not require thin pre-processed samples, and provides volumetric image data. In addition, OPT can also obtain intensity OCT images that can be used to construct mechanical loading induced local strain map as in optical coherence elastography [31,32]. Therefore, OPT provides an opportunity to fully reveal the precise zone-dependent structural and mechanical responses simultaneously, which has not been obtained previously.

In this study, an OPT based method was developed to study the complex relation of zonal architecture and mechanical properties in articular cartilage under unconfined compressive loading. By using simultaneous tractography and elastography imaging, we were able to image detailed zonal organization and their mechanical properties simultaneously in cartilage samples. Our results revealed clear zone-dependent structural and mechanical property changes in response to compressive loading. The data revealed that the conventional RZ contained two sub zones with distinct mechanical properties. The “upper” part of the RZ played a critical role in “absorbing” the mechanical compression. The study also provided strong evidence that zonal fiber organization can greatly affect the cartilage structural and mechanical responses.

2. Materials and methods

2.1 Cartilage sample

We tested 12 cartilage samples excised from six bovine femoral bones of different animals that were 24- to 30-month old. The samples were obtained from a local grocery store and the meat science laboratory at University of Missouri-Columbia. No significant group difference was observed in the results obtained using samples from the two different sources. From each femoral bone, two articular cartilage samples were excised from the central part of the lateral trochlear ridge (LTR). A cartilage section with ~7-10 mm thick remaining bone tissue was first removed from the LTR using a thin diamond saw. Then a small cubic sample (~5 × 5 × 5 mm³) was carefully excised manually by using a precise surgical blade. The 12 trimmed cubic LTR samples were fixed in 10% formalin solution at room temperature (23 C°) for at least 12 hours before the experiment.

2.2 OPT imaging

The fiber structure of the samples was imaged using a 0.85 μm wavelength frequency-domain PSOCT system that has been described in detail previously [27,33]. The image volume covered a 6 × 6 × 1 mm³ sample volume (in B × C × A scans) with an isotropic 6 µm voxel size. Three types of images were obtained: OCT intensity, local birefringence, and local optic axis. The OCT intensity image was formed using backscattered light intensity that revealed the inner structure of the sample. By computing the speckle changes in OCT intensity images before and after compression, the images of displacement and strain can be constructed as described in Section 2.3 below. The birefringence imaged in OPT represents the “apparent” optical birefringence within a measurement plane perpendicular to the incident light [4,34]. Such “apparent” birefringence is the true sample birefringence when the collagen fibers lie within this measurement plane. In general, a higher birefringence value may suggest a higher collagen fiber concentration or a more organized collagen fiber structure.

The image of the optic axis represented the collagen fiber orientation in the sample, which can be rendered in forms of tractography images as illustrated previously [35,36]. The birefringence and fiber orientation images obtained in OPT were similar to those obtained in quantitative PLM [37]. However, OPT works in the reflection mode and can provide 3D depth-resolved subsurface images. The organization of the fibrous organization can be quantified by calculating the fiber alignment index using circular statistics of the orientation data as described previously [33]. A 7 × 7 pixel window size was used in calculating the fiber alignment index which had values within [0, 1] with 1 indicating perfectly aligned fiber structure and 0 indicating randomly organized fibers.

2.3 Compressive mechanical loading

Figure 1(a) illustrates the experimental setup. The sample was located in a chamber filled with 10% formalin solution to avoid dehydration during the test. The compressive force was applied slowly using a translational stage to the cartilage surface along the x-axis to form a uniaxial loading configuration. The applied load was recorded using a force sensor “S” (DFS-BTA, Vernier Software & Technology, Beaverton, OR, USA) attached to the translational stage.

 

Fig. 1 The experimental setup, force-strain data, and cartilage thickness of all 12 samples. (a) A schematic drawing of the sample holder and the imaging geometry. (b) An example of a cartilage sample with labeled coordinates. (c) The tractography of the fiber orientation in the superficial zone as imaged from the cartilage synovial surface. The x, y, z coordinates were aligned with the B, C, and A scanning directions, respectively. (d) The force-strain data from 12 samples. (e) The measured cartilage thickness at different amount of compression.

Download Full Size | PDF

Figure 1(b) shows an example LTR cartilage carefully excised from a bovine femur bone. The x, y, z coordinate directions were aligned with the B, C, and A scanning directions, respectively. It is important to section the image plane (x-y) in parallel to the collagen fibers in the SZ so that the collagen fibers were located inside the x-y plane. Therefore, all samples were first imaged from the synovial surface (Fig. 1(b)) to determine the fiber orientation in the SZ (Fig. 1(c)), which was then used to define the y-axis. After mounting the sample to the setup, the x-y side of the cube (Fig. 1(b)) was faced toward the imaging system during the mechanical test. The “axial” compressive loading applied along the x-axis was perpendicular to the synovial surface or the fibers in the superficial zone. Therefore, the strains induced were primarily along the axial direction (x-axis); whereas the potential strains along the two lateral directions (y and z) were small as shown previously [38]. Between the two lateral directions, imaging perpendicularly to the superficial zone fibers is preferred and adopted in this study because imaging in parallel to the fiber may induce unnecessary artefacts [4].

A small preloading force (0.4 ± 0.2 N) was applied to ensure a good contact between the force stage and the cartilage surface. After the preloading was fully equilibrated, the sample was imaged at 0% strain first, and the compression strain was gradually increased at a nominal step of ~4% until reaching ~20% global compression (Fig. 1(d)). The changes in the cartilage thickness are shown in Fig. 1(e). The actual bulk compression in cartilage was calculated from the intensity images as 4.1 ± 0.54%, 8.0 ± 0.58%, 11.8 ± 0.77%, 16.2 ± 0.69%, and 19.9 ± 0.63% with the corresponding bulk force measured as 1.7 ± 0.5 N, 3.1 ± 0.9 N, 4.8 ± 1.2 N, 6.5 ± 1.8 N, and 8.2 ± 2.6 N, respectively. Two samples were only measured up to 12% and 16% bulk compression due to the sample mount malfunction. In all samples, a relaxation period of 10–40 minutes (increasing with strain) was allowed after each increment so that the force reached a steady state before the sample was imaged.

2.4 Images of local strain field

The non-rigid Demon registration algorithm [39] was utilized to calculate the displacement field induced by the compression. Using the two OCT intensity images I_i and I_{i + 1} obtained at two consecutive compression levels i and i + 1, respectively, the Matlab function imregdemons was applied to obtain the pixel-wise horizontal (u_i) and vertical (v_i) displacement between the i^th and (i + 1)^th compression. The displacement results (u_i,v_i) and (u_{i + 1}, v_{i + 1}) were superimposed to obtain the displacement between the i^th and (i + 2)^th compression.

To calculate the strain map from the displacement, a 1st order linear transformation was utilized to model the displacement within a 20 × 20 pixel (120 × 120 μm) evaluation window [40]:

2.5 Statistics

To accommodate the missing data points in two of the 12 samples, a linear mixed model was used with maximum likelihood method to determine the effects of zone and bulk compression on the structural (thickness, birefringence, and fiber alignment) and mechanical property (Young’s modulus) of the cartilage samples. The effect of the compression as represented in nominal bulk compression was treated as repeated measure and used as a covariate; whereas the zone was treated as a categorical factor. The fixed effects of zone, compression, and the zone × compression interaction were tested in the mixed model. In addition, the effect of compression was tested in each zone separately by using the Friedman test with the “Exact” method. All statistical analyses were conducted with the use of SPSS 25 (SPSS Inc., USA). Alpha was set at 0.05 for all tests.

3. Results

Figure 2(a) shows the example OPT images obtained in a cartilage sample during a sequence of compressions with global axial strains from 0% to 20%. These en face images (in the x-y plane) were extracted at 360 µm beneath the sample side surface. Similar results were obtained at other depths. The boundary between the cartilage and bone was readily distinguishable in the intensity images where the bone appeared significantly darker. Multiple structural features (circled) could be seen in the intensity images. Tracing the locations of these markers during the compression sequence revealed some interesting patterns of displacement over the cartilage depth. The markers close to the cartilage surface had larger displacement; whereas the features close to the calcified zone nearly maintained their positions throughout the compression.

 

Fig. 2 (a) An example of the OPT intensity, tractography, and birefringence images acquired under different amount of bulk compression (0% to 20%). The size bar in the 0% intensity image indicates 0.2 mm. The dashed boxes mark the original cartilage region before compression. The dashed lines illustrate the boundaries between superficial zone (SZ), transitional zones (TZ), and radial zone (RZ). Also shown are the 1D depth profiles of (b) image intensity, (c) birefringence (∆n), (d) fiber orientation, and (e) fiber alignment under different bulk compression of the sample. These 1D results were calculated by averaging over the lateral direction for each image.

Download Full Size | PDF

The tractography images revealed the classic zonal organization in the cartilage. The fibers were predominantly parallel to the surface at −8.3° ± 20.2° in the SZ, whereas the fibers in the RZ were highly organized and aligned at 89.6° ± 9.7° in relation to the cartilage surface. A very thin TZ layer was beneath the SZ and had a broad range of orientations 55.7° ± 41.6°. During the compression, the SZ thickness slightly decreased in some regions at 4% while it appeared thicker at 20% strain. The fiber orientation varied significantly in the TZ with interwoven parallel and perpendicular fibers.

The OPT birefringence images showed a depth-dependent pattern that was typically in PLM results [37]. A thin SZ close to the cartilage surface had high birefringence values of 9.7 × 10⁻⁴ ± 3.1 × 10⁻⁴. The TZ appeared as a band with small birefringence of 7.4 × 10⁻⁴ ± 2.2 × 10⁻⁴. The birefringence started to increase toward the RZ and reached 9.3 × 10⁻⁴ ± 2.5 × 10⁻⁴. The birefringence significantly increased to 1.8 × 10⁻³ ± 5.6 × 10⁻⁴ when approaching the calcified zone. At a greater level of compression, the subsurface region of low birefringence in the upper part of the RZ moved deeper into the lower part of the RZ.

Figures 2(b)-(e) show the quantified 1D depth profiles of intensity, orientation, birefringence (∆n), and fiber alignment for the sample in Fig. 2(a) under different amount of compression. These depth profiles were calculated by averaging over a 1-mm long section along the lateral direction (y-axis) for each image. Taking into consideration of the sample surface movement due to compression, all depth profiles appeared to be quite consistent under different compression. The 1D intensity profiles (Fig. 2(b)) showed minimal changes within the cartilage region and then became lower toward the bone areas, which was consistent with the appearance in the intensity images shown in Fig. 2(a).

The birefringence profiles started with a high value (~1.1 × 10⁻³) in the SZ, decreased to ~0.6 × 10⁻³ in the TZ, and then gradually increased to ~0.9 × 10⁻³ in the RZ. The birefringence increased (~1.7 × 10⁻³) into the calcified tissue. To avoid angle wrapping, the orientation values were converted to the range of [-180°, 180°] in Fig. 2(d). The orientations (Fig. 2(d)) showed a typical “step” pattern from ~8° in the SZ to~90° in the deeper RZ as observed in typical PLM results [15]. The alignment index profile gave additional information about fiber arrangement in the cartilage. An alignment index of ~0.93 indicated highly organized fibers in the SZ. The index decreased drastically in the TZ and then went back to ~0.98 and remained almost constant in the deeper RZ. The alignment index decreased again toward the calcified zone.

The boundary between cartilage and bone was automatically identified by fitting decreasing portion of the 1D intensity profile (Fig. 2(b)) using a hyperbolic tangent function y = p₁ × tanh((x-p₂)/p₃) + p₀. The bone boundary was selected at the location of the peak of the first derivative of the fitted curve. The boundaries of the SZ, TZ, and RZ were determined based on the 1D fiber orientation profile (Fig. 2(d)). The decreasing portion within the first 500 μm of the orientation profile was fitted using also a hyperbolic tangent function. The fitted function was then normalized to values within [0, 1]. The boundaries between superficial-transitional zones and transitional-radial zones were determined as the depths where the normalized hyperbolic tangent function had values of 0.9 and 0.1, respectively. A similar fitting approach has been previously applied in a PLM study to determine zonal boundaries [41].

Figure 3 shows the images of displacement (u and v) and strain (ε_y and ε_x) calculated from the intensity images shown in Fig. 2. Using the aforementioned method, the cartilage area and the boundaries between SZ, TZ, and RZ at the end of each bulk compression were marked using boxes and dashed lines. As a result of the fixed bone position, the axial displacement (v) was the highest at the sample surface and then gradually decreased with depth. However, Fig. 3(a) clearly illustrated that the amount of axial strain reached maximum within the TZ at all compression levels. The axial strain started to decrease in the RZ and approached zero toward the bone. In comparison with the axial displacement (v) and strain (ε_x), the lateral displacement field (u) was negligible for all bulk compression levels and the lateral strain (ε_y) was also very small with no clear pattern observed (Fig. 3(a)).

 

Fig. 3 Example results of displacement and strain measurements calculated from the intensity images shown in Fig. 2. (a) The images of displacement and strain along the depth (x) and lateral (y) directions. The dashed boxes mark the cartilage region at the end of each compression. The dashed lines illustrate the boundaries between superficial zone (SZ), transitional zones (TZ), and radial zone (RZ), as well as between the upper radial zone (UR) and the lower radial zone (LR). The size bar indicates 0.2 mm. Also shown are the corresponding 1D depth profiles of (b) strain and (c) Young’s modulus (E) under different bulk compression for the same sample shown in Fig. 2. The zonal boundaries labeled in (b) and (c) were obtained for the sample before compression.

Download Full Size | PDF

Figure 3(b) shows the 1D depth profiles of (ε_x) calculated by averaging all images over the lateral direction. These 1D curves confirmed that the axial strain was small in the most superficial zone, coincident with the region of a higher birefringence as shown in Fig. 2. It started to increase and reached the maximum in the TZ; and decreased thereafter into the RZ. The corresponding Young’s modulus (E) profiles (Fig. 3(c)) showed that the most superficial zone (within 0.1 mm) had slightly higher stiffness. The modulus decreased quickly thereafter, remained low in the TZ, and then gradually increased until 0.8~0.9 mm in the RZ.

A careful analysis of the 2D axial strain images in Fig. 3(a) suggested that the RZ may be divided into two different regions. The upper RZ (UR) showed a considerably large strain (appeared in light-blue color), whereas the lower region of the RZ (LR) maintained a very small strain (appeared in dark-blue color) throughout the compression. Using the 1D depth profiles of the axial strain, we estimated the boundary between UR and LR at the depth where the accumulated change of strain was less than 20% of the maximal change observed in TZ. The resulting boundary was labeled in Fig. 3.

The axial strain generally increased with compression as expected. It is interesting to note that the Young’s modulus increased with compression in the SZ and lower RZ; but decreased in the upper RZ. This implied that the SZ and lower RZ became stiffer, while the upper RZ became softer during the compression in this sample. The trend change inside the RZ occurred at about the boundary between the UR and LR, which supported the existence of two sub RZ regions with different mechanical properties. Such a “strain-softening” effect has been reported in a few previous studies, although the specific zone linked to this phenomenon has not been investigated [16,42]. This interesting behavior has been explained as the result of overpowering the osmotic pressure induced tension by the compressive loading.

Figure 4 shows the zonal structural parameters (thickness, birefringence, orientation, and fiber alignment) measured at each compression level. The birefringence, orientation, and fiber alignment were averaged values obtained from the central 70% portion of each zone at each compression level. Figures 4(a1)-(a4) show the absolute thickness measured for each zone. The zonal thickness showed different changes in each zone in response to the mechanical compression. To compensate the cartilage thickness variations in different samples, the zone thickness was represented as the percentage of the total cartilage thickness in each sample and at each compression as shown in Figs. 4(b1)-(b4). The SZ contributed between 7~20% of the total cartilage thickness before compression, which increased to 10~30% at 20% compression (Fig. 4(b1)). The TZ made up 2~15% of the total cartilage thickness, and showed little change after compression (Fig. 4(b2)). The upper RZ contributed 20~60% of the total cartilage thickness before compression, which decreased to 6~40% at 20% compression (Fig. 4(b3)). The lower RZ thickness showed a consistent trend of increase with compression: from 20~50% of the total cartilage at pre-compression to 30~64% at 20% compression (Fig. 4(b4)).

 

Fig. 4 The effects of compression on the following zone-specific structural properties: (a1-a4) absolute thickness, (b1-b4) percentage thickness, (c1-c4) optical birefringence, (d1-d4) fiber orientation, and (e1-e4) fiber alignment. The box boundaries indicate the 25^th/75^th percentiles; and the error bars indicate the 10^th/90^th percentiles. The circles indicated outlying data points.

Download Full Size | PDF

The linear mixed model indicated that the zone thickness (%) was significantly different in different zones (F = 53.16, p = 0.000). Although the fixed effect of compression was not significant (p = 1.000), the interaction of zone × compression was highly significant (F = 107.13, p = 0.000), which suggested different compression effects in different zones. The SZ changed significantly due to compressions (χ²(5) = 43.31, p = 0.000, Friedman test). The median SZ thickness increased from 13.4% at zero compression to 21.8% at 20% compression. No thickness change was observed in the TZ (χ²(5) = 1.77, p = 0.893). The Friedman test indicated that the thickness was significantly changed in both the upper RZ (χ²(5) = 48.97, p = 0.000) and the lower RZ (χ²(5) = 50.00, p = 0.000). The median upper RZ thickness decreased from 46.7% at zero compression to 27.0% at 20% compression, whereas the median lower RZ thickness increased from 34.1% at zero compression to 43.2% at 20% compression. However, a careful examination revealed that the actual physical thickness of the lower RZ had little change during the compression (Fig. 4(a4)). The apparent increase in the percentage thickness of the lower radial (Fig. 4(b4)) was the effect of the reduction in the total cartilage thickness.

As shown in Fig. 4(c1)-(c4), the optical birefringence was generally higher in the RZ than in the SZ and TZ. The TZ had the smallest birefringence, while the lower RZ had the highest birefringence. The linear mixed model using birefringence as the dependent variable confirmed a significant fixed effect of the zone (F = 11.17, p = 0.000). It also revealed a significant compression effect (F = 6.75, p = 0.010), but the zone × compression interaction was not significant (F = 1.54, p = 0.204). The Friedman test confirmed a significant effect of compression in the SZ (χ²(5) = 12.97, p = 0.018) and lower RZ ((χ²(5) = 24.46, p = 0.000), but not in the TZ and upper RZ. As a group, the birefringence had a small but consistent decrease with compression in the lower RZ.

As expected, the fibers in the SZ were largely orientated in parallel to the synovial surface, whereas they became perpendicular in the RZ (Fig. 4(d1)-(d4)). The fiber orientation in the TZ varied greatly. The linear mixed model confirmed a significant effect of zone (F = 6.75, p = 0.000); but neither the compression nor the zone × compression interaction was significant. The Friedman test confirmed that the compression effect was not significant in any of the four zones. The outlying points shown as circles in Fig. 4(d3) and (d4) were from the same sample. The fiber orientation in the lower RZ of this particular sample was ~140° that was significantly deviated from the vertical orientation in other samples. In the upper RZ of the same sample, the fiber orientation appeared to increase with compression and approached the fiber orientation in the lower RZ.

As shown in Fig. 4(e2), the fiber alignment was noticeably lower in the TZ, whereas the other three zones had close to 1.0 fiber alignment. The linear mixed model showed that fiber alignment was significantly affected by the zone (F = 83.59, p = 0.000), but not by the compression (F = 0.95, p = 0.332) nor the zone × compression interaction (F = 1.24, p = 0.298). However, the Friedman test revealed a significant effect of compression on fiber alignment (χ²(5) = 28.80, p = 0.000) in the upper RZ where the fiber structure became less organized with compression. A close examination showed that such a decreasing trend in the upper RZ was more dramatic in some samples. For example, in the outlying sample shown as circles in Fig. 4(e3), the fiber alignment decreased from 0.89 at zero compression to 0.69 at 12% compression.

Figure 5 shows the mechanical responses and Young’s modulus measured in different zones at different compression levels. Similar to the measurements in Fig. 4, those values were averaged within the central 70% portion of each zone. The negative sign in the strain values indicated compression. In general, the amount of strain increased with the level of the bulk compression in all zones. As expected, the strain was significantly smaller in the lower RZ, while remarkably increased in the upper RZ. The mixed model showed a significant effect from compression (F = 968.69, p = 0.000) and zone × compression (F = 114.56, p = 0.000).

 

Fig. 5 The average (a1-a4) zonal strain and (b1-b4) Young’s modulus measured at each compression in each zone. The box boundaries indicate the 25^th/75^th percentiles; and the error bars indicate the 10^th/90^th percentiles. The circles indicated outlying data points.

Download Full Size | PDF

Figure 5(b) shows the mean Young’s modulus in different zones at different compression. The Young’s modulus was significantly higher in the lower RZ than the other zones, which was consistent with the observation of small strain in this zone. Young’s modulus appeared to increase with compression in all zones except for the upper RZ. The linear mixed model confirmed that the Young’s modulus was significantly different in different zones (F = 6.07, p = 0.001) and at different level of compression (F = 62.20, p = 0.000). The interaction of zone × compression was also highly significant (F = 56.00, p = 0.000). In the SZ, the Friedman test confirmed that the Young’s modulus changed significantly with compression (χ²(4) = 21.36, p = 0.000). Similar results were observed in the TZ (χ²(4) = 29.76, p = 0.000) and lower RZ (χ²(4) = 32.56, p = 0.000). However, as a group, the change of Young’s modulus was not significant in the upper RZ (χ²(4) = 2.96, p = 0.589), suggesting an overall linear elastic behavior during the compression in the upper RZ.

4. Discussion

This study utilized OPT to image simultaneously the detailed fiber structure changes and mechanical responses in articular cartilage. OPT can image detailed fiber orientation in tissue during compression, which allows a precise determination of the cartilage zonal structural changes due to compression. This is a major improvement over previous mechanical studies where an arbitrary division of articular cartilage thickness into three [43] or even two [38] subdivisions made it hard to measure precisely each zone’s actual contribution.

Although our results indicated that the SZ could become slightly thinner at small compression (~4%), the overall results in Fig. 4 clearly revealed a trend of SZ thickening (in both absolute thickness and percentage thickness), especially at > 8% bulk compression. This result disagree with previous reports of unchanged [20] or decreased [17,21] upper zone (SZ plus TZ) thickness during modest compression; but agreed with other studies using µMRI [44], SEM [18], and PLM [19]. This phenomenon is likely due to the reorganization of collagen fibers when portion of the UR was converted into SZ and TZ during the compression. Such reorganization has been demonstrated in a previous SEM study [18]. The TZ contributed to less than 15% of the total cartilage thickness and no clear change in its thickness was observed.

The majority of the cartilage (70% ~90%) belonged to the RZ where the fibers were perpendicular to the cartilage-bone interface. The Young’s modulus in the RZ increased significantly with depth from the transitional-radial boundary (upper part) and peaked in a region (lower part) toward the bone/calcified region. The modulus in the lower RZ was ~4-12 times higher than that in the upper RZ. In the 12 samples tested in this study, the stiffer lower RZ occupies 27% ~76% of the total RZ and had very small (<6.5%) strain at the 20% compression. The existence of two different regions in the RZ was previously described by Kaab et al. in an SEM study [18] where different changes in structure and fiber orientation were observed under mechanical loading. The greater stiffness in the lower RZ may be due to the increase of calcified cartilage toward the underlying bones [45].

The SZ and TZ were more compliant during the compression and had larger strain than the RZ. However, their absolute contribution to accommodate the external compression was small due to their small thickness. Overall, the upper RZ played the most important role in accommodating the external compression due to its large thickness and compliance. Figure 6(a) shows the portion of the total bulk compression that was attributed to the compression in each zone by calculating the ${\epsilon }_i{t}_i/{{\displaystyle \sum }}^{\text{}}{\epsilon }_i{t}_i$ where ε_i and t_i are the averaged strain and thickness of zone i, respectively. It is clear that the lower RZ contributed the least to the total cartilage compression, whereas the upper RZ contributed the most of the compression. The contribution from the upper RZ increased from 42.5% (median value) at 4% compression to 55.3% (median value) at 20% compression.

 

Fig. 6 The contribution of each zone to the overall bulk compression (a1-a4) and the thickness gained in each zone during the compression (b1-b4). The box boundaries indicate the 25^th/75^th percentiles; and the error bars indicate the 10^th/90^th percentiles. The circles indicated outlying data points.

Download Full Size | PDF

While it was expected to see that the superficial, transitional and lower RZs all became stiffer with compression (Fig. 5) [46], it is intriguing to observe that the SZ became thicker while being compressed. Figure 6(b) further quantified this phenomenon by calculating the thickness gained in each zone using t_ic-(1 + ε_i) t_i0, where ε_i is the averaged strain in zone i; t_i0 and t_ic are the zone thickness before and after compression. It is clear that the SZ gradually gained thickness in response to compression, despite of the fact that compression reduced its original thickness. The lower RZ also gained slightly, whereas no consistent trend was observed in the TZ due to its small size. In sharp contrast, the upper RZ consistently lost thickness in all samples that was in addition to compression-induced reduction. This analysis showed unequivocally that part of the upper RZ gradually reorganized into new TZ and SZ due to compression.

The aforementioned zonal reorganization during compression have been previously explained as the result of “fiber bending” in the RZ [18,47,48], where the fiber orientation change increased strain in the lateral direction. However, our results did not support this explanation because as a group of 12 samples, the fiber orientation was not significantly altered in any of the zones and the lateral strain ε_y remained minimum (Figs. 2 & 4). One sample (the “circle” outlier in Fig. 4) did show large orientation change in the upper RZ by 26.4° at 16% compression in comparison with <5° in the other 11 samples. From the results shown in Fig. 7, the most obvious difference in this particular sample was the different fiber orientation in the upper and lower RZs at 0% compression. The fibers in the upper RZ were vertical (~90°); while the fibers in the lower RZ were ~50° deviated from 90°. In comparison, the angular deviation in the lower RZ was between 0° and 18° in the other 11 samples (Fig. 2). It is plausible that the fibers in the upper RZ bent to align with the fibers in the lower RZ [49]. However, such “bending” behavior would have induced significant deformation in the lateral direction, which was not observed in the strain results. Therefore, we speculate that the realignment of the upper radial fibers observed in this sample may be attributed to non-classical mechanical mechanisms for example, reorganization of collagen crosslinking [50]. More studies are necessary to clarify this observation.

 

Fig. 7 A sample with larger fiber orientation changes in the upper radial zone during compression. (a) intensity, (b) birefringence, (c) tractography, and (d) strain maps.

Download Full Size | PDF

Coincidently, the sample shown in Fig. 7 also had the smallest birefringence and Young’s modulus in the lower RZ among all samples. It had a birefringence of 5.1 × 10⁻⁴ and Young’s modulus of 3.6 MPa in the lower RZ at 16% compression. As a comparison, the sample in Fig. 2 had a birefringence of 8.7 × 10⁻⁴ and a Young’s modulus of 8.1 MPa at the same compression. There appeared to be a correlation between the optical birefringence and Young’s modulus, which was also implied in their corresponding 1D depth profiles. As shown in Fig. 2(c) and Fig. 3(c), both birefringence and Young’s modulus quickly descended within the SZ and then gradually increased in the RZ. Further quantitative analysis suggested a linear correlation existed between the averaged optical birefringence and the logarithm of Young’s modulus. Figure 8 illustrated such correlation calculated at 8% and 16% compression (r² = 0.45 and 0.43, respectively). Similar results were observed at other bulk compression levels. The association between bulk Young’s modulus and collagen content in cartilage has been previously reported using Fourier transform infrared imaging [2]. A higher optical birefringence may indicate a higher collagen fiber density, which in turn leads to a higher tissue stiffness. The results shown in Fig. 8 revealed that such linear relationship also existed in individual zones in particular the SZ, TZ, and upper RZ.

 

Fig. 8 Correlation between Young’s modulus (E) and optical birefringence (∆n) at (a) 8% and (b) 16% compression. The dashed lines are linear regressions: log₁₀E = 5.27 + 1.46 × 10³∆n with r² = 0.45 in (a) and log₁₀E = 5.37 + 1.51 × 10³∆n with r ² = 0.43 in (b).

Download Full Size | PDF

In this study, the Young’s modulus was calculated under the assumption of a constant sample cross-sectional area. In fact, the measured lateral strain was negligible in comparison to the axial strain (Fig. 3). This observation implied a small Poisson’s ratio in cartilage sample, which is in good agreement with previous studies [38,51]. Under such conditions, 2D elastography is sufficient in measuring sample mechanical responses and has been widely used in cartilage studies [42,52,53]. However, 3D elastography is necessary to precisely characterize any subtle effects of anisotropic deformation. It is also important to note that tissue fixation was applied in this study to avoid potential sample structural changes during the lengthy sample preparation and testing process. Under equilibrium conditions, the cartilage mechanical properties are mainly determined by its collagen network that is well preserved using formalin fixation [54]. We compared the sample structural properties before and after compressive loading in a few samples and did not find any obvious alternations. Nevertheless, testing in fresh samples is preferred if a feasible protocol is developed.

5. Conclusion

By using simultaneous tractography and elastography imaging, we were able to delineate the zone-specific structural and mechanical responses of articular cartilage under compression. Experimental data in bovine cartilage samples confirmed the traditional three-zone architecture in cartilage based on the fiber orientation obtained in tractography. In addition, the conventional RZ can be divided into upper and lower RZ due to their distinctly different mechanical properties. The SZ thickening was observed at larger mechanical compression. This phenomenon was attributed to the conversion of the part of the upper dial zone into new SZ and TZ. Moreover, the upper RZ absorbed the most compression especially at larger compression levels. In samples where the lower radial fiber orientation was perpendicular to the cartilage surface, no reorganization in fiber orientation was observed. However, in samples where the lower RZ fiber orientation was oblique to the cartilage-bone interface, the upper RZ fiber orientation was altered with compression to be more realigned with the lower radial organization. A linear correlation was observed between the optical birefringence and logarithm of the Young’s modulus. The results obtained in this study underscored the importance of the collagen fiber orientation and organization in modulating the articular cartilage’s mechanical responses. A detailed knowledge of such complicated structure-function relationship is critical to improve our understanding of disease progression and develop better tissue engineering based treatment strategies.