A. Data Collection
Hybrid SPECT/CT images (99mTc-HDMP using a Symbia T16 scanner, Siemens, GmbH, Germany) of 20 patients (Age: 60±10; BMI: 29±5 kg/m²) slated to receive a medial UKA prosthesis (Oxford partial knee, Zimmer-Biomet, Warsaw, USA) were acquired before surgery (20±12 days) and one year following surgery (369±16 days). All surgeries were carried out by the same surgical team. Acquisition protocols included SPECT voxel size of (x/y/z) 4.79/4.79/4.79 mm and CT voxel size of either 0.98/0.98/3 mm or 1.27/1.27/3 mm. Lower CT resolution was used for pre-operative scans of three patients and post-operative scans of six patients.
B. Data Processing
Throughout this paper, the coordinate system used, x, y and z are parallel to the patient’s frontal, sagittal and transverse axis respectively. The SPECT/CT scans were subjected to three operations within SCreg, all performed in MeVisLab (MeVis Medical Solutions AG, Germany)– pre-processing, registration and SPECT-based statistical analysis of osteoblastic activity (Figure 1). Each SPECT/CT scan had a common coordinate system for both imaging modalities; hence, same transformations were applied to SPECT scans during registration as those applied to the associated CT scan, without any disruptions in the alignment.
To account for resolution differences between image datasets, all CT and SPECT scans were resampled to a voxel size of 0.98/0.98/3 mm. Right legs were reflected about the sagittal plane to obtain a ‘left’ laterality for all specimens. To consider all possible changes in joint position, images were divided in three volumes of interest (VOIs)—femur, tibia and patella.
The cuboidal VOIs were defined in CT space using bony landmarks shown to be accurately identifiable (33) as well as other distinctive bony landmarks previously described (34) on the femur (the fibula head, the femoral trochlea proximal (FTP), the femoral lateral epicondyle (FLE) and the femoral medial epicondyle (FME)), tibia (femoral knee center (FKC) and tibia lateral and medial peaks) and patella (proximal, distal, lateral and medial patella margins). Within each VOI, an intensity threshold was used not only to isolate the bony tissue, but also to mask the prosthesis in post-operative scans thereby avoiding complications during registration.
The registration process was developed using the Elastix toolbox (25, 32, 35–36), developed specifically to perform image registration. However, despite including all the required parameters, their values must still be defined by the user. To account for differences in positioning inside the scanner, an affine registration (AR) was used for intra-subject registration. In contrast, inter-subject registration, which required the consideration of additional anthropomorphic differences, involved AR followed by non-rigid registration (NRR) with additional non-affine transformations. Registration parameters (Table 1) were tuned by registering five randomly selected femoral bone images with typical values from literature (25, 32, 35). The selected parameters were the one preserving the image integrity upon visual inspection and maximizing the Dice similarity coefficient (DSC) (25, 35, 37–38).
For intra-subject comparison, registration was performed between images of the same patient but at different timepoints. The reference image, further referred as fixed image, was the post-operative scan, with masked-out prosthesis, and the image transformed to the fixed image, further referred as the moving image, was the pre-operative one. For the inter-subject registration, the fixed image was the post-operative bone image of a reference subject randomly selected amongst the 14 post-operative scans with a higher resolution and the moving image was
the post-operative bone image of any other subject. Prosthesis masks were not required for inter-subject registrations, since the same prosthesis was used for all subjects, hence no intrinsic significant differences requiring a masking. In addition, it allows to take into account the extrinsic differences in terms of implantation. The transformations from inter-subject registration were concatenated with those from the intra-subject registration.
To increase the chances of successful registration, a multiresolution approach, aimed at simplifying the data by iterative image smoothing, was implemented (35). This approach registered large and dominant structures before moving on to progressively smaller structures. Three resolution stages were implemented for AR, and four for NRR, using a common gaussian pyramid blurring the image. Smoothing factors were defined for every stage and direction (Table 1), with each factor representing lowered resolution in the x, y and z directions of the images respectively.
A B-Spline representation (39), modeled as a weighted sum of B-Spline basis functions placed on a uniform control grid, was used to define NRR for inter-subject registration. The resolution of this control grid defined the flexibility of the transformations, starting with a larger grid for the first resolution stage and progressively decreasing to end with a thinner grid, allowing for more deformation of the image, for the last resolution stage.
Normalized mutual information (40–41) was chosen as the metric for image similarity. This metric required the computation of the joint histogram based on the dynamic range of the images. Following parameter tuning, 16 bins were defined. For NRR, since the presence of the prosthesis increased the dynamic intensity range of the images, the number of histogram bins was greater than that used for AR.
To find the optimum of the mutual information-based cost function, an iterative optimization was performed using a robust and adaptive stochastic gradient descent method (Robbins-Monro) (42–43). This method reduced the computation time by using a small subset of points from the fixed image at each iteration to compute the derivative of the cost function. Maximal number of iterations was used as a stop criteria for the optimization since it is the only criteria available in Elastix and its value for AR and NRR was defined during parameter tuning.
A sampler was used to select the subset of points required for the optimizer. It involved a random selection of points at each iteration for the stochastic optimization, here based on the Halton sampling, thereby improving the smoothness of the cost-function and avoiding two known problems related to the use of mutual information, “overlap problem” and “grid effect” (44).
Interpolation of the gray values between voxels was required owing to the selection of random points by the sampler, which did not necessarily ensure that the corresponding points in the moving image were at a voxel position. In case of the AR, nearest neighbor interpolation was selected, with a linear interpolator to generate the last deformed image (29), while for the NRR, a linear interpolator was selected, with a third degree B-spline interpolator to generate the last deformed image.
Once scans were registered, SPECT values were normalized to account for subject-specific differences in tracer metabolization. For each scan, average SPECT activity calculated within a rectangular volume of 6008.5 mm³ (2100 voxels) at a distance of 10 cm from the FKC along the z axis was used for normalization. Pre-operative SPECT activity was subtracted from post-operative values for each voxel, and an aggregate map was created by averaging differences within each voxel value over all co-registered subjects, thereby providing a unique representation of all the subjects.
C. Validation
Since no unique tool allows the complete validation of a registration (35), multiple methods were employed, including qualitative criteria for preliminary evaluation and quantitative criteria for objective evaluation. The registration was considered successful only if all criteria were fulfilled.
All images were first qualitatively evaluated for image integrity, including the presence of image blur, irregularities, such as holes, and loss of bone contour shape. This was followed by qualitative assessment of the conformance of bone contours on a two dimensional overlay of the two registered images over a few horizontal slices of CT selected in the region with the prosthesis, since this region was expected to be exposed to the most errors in registration.
Following qualitative assessment, volumetric agreement between the registered images was quantified by computing the DSC and used as quantitative validation. Considering the resolution of the CT and SPECT images, complexity of the structures, and the presence of an implant, a DSC over 80% was defined as the criterion for successful registrations (25). In addition, DSC were also acquired before and after intra-subject and inter-subject registration, in order to observe if registration significantly improved spatial alignment. Since the data did not follow a normal distribution, as assessed by a Kolmogorov-Smirnov test of normality, a Wilcoxon signed rank sum test (Matlab R2016, MathWorks Inc., Natick, USA) was performed to analyze the effect of registration on DSC, and a Mann Whitney U test was applied to evaluate the impact of varying resolution of pre-operative and post-operative CT images on DSC. P-values below 0.05 were considered statistically significant.
The difference between coordinates of anatomical landmarks in registered images was used as the second criterion for quantitative validation. Considering the ratio of the SPECT resolution over the CT resolution, a difference lower or equal to (2;2;0) voxels was considered as satisfactory to assure that the registered points were associated to the same SPECT voxel.
SPECT activity was used as the final criterion for quantitative validation. Apart from comparisons with literature (14), the location of main increased activity on the aggregate map was also compared with an ongoing study on the same population as this project using the regional classification system for analysis (18). In addition, since osteoblastic activity is directly linked to bone strain as explained by the Wolff’s law (45), SPECT activity was compared with the location of increased bone strain following UKA (46).