Hierarchical Rigid Registration of Femur Surface Model Based on Anatomical
Features

Existing model registration of individual bones does not have a high certainly of success due to the lack of anatomic semantic. In light of the surface anatomy and functional structure of bones, we hypothesized individual femur models would be aligned through feature points both in geometrical level and in anatomic level, and proposed a hierarchical approach for the rigid registration (HRR) of point cloud models of femur with high resolution. Firstly, a coarse registration between two simplified point cloud models was implemented based on the extraction of geometric feature points (GFPs); and then, according to the anatomic feature points (AFPs) in two level namely shape features and structure features, the fine weight-based registration was performed to achieve anatomical alignment; finally, the origin source model was automatically transformed by applying the obtained coarse matrix and fine one in sequence. Experimental results show that the hierarchical registration method can rapidly and accurately register point clouds of individual femurs, and achieves the medical semantic alignment, and provides a basic tool for the understanding and comparison of femur anatomy and structure.


Introduction
Three-dimensional (3D) models which describe the intuitive and precise anatomic morphology of human bone play a more and more important role in many biomedical fields, including orthopedic diagnosis, surgical planning and basic researches [1]. Based on the volume data acquired from computed tomography (CT) scan slices, point cloud models representing outer surface shape of skeleton can be generated with a set of 3D points. As is well known to all, the shape registration of one sample bone to another one is an essential prerequisite process for some comparative and quantitative applications of femoral surface morphologies, such as the construction of statistical shape model (SSM) [2] and parameter measurement [3]. Unfortunately, the point cloud or mesh models of femur samples can seldom be aligned well in a uniform coordinate system because of various viewpoints of scanner or different modeling platforms [4][5]. Therefore, transforming 3D models and determining the correspondence between them have been indispensable components for surface shape comparison and anatomy statistical analysis.
In the past two decades, many works so far have focused on point cloud registration (PCR) in the research field of computer graphics [4]. The key goal of PCR is to obtain the 3D correspondence of various models, so that the same local parts can be matched each other in a common spatial system. Existing PCR methods can be roughly divided into two major categories in terms of transform properties, namely rigid registration and non-rigid registration [5]. In the former, models were usually transformed with rigid operations including translation and rotation, and the distances and interrelationships of adjacent points could be reserved. While in the latter, 3D models were processed with affine transformations including scaling and shear mapping, so that the input model were eventually deformed to match the reference one. In practice, 3D rigid transformation of skeleton model was more widely used than the non-rigid one in many applications, because that human bones are incredibly strong and their shapes cannot be deformed easily. Accordingly, the rigid registration of human bone model is critical for many orthopedic researches and applications.
There are many excellent reviews in the literature dealing with various concepts of PCR [5][6]. The distance-based, filter-based, probability-based methods were widely applied in many studies [7][8][9][10][11][12][13][14][15]. Distance-based registrations are the most common method, which were achieved with two processes including correspondence finding and distance error minimizing. Besl et al. [6] proposed a classical algorithm named iterative closest point (ICP), in their study, the mean squared error of Euclidean distance of point pairs was iteratively minimized. Subsequently, to accelerate convergent speed and improve registration accuracy, many efficient ICP variants including FICP, EM-ICP [7] and LM-ICP [8], were subsequently proposed by using different strategies such as point selection, point matching, pair weighting, pair rejecting and error metric minimization [9]. ICP-based algorithms have the advantage of low computational complexity and support parallel computing. Nevertheless, these approaches are sensitive to the initial pose and are easy to be trapped into local minima for complex models considerably. To overcome major drawbacks of ICP-based methods, the robust point matching (RPM) method by employing deterministic annealing and soft-assign optimization [10][11] was proposed to obtain the correspondence and transformation parameter of cloud points. However, the applications of RPM method were limited in the complex environment with noise and outliers [12]. Later, some filterbased methods were suggested and achieved by using state space model. Ma et al. [13] proposed an unscented particle filter for rigid registration. Sandhu and his team [14] suggested a particle filter schema in which the local optimizer formulated based on point matching. Although filter-based methods were used to process a large number of point sets, their primary limitation is the pre-calculating of the corresponding relations of point sets. Myronenko et al. [15] proposed a novel approach named coherent point drift (CPD), which core idea was the solving of a maximum likelihood (ML) estimation problem by using Gaussian mixture model (GMM). Probability-based method has better performance, but the calculation cost was significantly higher than other methods.
Most PCR studies have achieved the automatic alignment of models constructed in different views or platforms. It is important to highlight that the source model and the reference are various poses of the same object, and they can overlap completely in theory. It is well known that the structure of one human femur is similar to others, and there are significant differences between individuals. In fact, it should be noted that very littler work has been carried out on the alignment of different individual bones. On one hand, the large number of vertices in high-precision point cloud model increases the computational complexity of geometric objects, and the improvement of registration algorithm efficiency is still a challenging work. On the other hand, existing PCR researches pay little attention to the anatomy alignment, and little evaluation criteria for registration result between different skeletons has been suggested. Therefore, how to align different bone models with medical anatomic semantics for biomedical applications is yet another hard issue to be researched. In summary, model registration of different bones does not have a high certainly of success because of the lack of anatomic semantic. In light of the surface anatomy and functionally structure of bone, we hypothesized different point clouds of femur could be aligned through feature points both in geometrical level and in anatomic level. In this study, a hierarchical approach for the rigid registration (HRR) of individual femur models with high resolution and was propose based on anatomical feature to achieve the rapid alignment of individual bones and to facilitate the understand and comparison of anatomy and structure of femur.
The paper is organized as follows. Section 2 presents the study materials. Section 3 describes an overview of the proposed method of model registration. Section 4 provides detailed algorithms for point cloud preprocessing, feature extraction and registrations in two levels. In Section 5, the approach was implemented, and experimental results were tested, then time-consuming, error and anatomic alignment were discussed. Section 6 concludes with a summary of the content of this study and proposes future research directions.

Materials
In the study, 20 heathy right femurs of Chinese female volunteers lived in Southern Jiangsu of China with no previous trauma and no history of diseases were selected and scanned by 64-slice CT (MSCT, Aquilion 64, Toshiba, Zoetermeer, Netherlands), and then 3D models were constructed in the form of point cloud. The average age and the average height of volunteers were 37.4 ± 6.2 years (range: 30 to 39) and 158.3 ± 5.3 mm (range: 156.7 to 163.2). Although the selected femur samples are very similar in overall structure, there are significant differences in length and detailed shapes between individuals.

Overview of Proposed Approach
Given preprocessed point cloud models of femur, our approach aims at registering the two models with rigid transformation and achieving the anatomic semantic alignment. Firstly, a coarse registration between two simplified point cloud models was implemented based on the extraction of geometric feature points (GFPs); and then, according to the anatomic feature points (AFPs) in two level namely shape features and structure features, the fine registration was performed by using the weight-based method to achieve anatomical alignment; finally, the origin source model was automatically transformed by applying the coarse matrix and the fine one in sequence.
The complete flow of hierarchical registration of femur model based on anatomic features is shown in Fig. 1. Steps II-IV are the most important of the steps and will be expounded in the following sections.

Point Cloud Preprocessing
In the original 3D point cloud model of individual samples reconstructed from CT scanning images, available points are usually accompanied by some noise and outliers. To remove the noise points, individual original models were preprocessed by using filtering and smoothing before subsequent 3D modelling. Recently, many effective methods have been proposed to remove the outliers such as based on partial differential equations [16], based on Laplacian signal processing [17], based on neighborhood filtering [18], based on projection [19] and based on statistical [20]. Considering the preservation of the original shape of bone, especially local details in the noise removal stage, the denoising smoothing algorithm was applied in this work based on the method proposed by Gu [21]. The principal pipeline of the algorithm and a preprocessed smooth model (abbreviated as M P ) were shown in Figs. 2 and 3 respectively.

Geometric Feature Extraction
High-precision model composed of enormous geometric points have the advantages in the detail representation of complex shapes, however, the process of mass geometric objects in low-level may consume large amounts of computational and memory resources. Considering different contributions of each point on the description of skeletal shape, reducing point number in the model under the premise of accurate representation may be an ideal approach to improve the processing efficiency of point cloud model. Therefore, the detection and extraction geometric feature points was an indispensable basic technology in this work.
In the geometric level, normal vector and curvature were typical factors and widely used in the shape feature description. The point which normal angle was larger than a threshold value can be regarded as geometric feature point [22]. In traditional curvature-based methods [23][24], the feature points were extracted according to the bending or smoothness of local surfaces reflected by curvatures of each point. To avoid the limitation, normal vector sensitive for noise and outliers, we proposed an optimized curvature-based strategy for geometric feature extraction in this study, and the feature points were selected from the curvature variation of neighbors. The detailed processes of curvature solution were processed as follows.
Step 1. Calculate the gravity center O cnt of point P i : where P nbi was the neighborhoods of P i , k was the number of the neighbors.  Step 2. Obtain the normal vector N i with the minimization of fitted function F in Eq. (2), where a, b, c, d and e were the coefficients.
Step 3. Calculate the curvatures of points: where K Gaus was Gaussian curvature, K Mean was mean curvature, K 1 and K 2 were principal curvatures, respectively.
Step 4. Select GFPs based on curvature variation: where C idx (P i ) was the curvature factor of P i which ranged from 0 to 1: P ccv ¼ fP i ;C idx ðP i Þ <¼ v 2 g (11) According to Eqs. (10) and (11), both the points which C idx (P i ) were larger than the threshold ω 1 in convexity and the ones which C idx (P i ) were less than the threshold ω 2 in concavity were extracted and defined as the point set P cvx and P cvv respectively, then P cvv and P cvx were combined into the geometric feature model (abbreviated as M G ) as show in Eq. (12).
Therefore, the M P shown in Fig. 3 was replaced by M G with a smaller number of GFPs as shown in Fig. 4.

Coarse Registration of Geometric Feature
At this stage, the sample model which length was the closest to the mean value of all M G models was selected as the reference model named with M GR , and the other GFP models were defined as source models namely M GS , and all M GS were registered to M GR coarsely in geometric level.
The overview registration process was depicted in Fig. 5. To meet the requirement of initial posts and avoid falling into local optimization, the centroid alignment between two models was achieved by translating the M GS with a matrix T C . As a result, M GS and M GR were locally overlapped. Then, the M GS was registered to the M GR with the ICP algorithm based on KD tree researching, and the rotation matrix R G and translation matrix T G were obtained. The key steps were explained as follows: Step 1. Create the centroids model, named as P c shown in Fig. 6(a) by computing the average coordinate values of each point cloud according to in Eq. (13).
where n was the number of the GFPs in a cloud model M G .
Step 2. Translate the M GS to realize the centroid alignment between M GS and M GR as shown in Fig. 6(b).
Step 3. Constrain the centroid of M GR , and the M GS was iteratively transformed until either the sum of distance between two models was smaller than a threshold value or the iteration times was over a specific number to obtain the coast registered model (M GC ) as shown in Fig. 7.

Anatomic Feature Extraction
Unlike other 3D models, there are always special structures and anatomic shapes in femur surface model, such as the head, neck, condyles and trochanters, and these anatomic features have distinct function roles due to being a part of hip joint or knee joint. It is of great significance for many clinical applications to combine the anatomical features of bones in the alignment of bone. Therefore, using anatomic feature points (AFPs) to describe the functional structure and local shapes in point cloud model was a premise task of 3D alignment. In this work, AFPs were defined and extracted from two aspects. Feature points of shape (P s ) were extracted from the surface based on the prior medical knowledge to represent detailed anatomy in local regions. In addition, Feature points of function (P f ) were created by discrete the fitted axes to describe the anatomic structure of a whole bone.  The definition scheme of anatomy model (M A ) with of AFPs was shown in Fig. 8, and M A consists of the set P s and P f . The detailed methods for feature definition and extraction can be found in our previous work [25], and the formal definition was given below: s:t: P s ¼ fP s1 ; P s2 ; P s3 ; . . . ; P sm g s:t: P f ¼ fP f 1 ; P f 2 ; P f 3 ; . . . ; P fn g

Fine Registration of Anatomic Feature
The aims of fine registration of anatomic feature were to accomplish the anatomic alignment and improve the medical significance based on the extracted AFPs, so that not only the detailed regions with significant medical semantics but also the functional structure of the whole skeleton may be aligned as possible.
The coast registered model (M GC ) obtained in Subsection 4.3 was transformed again by employing a weighted rigid registration. To achieve the goal of anatomy alignment, and the point clouds were processed with the following steps.
Step 1. Group the M GC and M A into a new point cloud named with feature model M F .
Step 2. Transform the new source model M FS to match the reference M FR by using ICP rigid registration.
Step 3. Calculate the rotation matrix R A and translation matrix T A . Therefore, the fine registration problem was converted to solve the minimization problem of the sum of Euclidean distance between point pairs of M F s as shown in Eq. (16).
where, E F (R A ,T A ) was the object function, M GCR and M GCS were the reference model and source model of M GC , while M AR and M AS were the reference and source model of M A , k and s were the point numbers in M GC and M A . β was the influence weight, which value were defined in the range of 1 to 100, and the default value was set to 20.

Final Transformation of Origin Model
Based on the mentioned transformation in two levels, namely coarse registration of geometric feature and fine registration of anatomic feature, the simplified model M F was aligned to the reference one. However, the description accuracy of the simplified model (M F ) was obviously less than the one of original preprocessed model (M P ) which composed of a large number of point sets. In consequence, the original preprocessed model (M P ) should be transformed according to the hierarchical matrixes obtained in the two registration processes to match the reference model. The final registered source model was (M R ) acquired with the following steps.
Step 1. Translate the source model M FC according to the vector T C obtained in the coarse registration, and create the intermediate model M FC1 .
Step 2. Transform M FC1 with rotation and translation in geometric level based on the matrix R G and T G respectively, and create the intermediate model M FC2 .
Step 3. Readjustment M FC2 in anatomic level based on the matrix R A and T A respectively, and obtain the final model M R .

Results and Discussion
To evaluate and verify the feasibility and effectiveness of the methodology and algorithms proposed in this study, 20 point-clouds of right femur were tested, and the error analysis, semantic alignment and timeconsuming were carried out to illustrate the registration efficiency and the accuracy in anatomic semantic.

Experiment Results
Using a 2.2 GHz i5-5200u CPU with 4G RAM, the proposed methodology and algorithms were implemented in the software environment of MATLAB 7.11.0.

Registration Results
To compare the align results conveniently, the same reference model was used in each registration experiment, while the source femurs were in various poses as shown in Figs. 9(a) to 9(c). The registration results of three samples selected were shown in Figs. 9(d) to 9(f).
To improve the medical semantics of registration, the value of the influence weight β was set to 10, 50 and 100, respectively, and the results of optimized experiments were tested and shown in Fig. 10.

Method Comparisons
Four representative point set registration algorithms including ICP, FICP, CPD, KC were selected to conduct the similar experiments. The performance of these common algorithms was validated on the same point clouds. First, the Root Mean Square (RMS) error was used to evaluate the performance. The registration results and the total processing times were calculated, compared and listed in Fig. 11. Then, the shaft center distance, namely the Euclidean distance between the shaft-centers of two registered samples proposed in our previous works [26][27] was used for the semantics effect assessment of registration. The average values of shaft center distance of five methodologies were compared and listed in Tab. 1.

Time-Consuming
Under the premise of preserving the original shape and anatomic structure information of femur model, using a small number of feature points to simplify the skeleton model and reducing the scale of the original  Fig. 11, it was found that not only geometric feature extraction but also anatomic feature extraction are important stages to improve the efficiency of the algorithm.

Error Analysis
As shown in the right column of Fig. 11, the final RMS errors of former four algorithms were found to be similar and were limited within 10 mm under the same experiment conditions. Nevertheless, the registered models (the red) in former two rows were inverted, and neither the result of ICP nor the one of FICP was reasonable. The main reason for the inaccuracy was that the ICP-based processes fell into local minimization due to the lack of overlap between the two models [4]. Compared with results in the third and fourth rows, it was shown that the final result our HRR can have a better match, especially in the parts of head and the condyle, and the alignments of these special functional parts are very important for medical applications. Therefore, on the basis of comparison results, it can be concluded that the HRR has certain advantages in anatomical alignment.

Anatomic Alignment
Considerable research efforts have been devoted to 3D point cloud registration, and much work so far has focused on the error minimization of geometric objects [5]. Anatomic features of the femur, such as the head, trochanters and condyles have distinct functional roles in the joint connections and daily activities [28]. Therefore, the bone models resulted from existing typical methods were seldomly accepted by medical professionals because of the lack of attention to anatomic semantics. In brief, not only local shapes and typical structure are of great significance for many medical researches and related applications. As far as we know, there is no report on anatomic registration with medical semantics.
In this study, anatomic features (AFPs) were defined and extracted in two levels (P s and P f ). On one hand, P s s were extracted from the point cloud with the combination of prior medical knowledges to  describe surface features. On the other hand, P f s were generated with the discretization of the fitted anatomic axes to represent structure features. Based on the registration of geometric features, the feature registration of anatomic and structure features were implemented by using a linear weighing method. As a result, the RMS error value of HRR was less than anyone of the other methods as shown in Fig. 11. What's more, the comparison results in Fig. 10 show that the influence of the AFPs in the entire feature points (the pint set of GFPs and AFPs) could be adjusted by the weight factor β, thereby, anatomical semantics of the model can be improved according to the ratio of numbers between GFPs and AFPs.
So far, there are no clear criteria for assessing anatomical alignment between femur models. Considering the significant individual differences in skeletal length and local shape, a simplified scheme to evaluate the alignment error was proposed in this paper. The center-offset, namely the distance between the central point of the shaft isthmus of registered source model and the one of the reference model. In this study, the overlap between centroid M G of source and the one of reference model was defined as a constraint in all registrations, therefore, the center-offset of HRR is less than the ones of other methods as shown in Tab. 1. Although the centroid M G may not match the center of shaft isthmus, the constraint of centroids overlap can facilitate the anatomical alignment.
The evaluation of anatomic and semantic alignments between different bones cannot only depend on the geometric objects, and the prior knowledge, medical anatomy and bone characteristics should be taken into consideration. Therefore, the study of describing the differences between bones with typical effective method, such as probability statistics and feature parameterization, which will be one of the research directions in the future.

Conclusion
To facilitate the anatomic alignment of different femur point clouds, this study proposed a novel mechanism for hierarchical registration of different models based on anatomic features of femur. The femur model consisting of a large number of point clouds was simplified with a small amount of feature points in geometric level to deduce the compute cost of model registration, and can improve the reconstruction efficiency. In addition, the features of anatomic surface and functional structure were defined and extracted according to prior medical knowledge, and were taken part in the rigid registration to achieve the anatomic alignment. This research indicates that for the first time we have proposed a novel method for registering individual femur models in anatomic level rather than only in geometric level. The main characteristics of this approach are as follows.
It simplifies the original femur model consists of large number of point clouds with the feature points in two-level, namely geometric feature points (GFPs) and anatomical feature points (AFPs), and improves the efficiency of registration algorithm significantly. It defines hierarchical anatomic features from two aspects of shape and function, and achieves the hi representation of anatomic features and the alignment of medical semantic It extends the existing methodology for 3D rigid registration, and provides the alignment of medical semantic for further medical applications However, some directions must be researched in future work, including the following avenues of investigation: