3D T1 relaxation time measurements in an equine model of subtle post‐traumatic osteoarthritis using MB‐SWIFT

The aim of this study is to assess whether articular cartilage changes in an equine model of post‐traumatic osteoarthritis (PTOA), induced by surgical creation of standard (blunt) grooves, and very subtle sharp grooves, could be detected with ex vivo T1 relaxation time mapping utilizing three‐dimensional (3D) readout sequence with zero echo time. Grooves were made on the articular surfaces of the middle carpal and radiocarpal joints of nine mature Shetland ponies and osteochondral samples were harvested at 39 weeks after being euthanized under respective ethical permissions. T1 relaxation times of the samples (n = 8 + 8 for experimental and n = 12 for contralateral controls) were measured with a variable flip angle 3D multiband‐sweep imaging with Fourier transform sequence. Equilibrium and instantaneous Young's moduli and proteoglycan (PG) content from OD of Safranin‐O‐stained histological sections were measured and utilized as reference parameters for the T1 relaxation times. T1 relaxation time was significantly (p < 0.05) increased in both groove areas, particularly in the blunt grooves, compared with control samples, with the largest changes observed in the superficial half of the cartilage. T1 relaxation times correlated weakly (Rs ≈ 0.33) with equilibrium modulus and PG content (Rs ≈ 0.21). T1 relaxation time in the superficial articular cartilage is sensitive to changes induced by the blunt grooves but not to the much subtler sharp grooves, at the 39‐week timepoint post‐injury. These findings support that T1 relaxation time has potential in detection of mild PTOA, albeit the most subtle changes could not be detected.

Post-traumatic osteoarthritis (PTOA) develops as a result of joint injuries, such as fractures, articular cartilage lesions, cruciate or collateral ligament rupture, acute meniscal tears, or a combination of these.The resulting degenerative changes in the articular cartilage persist long after the initial injury 1,2 and may trigger progressive cartilage loss. 3In the early stages of the degeneration, macromolecules of the cartilage are affected: the proteoglycan (PG) content decreases, accompanied by an increase in water content.These changes are followed by disruption of the collagen network. 3,4nsequently, changes in the biomechanical properties of articular cartilage are associated with the onset of osteoarthritis (OA). 5rly interventions, such as repairing the damaged tissue, or stabilizing the joint, can limit the progression of PTOA 6 and hence early detection of articular cartilage damage is a great asset in the prevention and treatment of joint disease. 7However, current noninvasive imaging methods 8,9 either lack sensitivity to detect subtle, early damage in PTOA, or are not yet available for clinical use.
Quantitative magnetic resonance imaging (qMRI) techniques can evaluate the compositional and structural changes in cartilage.
A common qMRI parameter, T 1 relaxation time, is dictated by the local field fluctuations at frequencies nearby the MRI resonance frequency of the water protons.The high concentration of macromolecules in healthy cartilage increases the probability of having frequencies of molecular motions near the resonant frequency, shortening the T 1 relaxation time.1][12] Native T 1 relaxation time has shown promise in the detection of articular cartilage damage in animal models of PTOA, 11,13,14 correlating well with both biomechanical properties 11,13,14 and water content, 10 and is one of the most sensitive univariate parameters for assessing cartilage degeneration. 11eep imaging with Fourier transform (SWIFT) is a technique capable of capturing signals from the most rapidly relaxing spins. 15,168][19] Multiband SWIFT (MB-SWIFT) 20 is based on multiband excitation, which allows reaching a very high bandwidth with relatively low radiofrequency (RF) power.Although both are inherently ultrashort echo time (TE) three-dimensional (3D) sequences, the use of multiple excitation bands and thus higher bandwidth in MB-SWIFT distinguishes it from regular SWIFT. 20MB-SWIFT's ability to achieve a zero TE coupled with reduced sensitivity to susceptibility artifacts and motion further support the use of the technique for quantification of T 1 in both ex vivo and in vivo conditions.Taken together, VFA-based T 1 measurement utilizing MB-SWIFT readout has a lot of potential in the diagnostics of joint-diseases such as PTOA.
However, the sensitivity of this methodology for detection of subtle changes in articular cartilage caused by mild PTOA is not known and requires testing.
The articular groove model is an animal model for PTOA that was initially used in dogs, [21][22][23] and later expanded to sheep and rats. 24,25ter, it has also been applied to the metacarpophalangeal joint of the horse. 26The close similarity between the human and equine articular cartilage makes it an interesting model to study PTOA changes in cartilage. 27Recently, the groove model was used in a 9-month study in the horse, comparing classic bluntly made grooves with sharp grooves that provoked hardly any tissue loss. 28Material from this study was deemed optimal for testing our hypothesis that T 1 relaxation time measured with VFA-MB-SWIFT is sensitive to structural and/or compositional changes caused by mild PTOA.We therefore studied the bluntly and sharply grooved samples together with controls from the horse study 28 and correlated the outcome with biomechanical and histological properties of the tissue.

| Sample preparation
In a comprehensive study exploring the groove model for equine PTOA, nine healthy adult female Shetland ponies (aged between 4 and 13 years) underwent a surgical procedure in which each pony had one radiocarpal and one middle carpal joint bluntly or sharply grooved in one randomly assigned front limb per animal. 28Three grooves in total (two in parallel in palmaro-dorsal and one in mediolateral directions) were created on the articular surfaces (Figure 1A,B) of the radial facet of the third carpal bone (middle carpal joint) and the intermediate carpal bone (radiocarpal joint) of equine subjects (Figure 1).The surgical procedure was performed under approval of the Utrecht University Animal Experiments Committee and the Central Committee for Animal Experiments (permit AVD108002015307) and in compliance with the Dutch Act on Animal Experimentation.The ponies were subjected to an 8-week incremental exercise program on a treadmill starting 3 weeks after surgery.After this, they were given free pasture exercise until euthanasia at 39 weeks after surgery.A more comprehensive description of experimental design and procedures is given elsewhere. 28subset of samples used in a previous study 29 (8 bluntly and 8 sharply grooved osteochondral samples along with 12 samples from contralateral control limbs) were imaged in this study.

| MRI
MRI experiments were carried out on two different 9.4 Tesla preclinical Varian/Agilent scanners (Vnmrj DirectDrive console v. 3.10) using the same 19 mm quadrature RF volume transceiver (Rapid Biomedical GmbH) on both.The samples were mounted on custom-made sample holders and immersed in a test tube filled with phosphate-buffered saline (PBS) containing enzyme inhibitors, with the normal of the articular surface broadly perpendicular to the main field (B 0 ), at the center of the RF coil.T 1 mapping was carried out by means of VFA measurement using 11 flip angles spanning 1°-20°( 1°-8°, 10°, 14°, and 20°).A 3D radial MB-SWIFT sequence with bandwidth of 385 kHz, TR of 2.97 ms, 16,384 spokes per flip angle, field-of-view (3 cm), matrix size of 256, 3 and 32 dummy projections before spatial encoding was utilized as the read-out sequence.VFA was selected as the T 1 relaxation time measurement method over another common technique termed Look-Locker, 30 as it has a higher signal-to-noise ratio. 31As the imaging was conducted using small tissue samples, impactful B 1 + variations (that could pose problems for the VFA technique) were not expected with the samples.

| Reference methods
Before the MRI, the osteochondral samples underwent biomechanical indentation testing, as reported previously. 29For each sample, six measurement locations (three points each on the dorsal and palmar parts of the cartilage surface) were chosen.
The calculated Young's moduli from equilibrium and instantaneous (peak force after displacement) stress/strain ratios were corrected using Hayes equation 32 and Poisson ratios of ν = 0.2 33 and ν = 0.5, 34,35 respectively.
After the MRI, PG content was measured as OD, obtained from digital densitometry (DD) images of Safranin-O-Fast-Greenstained histological sections, 29 matching the locations subjected to biomechanical testing and MRI.Similar to the procedure reported by Mohammadi et al., 29 OD analysis in this study was conducted for 12 regions of interest (ROIs) per sample, to estimate the depth-wise profiles through cartilage thickness, with a width of 1 mm using custom-made code in Matlab software (R2019, Mathworks).These 12 ROIs included the 6 biomechanical indentation testing locations, 4 locations at the grooves running in palmarodorsal direction and 2 at the intersection of grooves (Figure 1C).
The depth-wise profiles were obtained for the superficial half and for the full cartilage thickness.

| MRI data processing
VFA-MB-SWIFT images (Figure 1D) were reconstructed using an iterative gridding algorithm. 36For imaging experiments, T 1 relaxation time maps were fitted using nonlinear Gauss-Newton minimization in a voxel-wise manner (Equation 1).

∑ S T S T Sα
where S 0 is the initial signal, α j is the individual flip angle, and N is the number of flip angles.cartilage profile while leaving out the bone and PBS and then interpolated to 10 depth-wise points.For further quantitative analysis, the mean T 1 relaxation times were computed from these trimmed profiles.VOI calculations and analyses were performed using an in-house written tool in Matlab.

| MRI data analysis
Parametric surface maps were created to visualize the differences between the grooved regions and their adjacent regions.As a part of this, a triangulated 3D mesh was generated for the articular cartilage using a threshold-based multislice segmentation technique on the 6°f lip angle magnitude images.Surface normals for each triangle were then calculated for the mesh.After this, the 3D T 1 relaxation time maps were sampled from within the cartilage volume along the surface normals.The sampled volume was then eroded from the surface side to remove partial volume effects from the surrounding PBS.Finally, mean T 1 values calculated along the surface normals were displayed as a surface map on the mesh.The 3D triangular meshes for cartilage surface were computed using Materialize Mimics (Materialise NV).

| Statistical analysis
The statistical analyses were carried out using a linear mixed-effects (LME) model, with T 1 , OD (both from full-depth and half-depth regions), and equilibrium and instantaneous Young's moduli as the dependent variables.Samples from adjacent, grooved, and their intersection regions were considered all together for controls in statistical analysis, and samples from grooves and their intersection regions were together considered as grooves in LME analysis (to limit the number of possible comparisons).The ponies were selected as subject for random effects, accounting for additional covariance caused by taking observations from multiple locations within the same animal.The types of cartilage surface (control tissue, tissue at the grooves and tissue adjacent to the grooves for both groove types) were set as fixed variables, and the scanner information was considered as a covariate in the analysis of T 1 .
As only a subset of the specimens from the previous studies 28,29 were used in this study, all the data-analyses were redone.In addition, partial correlation analysis was performed to evaluate the relationships between the OD, T 1 relaxation times and the equilibrium and instantaneous Young's moduli.A p < 0.05 was considered as the limit of statistical significance.The normality test (Shapiro-Wilk) was used to check for normal distribution of T 1 and reference data, and, depending on this either Pearsons' (R p ) or Spearmans' (R s ) partial correlations were used, and approximate linearity of the dependent variables was achieved.All statistical analyses were carried out using IBM SPSS statistics (v.27 SPSS, IBM Company).

| Qualitative analysis of MRI data
Surface maps of the relaxation times in cartilage showed visual differences between the lesioned cartilage and adjacent tissue in the bluntly grooved samples.The damage was, however, not visually detectable in the samples with sharp grooves (Figure 2).
T 1 relaxation time maps revealed an increase in T 1 values between the grooved groups with respect to controls, especially at the grooved regions of bluntly grooved samples (Figure 3A).Changes at groove regions in T 1 followed the trend of changes in OD maps reflecting PG content, showing decreased PG content at the grooved regions especially with the blunt groove type (Figure 3B).
Visually evaluated profiles showed elevated T 1 values at the grooves and at groove intersections with respect to controls (Figure 4A).
Compared with controls, T 1 values of grooved samples appeared higher towards the superficial half of the cartilage at the grooves and their intersections.The T 1 relaxation times from both groove types overlap with those of the controls towards the deeper zones in depth-wise profiles.Samples with blunt groove injury showed higher T 1 values than sharply grooved samples in the depth-wise profiles at the grooves and intersections of those, for the superficial half of the cartilage (Figure 4A).
However, the OD profiles showed a decrease of PG content for each anatomical location (adjacent to grooves, at grooves and groove intersections) with respect to controls (Figure 4B and Table 1).

| LME model analysis
Higher T 1 values were noted in bluntly and sharply grooved regions than in the respective adjacent areas in the following interaction groups: between the blunt grooves versus the adjacent regions, sharply grooved Photographs of the specimens (A) and the corresponding surface T 1 relaxation time maps (B).Cartilage surface maps calculated with a segmented mesh could differentiate between the defective regions from bluntly grooved regions and their adjacent regions, highlighted with black arrows.Sharp grooves were visually almost indistinguishable from the adjacent regions.
versus adjacent regions of bluntly grooved groups, adjacent regions of sharp grooves versus blunt grooves, and between adjacent regions of both bluntly and sharply grooved groups (Table 2).With respect to controls, significant differences with increased T 1 values were noted with both groove types: with the adjacent regions of bluntly grooved groups, with blunt grooves, and sharp grooves.In the superficial half of the tissue, there were larger and more significant differences in T 1 relaxation times than in full-depth cartilage, that is, between the bluntly grooved group and the adjacent regions of the sharply grooved group.With respect to controls, different significance levels were noted with adjacent regions to bluntly grooved group and with sharply grooved.In addition, significant differences were noted between sharp grooves and their adjacent regions (Figure 5A and Table 2).
The differences in T 1 relaxation times were mostly in line with PG content as measured by OD.However, with OD, the differences were noted in all regions of the grooved (blunt and sharp) groups compared with the corresponding regions in the controls, in both full-depth and superficial half-depth tissue.In addition to the differences noted between groups with full-depth T 1 , OD showed also significant differences between control versus adjacent regions of sharp grooves, sharp grooves versus their adjacent regions, and blunt grooves versus sharp grooves.Compared with superficial halfdepth T 1 , there were additional significant differences between controls versus adjacent region of sharp grooves and sharp grooves versus blunt grooves (Figure 5B and Table 2).Biomechanical properties differ only between control and adjacent regions of sharp grooves (Table 2).

| DISCUSSION
In this study, we investigated the potential of T F I G U R E 4 T 1 relaxation times (A) and OD profiles (B) sampled and depth-normalized from the adjacent locations that also served as biomechanics testing locations, groove locations, and groove intersections.T 1 profiles from groove and their intersection locations indicate increased T 1 values in both groove types in the superficial half-depth with respect to controls, but comparatively higher in blunt grooves in the superficial part of the cartilage.OD profiles showed larger variations in the superficial half of the cartilage layer.
T A B L E 1 Mean and 95% confidence intervals of each group for T 1 relaxation time, OD, and mechanical properties.T A B L E 2 Mean differences and their significances.

Comparisons between groups
Full-depth T 1 (ms)  Higher T 1 values are broadly indicative of higher hydration of the cartilage tissue 10 due to, for example, more significant degenerative changes in the respective locations.In the adjacent regions, the changes were not large enough to cause a detectable increase in the T 1 relaxation time.For full-cartilage thickness, slightly smaller differences were noticed F I G U R E 5 Boxplots of T 1 relaxation times (A) and OD (B) from full-depth and from superficial zone, from all 12 measurement locations (6 adjacent to grooves, 6 on the grooves + intersections).Red lines are the median values of the measurements and red crosses indicate the outliers within the data limits.Star diagrams above and below illustrate the significant differences identified between the compared groups.BA, blunt adjacent; BG, blunt groove; C, controls; SA, sharp adjacent; SG, sharp groove.
T A B L E 3 Correlation coefficients obtained from the partial correlation analysis.

Correlation with equilibrium modulus Correlation with instantaneous modulus Correlation with OD
Anatomical location and scanner information as covariate | 2663 between the compared groups than in the superficial 50% of the tissue.
This was not unexpected since the damage caused by the grooves was most prominent in the superficial half of the cartilage and it should be noted that the sharp grooves were limited to the upper 400 μm of the cartilage.8][39] The correlation of T 1 relaxation times with the DD-measured PG content of cartilage and the detected differences in T 1 between the grooved regions and controls or adjacent-to-damaged regions indicate that PG loss of cartilage or cartilage in general had an influence on T 1 values.These correlations of T 1 with the reference methods were weaker than in previous studies on naturally degenerated human cartilage, 13 enzymatically treated articular cartilage, 11 and in an equine defect model, 14 which showed correlations between R = 0.5 and 0.8.In this study, the correlations of T 1 with the equilibrium moduli were higher than those with OD, which is consistent with previous studies. 11,13,14This outcome potentially indicates that the T 1 relaxation time in this groove model was also influenced by other degenerative processes, not just by the increase in free water content caused by PG loss at groove regions, although with such weak correlations, strong conclusions cannot be made.The stronger correlations for the T 1 relaxation times in the superficial half of the cartilage than in the fulldepth are in line with the superficial location of the affected region.
Previous studies on articular cartilage 18,19 have shown the feasibility of measurement of T 1 relaxation times for a broader range of spin populations with VFA-SWIFT, and demonstrated increased T 1 values in trypsin-treated cartilage samples.As the trypsin treatment induces PG loss in articular cartilage, the increase in T 1 relaxation at the grooved regions with reduced PG content in this study is in line with the previous results.39,40 Results of the current study are consistent with the previous studies in the sense that T 1 increase coincided with decreased PG content.However, the correlation found between the PG content and T 1 relaxation time was not as strong as reported in most of the previous studies.This could be due to the different methods of measuring T 1 but is more likely due to the very modest degree of degeneration in the present sample pool (especially in the sharply grooved group), with insufficient data spread for reliable correlation analysis.In a previous report, 19 it was noted that VFA-SWIFT has lower apparent T 1 relaxation time values (T 1 ranging between ~800 and 1800 ms) compared with measuring T 1 relaxation by inversion or saturation recovery sequences (T 1 between ~600 and 2400 ms).This shift in T 1 values and compression of the range of expected T 1 values could be one of the reasons for the observed low contrasts in the cartilage between the grooved and adjacent regions.
As reported in the previous study on the same equine model 29 with a slightly larger sample set than in this study, DD revealed degenerative changes in grooved and the adjacent regions.However, the biomechanical indentation testing in the previous study, conducted at the regions adjacent to the grooved tissue, showed a larger decrease in the sharply grooved than in the bluntly grooved group with respect to controls.The differences in observations between this and the previous study likely result from the different number of samples used in the study.
The current study is not without limitations.First, samples were immersed in PBS during imaging and thus the differentiation between the cartilage and background was not straightforward and may have resulted in inclusion of PBS or exclusion of the most superficial part of cartilage in the analyzed data.However, for cartilage segmentation from the PBS, careful evaluation of the depth-wise profiles was undertaken.
Second, the exact biomechanical indentation testing sites and sites used for histology were not coregistered with the sites used for T 1 .However, careful and consistent visual evaluations were followed to identify the locations of the biomechanical indentation testing and histology, which were then used for T 1 relaxation time analysis.Third, the sham-operated contralateral joints were considered as control samples, but there was no separate control group from non-operated horses.The sham-operated contralateral joints cannot be considered entirely free from structural alterations of the cartilage. 28,29In particular bilateral inflammatory responses might have been induced because of repeated arthroscopies and arthrocentesis. 41Fourth, due to technical constraints the full set of the samples from previous study could not be included and was not balanced between the control (n = 12 samples) and the grooved groups (sharp = 8, blunt = 8), which might have reduced the statistical power due to the asymmetric comparison.However, the relatively small sample sizes in the grooved groups did not significantly impair this study, as the goal was to investigate the feasibility of T 1 relaxation time mapping in detecting chronic changes due to lesions in grooved groups with respect to controls.Fifth, two different MRI scanners were utilized in this study due to hardware issues, which imposed an additional variable to be considered in the statistical analysis.Lastly, no B 1 + correction was applied in the VFA-based computation of T 1 relaxation time.However, the B 1 + field of the volume coil used in the study was carefully calibrated and was highly uniform within the relatively small specimens.Thus, the effects of B 1 + deviations on the T 1 measurements are expected to be minor and consistent between the samples.
In conclusion, this study reports the potential of T

Depth-wise mean
profiles of the T 1 relaxation time values from the grooves, regions adjacent to the grooves and at their intersections were calculated in cylindrical 3D volumes-of-interest (VOIs) (width ≈ 4 voxels and height ≈ 20 voxels (depending on cartilage thickness), covering the whole cartilage layer and some parts of background and bone).Similar profiles were also obtained from the control samples.A total of n = 12 VOIs were defined per sample, of which n = 6 VOIs were placed in the adjacent regions (carefully matched with the biomechanical testing points) and n = 6 VOIs were placed directly on the grooved regions (four on groove lines between the mechanical testing points and two on the intersections of the groove lines) (Figure1C).A similar spatial pattern of VOI placement was used for the contralateral control samples.Before the analysis, the depth-wise T 1 profiles from the cylindrical VOIs were trimmed to obtain just the F I G U R E 1 Schematic drawing of the right equine carpal joint (A), showing the grooved sites on the intermediate carpal bone in the radiocarpal joint and on the radial facet of the third carpal bone in the middle carpal bone with pictures of the corresponding surfaces where the grooved areas are indicated with dashed circles (B), drawing of a grooved joint surface showing the groove pattern (red lines) with biomechanical testing locations (red dots), and all 12 measurement points (six adjacent to grooves [red dots] and six on grooves [blue and green dots]) and pictures of typical examples of osteochondral samples with grooved surfaces (C), orthogonal slices through a three-dimensional (3D) volume averaged from 3D images acquired at flip angles 5°, 6°, 7°, 8°, 10°, and 14°.The arrows point at groove locations (D).
1 relaxation time, measured with VFA-MB-SWIFT, for the assessment of PTOA induced by surgically created grooves on the articular surfaces of the carpal joints in ponies.T 1 relaxation times were analyzed for changes, along with the reference properties of the articular cartilage measured via mechanical indentation testing (equilibrium modulus and instantaneous modulus) and DD (PG content).The findings indicated that the T 1 relaxation time of cartilage was altered due to the changes caused by the injuries and was dependent on the severity of the induced damage.Particularly, T 1 was sensitive to the F I G U R E 3 Coronal axis T 1 relaxation time maps (A) of seven consecutive slices averaged through the groove points and the corresponding optical densitometry measurement maps (B) from blunt and sharp grooved articular cartilage (groove locations marked with black and white arrows) and a contralateral control joint (with intact articular cartilage surface).Between the groove types, blunt grooves are well detectable in comparison with sharp grooves on the articular surface of the carpal joint.
more advanced post-traumatic degeneration caused by the blunt grooves in the articular cartilage.The very early post-traumatic degeneration in the adjacent tissue to the sharp grooves was hardly distinguishable from controls based on the T 1 relaxation times.Tissue damage was easier to detect in the superficial 50% of the cartilage compared with full-depth cartilage.The study revealed weak negative correlations between the T 1 relaxation times and cartilage PG content and mechanical properties.Blunt grooves were well detectable in comparison to the sharp grooves from the average T 1 relaxation time maps in slices through the groove points and in the surface map visualizations by their higher T 1 relaxation times.
−0.218* Note: Correlation coefficients obtained from the partial correlation analysis between T 1 relaxation time and the reference properties of articular cartilage for the two depths analyzed.*Statistically significant (p < 0.05).PALA ET AL.
1 relaxation time mapping in tracing the progression of modest PTOA in articular cartilage in surgically grooved tissue and its adjacent regions.The model utilized two different types of grooves featuring considerable differences in local tissue damage and appeared very suitable for the purpose of the study.The changes induced by the surgical model varied from modest in the blunt grooves to minimal in the sharp grooves.The T 1 relaxation times evaluated in small VOIs were highly sensitive to damage in the superficial half of the cartilage, especially at grooved regions in the bluntly grooved samples compared with healthy cartilage.The findings of this study suggest that T 1 as measured by VFA-MB-SWIFT is sensitive to modest tissue changes caused by PTOA, but cannot distinguish the mildest forms of damage in the adjacent tissues, which are also undetectable via biomechanical testing.Nissi, and Olli Nykänen were involved in MRI measurements, Ali Mohammadi, Nikae C. R. te Moller, Harold Brommer, P. René van Weeren, and Janne T. A. Mäkelä were responsible for indentation testing.Mohammadhossein Ebrahimi and Isaac O. Afara were responsible for OD measurements.Nikae C. R. te Moller, Harold Brommer, P.René van Weeren, and Janne T. A. Mäkelä were involved in animal experiments.Swetha Pala was responsible for data analysis and manuscript drafting.Swetha Pala, Mikko J. Nissi, and Olli Nykänen were involved in interpretation of results and critically revised the manuscript up to final version.All authors were involved in revision for intellectual content and final approval.