Statistical shape model of the talus bone morphology: A comparison between impinged and nonimpinged ankles

Abstract Diagnosis of ankle impingement is performed primarily by clinical examination, whereas medical imaging is used for severity staging and treatment guidance. The association of impingement symptoms with regional three‐dimensional (3D) bone shape variaties visible in medical images has not been systematically explored, nor do we know the type and magnitude of this relation. In this cross‐sectional case‐control study, we hypothesized that 3D talus bone shape could be used to quantitatively formulate the discriminating shape variations between ankles with impingement from ankles without impingement, and we aimed to characterize and quantify these variations. We used statistical shape modeling (SSM) methods to determine the most prevalent modes of shape variations that discriminate between the impinged and nonimpinged ankles. Results of the compactness and parallel analysis test on the statistical shape model identify 8 prominent shape modes of variations (MoVs) representing approximately 78% of the total 3D variations in the population of shapes, among which two modes captured discriminating features between impinged and nonimpinged ankles (p value of 0.023 and 0.042). Visual inspection confirms that these two shape modes, capturing abnormalities in the anterior and posterior parts of talus, represent the two main bony risk factors in anterior and posterior ankle impingement. In conclusion, in this research using SSM we have identified shape MoVs that were found to correlate significantly with bony ankle impingement. We also illustrated potential guidance from SSMs for surgical planning.


| INTRODUCTION
Ankle impingement syndrome can be described as abnormal entrapment or contact of structures in anterior or posterior areas of ankle joint during dorsiflexion or plantarflexion of the foot and can cause pain, stiffness, and limited range of motion. [1][2][3] Based on pathoanatomic features developing ankle impingement, it is divided into bony ankle impingement (also known as ankle impingement exostosis 4,5 ) and soft tissue ankle impingement. 6,7 The mechanisms that create obstruction and give rise to pain in both anterior and posterior bony ankle impingements are driven by abnormalities in the morphology of ankle bones, mainly osteophytes (bone spurs) in the tibia and talus. 6,8,9 In the anterior bony ankle impingement, the pain occurs during dorsiflexion of the foot, when osteophytes on the talar neck or distal tibia entrap the soft tissues in the anterior part of the ankle joint or when these osteophytes undergo irritation because they collide with each other or receive an external force. 2,10 In posterior bony ankle impingement, on the other hand, obstruction and pain happen when during plantarflexion of the foot some soft tissues get squeezed in the area between the tibial rim or calcaneus and the bony excess of the posterior side of talus, including os trigonum and Steida's process, or when an action like a plantarflexion or an external force affects the posterior process of the talus. 2,9,11 Diagnosis of bony ankle impingement is focused on the size and location of these bony processes or osteophytes. 2,10 In addition surgical treatments of ankle impingement aim to resect these osteophytes and bony excesses to improve functionality and relieve pain. 6,12 Therefore, having a good understanding of the morphology of these bony abnormalities is vital for severity detection as well as preoperative planning of bony ankle impingement.
In clinical examinations, simple palpation may provide detection of prominent osteophytes or bony excess, but it engages a high rate of false-negative results. 2,13 This necessitates further support from medical imaging modalities like radiography and computed tomography (CT). 2 Radiographs, while useful for diagnosis of bony ankle impingement, provide a reduced understanding of bone anatomy due to their two-dimensional (2D) nature and also might not always be sufficient for the assessment of bone abnormalities in their standard views due to the overlapping phenomenon of tibia and talus. 2,14 CT scans, by comparison, allow us to render the image in three dimensions (3D), restrict the visibility to only bones, look at bones from any angle, and so CT scans can help clearly distinguish osteophytes in bony ankle impingement. 10,15 Classification systems have been proposed for bony ankle impingement that work mainly based on the location and size of bone spurs in talus and tibia using radiographs. 6,12 Quantification of specific visible bone shape abnormalities are also usually performed by means of some predefined measures. Berberian and colleagues documented size measurements of tibial and talar osteophytes on plain radiographs and axial CT images of nine anterior ankle impingement patients. 7 A more recent study documented the coordinate location and size of the osteophytes from the cadaveric talus and tibia of a population with bony ankle impingement syndrome. 8 Diagnosis and quantification using above methods is focused on predefined categories and measures, 16,17 hence these measurements call for substantial manual work, are prone to errors and also provide a limited description of total bone anatomy, ignoring minor morphological variations that are difficult to capture but might play an important role. 18,19 This motivates the use of more automated and exhaustive quantification methods that can take into account the whole 3D anatomy of the bone instead of a set of limited predefined measures.
Using statistical shape modeling (SSM), a sophisticated technique that facilitates grouping and characterization of coincident shape variants into a set of discrete components, 20 it is possible to quantitatively describe the total morphology of bones. [21][22][23][24] An SSM describes all variations in shape that are present in the study population, so it can extract and quantify bony shape variations that potentially give rise to impingement.
Surgical treatment planning of bony ankle impingement also can benefit from SSMs, by gaining insights into the area and the amount of osteophytes that might need resection. Inadequate resection or overcorrection of osteophytes have been discussed as factors causing surgery inefficiency, the need for revision surgery, and painful debilitating reactions. 15,[25][26][27] Branca and co-authors 28 clarified the results of arthroscopic treatment of anterior ankle impingement and reported that in patients with more severe ankle impingement there were more cases of osteophyte recurrence following the surgery, probably as a result of under-resection of talar or tibial osteophytes and degeneration. 15 Evaluation of the area for resection is sometimes based on predefined measurements or idealized assumptions about the bone shape. 15,29 Frigg and co-authors described a surgical procedure of posterior ankle impingement. They measured the preoperative size of talus os trigonum and posterior talar process from radiograph by assuming it as the bony area lying posteriorly outside the talar circle in radiographs. They evidenced several painful stress reactions after surgery. 27,30 SSMs, in contrast, can be applied to extract the bone morphology and estimate the normal bone shape for a specific patient without idealized assumptions about underlying geometry. 26 SSMs have also been recently investigated for talus implant design in surgeries like total ankle replacement. 19 According to a very recent systematic review there is no published research comparing the imaging finding of people diagnosed with ankle impingement to asymptomatic people. 31 In this study we used SSM on 3D reconstructions of talus bones from CT images to quantitatively formulate the discriminating talus bone shape variations of a group of patients with bony ankle impingement from a nonimpinged control group. Finding such discriminating imaging features can be useful in quantitative diagnosis of bony ankle impingement. We have also illustrated possibility of SSMs assisting surgeons in preoperative planning of bony ankle impingement.

| METHODS
The main steps involved in this cross-sectional case-control study for the generation of 3D SSM of the talus are presented in Figure 1 Scalismo uses Gaussian Process models (GP) for modelling shape variations. 33

| Segmentation of bones
Using the commercial segmentation software (Mimics 20; Materialise), talus bones have been segmented in three dimensions from each CT scan. Some of the steps in semi-automatic segmentation are thresholding, labeling bones, region growing, filling holes inside segmentation masks and smoothing, among others. Using the same software, triangulated bone surfaces were extracted from the segmentation results. Next, triangulated bone surfaces of segmentations got exported as volumetric meshes. As the dataset used in this study was a mixture of left and right side ankles, right side ankles were mirrored in the sagittal plane. The dataset consists of 33 bone surface shapes, including 16 talus shapes with impingement diagnosis and 17 control talus shapes with no impingement diagnosed. In later steps SSM of talus bones will be generated from this set of 33 training shapes.

| Reference selection
To capture the shape variations, in the process of SSM generation ( Figure 1), first a "reference" shape is defined, to which all the other shapes of the training set will be registered. In case the reference is one of the training shapes, the SSM may become biased to that specific shape. To avoid such bias, an iterative method 35 is used to generate a shape that can act as an objective reference. In this method one of the shapes has been used as the reference in the first step, and then the shape modelling pipeline is undergone. The mean of the generated shape model is used as the reference for the next step. This process has been performed until the there is no significant change in the generated talus mean (less than 0.01 mm).

| Alignment of bones to the reference
In order to capture the true anatomical shape variations inside the shape population, it's necessary to eliminate the non-anatomical variations in the training set. Hence, all variations that are the effect of rigid body rotations or translations as well as size differences among bones should be eliminated as much as possible. This process is called alignment.
In this process, for a list of landmarks placed on all shapes, the best transformation T X ( ) parameters that can transform the reference landmarks (Ω ref ) to each target shape landmarks (Ω target ) are found using a least square estimation method in, 36 which leads to an aligned dataset of shapes.
For this study, first, 6 anatomical landmarks ( Figure 2) have been manually placed on each talus at lateral tubercle, medial tubercle, lateral and medial aspect of talar dome, lateral process, and head. Then the landmark alignment step is undergone to generate a rough alignment of all the bones to reference. Without this first landmark alignment step, chances are that the next step alignment gets trapped into local minimum.
In the second round of alignment, 1000 points are taken uniformly from the surface of each talus shape as landmarks, and used to improve the previously achieved alignment.

| Correspondence establishment
Setting up correspondence between points of the reference shape and each target shape one by one such that the morphology of the target shape can be described in terms of the reference, is called the correspondence establishment problem. 37,38 The goal is to find a deformation field from a space of deformation fields, which map each point of reference to its corresponding point on target surface: where Γ R is the reference, Γ T is the target surface, u is a spatial transform that maps points from surface Γ R to homologous points on Γ T from a space of deformations, F. In this equation the D is a similarity measure and R is a regularizer. 37 In this study, a shape model fitting method 37 is implemented using Scalismo, in which admissible deformations are modelled as a GP, then a F I G U R E 2 Surface of a talus shape from training set with six landmarks used for alignment low-rank approximation is computed, and the correspondence establishement problem is formulated as a parametric optimization problem.

| Extraction of bone shape MoVs
The covariance matrix of the data vectors containing the ordered coordinates of the corresponding points is established across all talus shapes. Principal component analysis was then performed on the covariance function, resulting in the shape MoVs with the direction of change.
New sample or existing talar shape (x) can be represented, using the mean shape of the talus (x ) and a weighted sum for MoVs of the talus (Φ) and their corresponding coefficient, representing the magnitude of the MoV for that specific shape given by the vector b (Equation (3)), as shown in step 5 of Figure 1: Where the b values describe the contributions of the first c modes of shape variation to the bone shape x. In other words, a shape parameter (e.g., b 1 ) for a given mode of shape variation (e.g., Φ 1 ) describes how far a shape x is away from the mean shape x in the specified direction of Φ 1 . For the tests in the following sections we sampled the shape model to create new shapes or estimate shape parameters for an existing shape as in Equation (3).

| Comparison of bone shapes between groups
While the mean talus shape, x , and the MoV, Φ, are identical by definition for all taluses of the training set regardless of their group (patients with impingement vs. controls), the shape parameters (i.e., b values in formula 3) are different. Therefore, the shape parameters for the nonspurious modes of shape variations (suggested by parallel analysis) were compared between the two groups using a t test and thereafter p values were reported to identify the most distinguishing shape modes between groups. For this comparison, we used results from compactness test (explained in the next section) to select the most significant MoVs. We also used Horn's parallel analysis 39,40 as a more objective way to identify the number of significant (nonspurious) PCA modes of shape variations, to include in later analysis. This method compares the eigenvalues of the training shapes with eigenvalues from a random dataset of the same size and distribution of training set to identify the number of MoVs that capture true variations than noise.

| Evaluating performance of model
To evaluate the quality of generated shape model, four wellknown shape model evaluations tests have been performed. [41][42][43][44] Generality test evaluates how well the constructed shape model can represent an unseen shape of the same class as the training set. Specificity test on the other hand, assesses if any generated shape of the model belongs to the same class as the training set.
Compactness test of the shape model by measuring the amount of variability that the model captures within each MoV tries to capture how compact the SSM can be. We have also used a test as discussed by Van de Giessen et al. 21,42,44 to evaluate if the number of used shapes is sufficient for creating an accurate shape model for the purpose of this study.
Algorithm of these tests are depicted in Table 1  In this step, we create an SSM from the group of control talus bones to capture the overall shape variations of nonimpinged taluses in samples.
By fitting this control group SSM to a talus bone from impingement group and representing it as a distance color map we will assess if this could facilitate for the purpose of illustrative and quantitative surgical planning.

| Web-based interactive shape model visualization tool
In order for the readers of the paper to easily interact with the In the implementation of this tool, we've made extensive use of the Javascript, Python and Scala programming languages, and especially Node.js, VTK.js, Flask and Scalismo libraries.

| RESULTS
Following the steps in the methods section we generated the SSM from the group of 33 talus shapes as training set. In the reference selection process, distance between the two consecutive means was changing less than 0.01 mm at around 0.16 mm. The average distance of the shapes from the mean of the generated SSM was 0.62 (±0.1) mm, in terms of mean and standard deviation.

| Model evaluation
The result of applying the four model evaluations is illustrated in Figure 4.
Generality test results show the improvement of shape reconstruction accuracy with increasing the number of shape modes.
After around 18 shape MoVs, the surface reconstruction error reaches below 0.57 mm which is basically the average resolution of the CT images, while including at least 95% of data variance.
The compactness test shows the amount of data variance and cumulative variance captured in the model, by increasing the number of principle components. In compactness chart we also depicted the eigenvalues for simulated dataset from Scree plot of parallel analysis.
F I G U R E 3 Basic architecture of the WebSSM 1 : Parameters of the shape model generated from a population of shapes are 2 extracted into CSV files. 3 CSV files are read in the backend of the web app, which serves it for 4 the client that is using the app from her web browser. SSM, statistical shape modeling According to parallel analysis the number of significant shape modes to retain for later analysis is 8, since for all modes above mode 8 the variance taken by randomized data is more than the variance taken by real data. The model specificity shows a value of 0.57 (±0.04).
The results of the test for sufficiency of training set size, show that as the number of training shapes increased the SSM of talus could more accurately describe the shape of bones that were not included in the training set. The mean estimation error is below average voxel size of the images (0.57 mm) towards the end of the curve (using more than 30 shapes). We can conclude that the amount of training shapes in our database seems sufficient to represent the population described by the model.

| General bone dominant shape variations
As shown in Section 3.1, the results of parallel analysis suggests the first 8 shape modes to be the most nonspurious shape MoVs according to their eigenvalues. These 8 shape modes capture around 78% of total shape variations. The amount of shape variations captured by each of this set along with the p values of t test examining the discrimination of the specific shape mode between the control and impingement group of bones is shown in Table 2.
As Table 2 shows, the most prominent shape modes that distinguish between the groups of control and impingement are shape Modes 1 and 3. Shape Mode 8 also shows a lower p value, but it is only responsible for 3% of the whole variations in the population of shapes. For the rest of the analysis we focus on the two discriminating shape modes.

| Physical interpretation of bone shape variations
The two prominent shape modes (Modes 1 and 3

| Comparison of bone shape variations between nonimpinged and impinged groups
Hypothesis testing between the two groups of nonimpinged and impinged groups with the null hypothesis of "the distribution of this shape mode is not significantly different between the two groups," is depicted in Figure 6 with p values. For shape Modes 1 and 3 a   F I G U R E 5 Visual description of bone shape variations for the first and third shape modes of variation. Red arrows show the areas that are mainly affected in each MoV. Red curves also are placed to better visualize the areas of changes between 3 SD distances from mean. MoV, modes of variation respectively, which is a key factor that can give rise to impingement.
Examples of 3 talus bone shapes of the population in this study, with shape Mode 1 around −3 standard deviation, 0 and +3 standard deviation from the mean shape are illustrated in Figure 7.
The other prominent observation was in shape Mode 3. In this shape mode while the talus of the nonimpinged controls deviated positively from the mean talus shape, those of the impingement patients had negative deviations (Figure 8). In the impingement group, the anterior process of the talus bone shows some form of a boney excessive. On the contrary, the healthy subjects had less anterior boney extension and were showing more smooth anterior parts i.e., positive deviation from the mean talus shape Mode 3, in Figure 8). Figure 9 shows one shape of the impingement group and its corresponding shape generated by fitting the SSM of the control group to this shape.

| Surgical planning
A color map of distance between a shape of impinged group from its fit of the control group SSM, is also illustrated in Figure 10.
Plots showing the distribution of MoV 1 and 3, each separated by the groups of control (blue) and impinged subjects (orange). Box and whiskers represent the median and interquartile range for each MoV of the two groups. Vertical axis is standard deviation. MoV, modes of variation F I G U R E 7 Talus bone of three subjects in this study: middle, subject whose value for coefficient of first MoV is the closest to the mean talus, right: a shape from control group, and left, shape from impingement group. MoV, modes of variation ARBABI ET AL.

| 191
As Figure 10 depicts, the areas colored in red by color map or deformation field are a potential site for debridement during surgery.

| DISCUSSION
In this study, a 3D Statistical Shape Model was generated to compare the shape of the talus bone between a group of patients with ankle impingement syndrome and a group of nonimpinged controls.
This also made it possible to quantitatively describe the complex 3D geometry of talus bone in terms of modes of shape variations and to identify the shape modes that show significant correlation to increased risk of ankle impingement in the population. We also illustrated how making an SSM of a group of healthy bones can be used to quantify disease severity and surgical planning.  shape modes, as depicted in our generality test. The first prominent shape mode in their shape model captures the effect of size of shapes since this F I G U R E 8 Talus bone of three subjects in this study: middle, shape whose value for the third MoV coefficient is the closest to the mean talus shape, right: a control shape, left: a impingement patient F I G U R E 9 left: a single talus from the impingement group, right: the shape generated by fitting the control group SSM to the abovementioned impinged shape. SSM, statistical shape modeling has been one aim of their study to capture the effect of size, yet as in our study we were not interested in the effect of size, we had a scaling step as mentioned in the alignment phase taking size our of our model. Related to the fact that this is a cross-sectional study meaning that only a single CT scan was studied for each ankle, makes it impossible to study the injury mechanisms and causality of ankle impingement, as well as the progression or changes of disease through time.
Currently, the diagnostic workup of ankle impingement is based on visual analysis of 2D X-rays, which were shown to lack the reliability and information content to diagnose and characterize the syndrome with high sensitivity and in 3D. In this study, we provided a tool that considers the 3D shapes of the bone in a quantitative F I G U R E 10 left: shown as surface in white is bone shape generated by fitting SSM of control group to the bone of one impingement patient. The arrows of deformation field show the direction and amount of distance from each point of the fit to its corresponding point on the impingement shape. Right: the color plot on the surface of the impingement shape showing the same concept of distance between control group SSM fit and the impingement talus. SSM, statistical shape modeling ARBABI ET AL.
| 193 manner, instead of only a 2D projection. In future, SSM can be used to improve diagnoses of impingement and characterization of the bone anatomy in 3D, which may aid prediction and prevention measures as well as adding to surgical planning in clinical practice.

| CONCLUSION
In conclusion, the current study demonstrated the feasibility of 3D SSMs to quantitatively characterize the specific 3D shape variations of talus bones, capturing of modes of shape variations that correlate with ankle impingement as objectified in surgical reports. We also demonstrated the use of SSMs for surgical planning of ankle impingement. We found that two shape modes of the talus significantly differed between the impingement and control groups.