- Email: [email protected]
Predictive geometrical model of the upper extremity of human fibula
biocybernetics and biomedical engineering 36 (2016) 172–181
Available online at www.sciencedirect.com
ScienceDirect journal homepage: www.elsevier.com/locate/bbe
Original Research Article
Predictive geometrical model of the upper extremity of human fibula Milica Tufegdzic a,*, Stojanka Arsic b, Miroslav Trajanovic c a
Department for Mechanical Engineering, Machine-electrotechnical School, Krusevac, Serbia Department of Anatomy, Faculty of Medicine Nis, University of Nis, Nis, Serbia c Department for Production, IT and Management, Faculty of Mechanical Engineering, University of Nis, Nis, Serbia b
article info
abstract
Article history:
Computer assisted preoperative planning in orthopedic surgery, as well as designing and
Received 27 June 2015
manufacturing of personalized fixators, implants and scaffolds requires a good three-
Received in revised form
dimensional model of bone(s) of the treated patients. Existing methods that convert the
30 November 2015
Computer Tomography (CT) images into the polygonal three-dimensional models are time-
Accepted 11 December 2015
consuming and inefficient. Therefore, we propose a predictive model that allows quick
Available online 29 December 2015
creation of three-dimensional (3D) surface model of a particular bone by measuring the relevant parameters from an X-ray or CT image.
Keywords: Geometrical model
In this paper, we present the process of creating a predictive geometrical model using the case of proximal end of fibula as an example. The predictive model is built by defining the
Parameters
referential geometric entities that correspond to anatomical features, based on which
Fibula
appropriate points, axes, planes and curves are created. Using the method of linear and
Prediction
nonlinear regression with four different parameters, which can be measured from X-ray
Regression model
images or anterior-posterior projection of fibula at CT scans, the equations for X, Y and Z coordinates of the selected 168 points are obtained and their predictive values are calculated. These values are used for creating 3D surface model with the aim of two different methods: using loft function and converting these coordinates into point cloud. These models were compared and verified through analysis of deviations and distances between initial model and predictive models. The resulting 3D model has satisfactory accuracy, and the process of its building is much shorter. # 2015 Nałęcz Institute of Biocybernetics and Biomedical Engineering. Published by Elsevier Sp. z o.o. All rights reserved.
* Corresponding author at: Department for Mechanical Engineering, Machine-electrotechnical School, Cirila i Metodija 26, 37000 Krusevac, Serbia. E-mail addresses: [email protected] (M. Tufegdzic), [email protected] (S. Arsic), [email protected] (M. Trajanovic). http://dx.doi.org/10.1016/j.bbe.2015.12.003 0208-5216/# 2015 Nałęcz Institute of Biocybernetics and Biomedical Engineering. Published by Elsevier Sp. z o.o. All rights reserved.
biocybernetics and biomedical engineering 36 (2016) 172–181
1.
Introduction
Three-dimensional (3D) models with accurate geometric representation of the long bones are being increasingly used for various aspects of clinical practice and research. They are an indispensable basis for preoperative planning in orthopedic surgery, designing and manufacturing of personalized fixators, implants and scaffolds. There are many programs that allow obtaining the polygonal 3D models of bones based on Computer Tomography (CT) images in DICOM format. Independently from the program, building a 3D surface model of a bone from CT scans must be conducted through several steps: data management, image enhancement, segmentation and 3D reconstruction [1]. Defining the boundaries of the bone is necessary for healing of the model, but it often cannot be done accurately, the process is time consuming and requires a great deal of patience. It is well known that in clinical practice the patient cannot undergo a CT scan because of previous exposure to very high levels of ionizing radiation. In such cases it is impossible to build a 3D model of a bone. Our idea is to create a method for obtaining a predictive model of long bones from X-ray or CT images, using Computer Aided Design (CAD) program and statistical tools. In this paper, we are presenting the method using the upper extremity of human fibula as an example. According to the anatomical description, the human fibula is a long, slender bone, located at the lateral aspect of the leg. It has three main parts: a proximal head, a narrow neck, a long shaft and a distal lateral malleolus [2]. The upper extremity of the fibula is composed of two parts, head and neck. The head of the fibula has a most prominent part called the apex and a flattened articular surface, which forms the proximal tibiofibular joint with a corresponding fibular articular surface on the lateral condyle of the tibia. Articular capsule is reinforced by the anterior and posterior tibio-fibular ligament of the fibular head. Main parts of human fibula are presented in Fig. 1 (down below). It is important to point out, as an additional biomechanical and functional description of the fibular upper end, that there are two types of superior tibio-fibular joints, horizontal and oblique, due to the anatomical variations [3]. They provide a compensatory motion in internal and external rotational movements of the tibia. The primary functions of the superior tibio-fibular joint are: (1) the dissipation of torsional stresses which occur at the ankle joint, (2) wasting of the lateral tibial flexing, and (3) tensile weight bearing [3]. The previous work on the 3D modeling of the fibula or the parts of the fibula has resulted in generating a 3D model of fibula bone ankle joint (including the lower end of the fibula) by using software for processing the medical images and the CAD software. The process of the 3D modeling of human fibula using medical image software is described in [4]. IGES surface of the human fibula using MIMICS software and a 3D surface model in Solid Works, in order to test a newly developed and adapted osteosynthesis implant, used for transsindesmotic fibula's fracture is obtained in [5]. In another paper a finite element model of the knee joint, including distal femur, tibia, fibula, menisci, articular cartilages and ligaments is obtained from MR images of a healthy male [6]. The models of different
173
bones such as distal femur, patella and upper extremities of fibula and tibia are obtained from the CT data. These CT slices were converted into the point clouds (using Mimics) and a CAD surface model was obtained in Geomagic. These CAD models were further processed for the purpose of Finite Element Analysis [7]. In order to study the influence of the fibula and talus on the distribution of stress on the tibia, a 3D solid model of the fibula was created using CT images for generating 3D solid models of tibia, fibula and talus, with the help of Mimics in the form of point clouds [8]. These point clouds were exported in Geomagic where the failure fixing and smoothing on the surfaces was performed. The 3D models were saved in STL format and later assembled using Mimics FEA module. A case of reconstruction of a maxillary defect and orbital floor with a micro vascular fibula graft is shown in [9]. An individualized titanium mesh, as a rapidly prototyped template for the missing part of fibula, is designed from a CT scan using Materialize. 3D surface model of the fibula was obtained in the process of Reverse Engineering (RE) [10], but this process was timeconsuming, as all processes mentioned above.
2.
Materials and methods
In the process of generating a predictive model of the proximal end of human fibula, a method of anatomical features (MAF) was used. This method was previously described in [11] and has been implemented on femur and tibia bones. In the case of fibula bone, first few steps of this method were applied in our research. But, because of the fact that linear regression was not suitable for fibula bone, the method was expanded by introducing different nonlinear regression models. These models were tested by statistical tools and for the level of significance a value less then 0.05 was adopted. In order to determine the most suitable model, from those which were statically significance, we have chosen one bone from the sample for which the value of the length of mechanical axis is the closest to the mean value of the sample. Regression values of coordinates for corresponding points on this bone were calculated for each model and compared with initial values. So, for the models were adopted the ones for which the difference from regression values was the least. Finally, we have created predictive geometrical models from the regression coordinates, with the aim to compare and verify our results using the deviation and distance analysis. So, our research was conducted through the following steps: CT scanning and digitalization of data, definition and determination of anatomical landmarks on the upper and lower end of fibula, creation of the referential geometric entities (RGEs) (points, axes, planes, surfaces) according to the anatomical and morphological characteristics of the bone [12,13], selection and measurement of the parameters on the upper end of the fibula, identifying and measuring the coordinates of the points for which the regression models are formed, mathematical models determination of the regression equations,
174
biocybernetics and biomedical engineering 36 (2016) 172–181
calculating the regression values for X, Y and Z coordinates, calculating the difference between regression and initial values for ‘‘control’’ bone and adoption of the best regression model, creation of the models using loft function and converting point coordinates into point cloud, and verification and comparison between the initial and obtained models, by using the deviation and distance analysis.
2.1.
CT scanning and digitalization of data
The two available scanning methods for quantitative 3D imaging are CT and MRI. Both CT and MRI provide accurate information for quantifying of the anatomical structures in a 3D environment. Both imaging methods have considerable application in clinical praxis, but CT is mainly used for the bone imaging and MRI for the soft tissue imaging. To create a predictive model of the upper extremity of the human fibula, a series of traverse CT images was taken from 11 patients (6 men and 5 women), 43–70 years old (58 years on average). These images were obtained on a Toshiba MSCT scanner Aquillion 64 (120 kV, 150 mA, thickness 1 mm, inplane resolution 0.781 0.781 (pixel size), acquisition matrix 512 512; field of view (FOV) 400 mm 400 mm). All the CT data were saved in DICOM format. After digital processing, these 2D data are converted to 3D STL (STereoLitography) models [14], considering that 3D modeling is more precise than 2D methods because it is less influenced by body position [15]. These models were imported in the CAD software and then further processed by eliminating the model errors (isolated triangles, non-manifold vertices, edges, etc.) and removing all unnecessary entities, healing, mesh creation, optimization and mesh smoothing. Model obtained in this manner was used for determining the RGEs in a CAD program.
2.2.
A-P (anterior–posterior) plane – defined by mechanical axis and ACF. ACF, AML, mechanical axis and the main parts of human fibula are presented in Fig. 1. The next stage was determining and measuring the 4 parameters in A-P plane (d1, d2, d3 and d4). These parameters are shown in Fig. 2: d1 – distance from the highest point of FACF in the A-P plane to the most distal point on the upper end of the fibula, measured in the direction perpendicular to the mechanical axis,
Methods
For the purpose of our research, RGEs at all fibula bones were identified. It was done according to the definitions described in [10,16]: (1) at the upper extremity (lat. epiphysis proximalis s. extremitas proximalis): ACF (lat. apex capitis fibulae) – apex of the fibular head, FACF (lat. facies articularis capitis fibulae) – surface for articulation with the lateral condyle of tibia, (2) at the lower extremity (epiphysis distalis s. extremitas distalis): FAML (lat. facies articularis malleoli lateralis) – surface for articulation with the talus faced medially, AML (lat. apex malleoli lateralis) – top of the lateral malleolus, as the most distal point on the lateral malleolus. At the basis of the previously defined RGEs the following entities were defined: the mechanical axis – which connects the centers of the articular surfaces at the upper (FACF) and the lower extremity (FAML) of fibula, and
Fig. 1 – Main parts and RGEs at the right human fibula, medial aspect.
biocybernetics and biomedical engineering 36 (2016) 172–181
Fig. 2 – Parameters on the upper extremity of the right human fibula.
d2 – distance from the center of FACF and the most distal point on the upper end of the fibula, measured in the direction perpendicular to the mechanical axis, d3 – distance from the lowest point on FACF in the A-P plane to the most distal point on the upper end of the fibula, measured in the direction perpendicular to the mechanical axis, d4 – the distance at the level (the intersection) where the head goes into the neck of the fibula, measured in the direction perpendicular to the mechanical axis. With the aim to get coordinates of points and their values, we have used two sets of curves, obtained with two sets of planes. These planes were utilized for intersections of polygonal model of the upper end of fibula. Measured values were applied for obtaining predictive values of X, Y and Z coordinates. The first set of planes is perpendicular to mechanical axis and set at places of the selected parameters (d1, d2, d3 and d4). So, the polygonal models of the upper end of fibula are intersected with these planes. As result, we got a corresponding set of curved lines (Curve 1, 2, 3 and 4, Fig. 3). The second set of planes was constructed at certain angles in relation to the A-P plane (158, 308, 458, 608, 758, 908, . . .), where a line parallel to the mechanical axis which passes through the ACF was used as the axis of rotation, according to our previous research [10]. The sections of the polygonal model of the upper end of fibula with these planes (including A-P plane), have resulted in corresponding set of 12 rotational curves, with ACF as common point. These two sets of curves always have intersections in two points and we have got 96 points. Because of the fact that the distance between the parameters d3 and d4 was great and we could not get satisfactory accuracy of the model, we have added 6 points on every curve from the second
175
Fig. 3 – Curves and points of intersections at the upper end of fibula.
set of planes. These points were chosen at the same distance (5, 10 and 15 mm) from the points of intersections of the curve 3 (Fig. 3) and the curves from the second sets of planes. At Fig. 3 are presented points of intersections between curves obtained at the section of polygonal model of the upper end with the A-P plane and the planes which are rotated for the angles of 158, 308, 458, 608, as an examples (respectively identified as A-P curve, 158 curve, 308 curve, 458 curve, 608 curve), and Curves 1, 2, 3 and 4. Points of intersection of rotational curves are marked as A-P-1, A-P-2 on Curves 1, 2, 3 and 4. The same methodology is applied for labeling corresponding points at the other rotational curves (15-1, 15-2, etc.), while additional points between Curves 3 and 4 are labeled as A-P-1-1, A-P-1-2, A-P-1-3, A-P-2-1, A-P-2-2, A-P-2-3, 15-1-1, and so on (for the simplicity of Fig. 3 only these points are marked). The coordinates of points of intersections in the newly formed coordinate system of the right orientation were measured. Because of the fact that ACF is the common point for the second set of curves, it was adopted for the origin, the X-axis is perpendicular to the A-P plane, Y-axis ‘‘lies’’ in the A-P plane (perpendicular to the mechanical axis), while the Z-axis is parallel to the mechanical axis. So, the sets of X, Y and Z coordinates were obtained and their values were measured. These sets of coordinates represent specific locations on the upper end of the fibula. The values of these coordinates for 10 fibula bones are the input data for the formation of the predictive model, while the 11th bone is separated and used as a ‘‘control’’ fibula for determination of the most suitable regression model. The main criterion for the selection of ‘‘control’’ fibula was mechanical axis length. It was selected in such manner that its length has the smallest deviation from the mean value of mechanical axis length of the sample.
176
biocybernetics and biomedical engineering 36 (2016) 172–181
Values of measured coordinates (Xi, Yi, Zi) were the input data for the statistical analysis, treated as dependent variables, while the values of the parameters d1, d2, d3 and d4 were treated as independent variables. The relationships between these two sets of variables were established by means of correlation coefficients and scatter plots were used for choosing appropriate regression models. Based on this analysis, 6 regression models (1 linear and 5 non-linear) were chosen and tested for every predicted value of the coordinates: (1) Linear model represented by the expression Xi, Yi, Zi = a + b di, (2) Non-linear models: (a) quadratic models, represented by the expressions Xi, Yi, Zi = a + b di2 and Xi, Yi, Zi = a + b di + c di2; (b) 2 logarithmic models represented by the expressions Xi, Yi, Zi = a + b ln(di) and Xi, Yi, Zi = a ln(di); (c) exponential model represented by the expressions Xi, Yi, Zi = a exp(b di); i = 1, 2, 3, 4; a, b – coefficients of linear regression, while di stands for the values of the parameters d1–d4. Coefficients in regression equations were calculated using statistical software and taking into account that the value of statistical significance p is p < 0.05. Models that have p-values greater than 0.05 for any of the coefficients in the regression equations were discarded. From the models in which the p-value was less than 0.05 was elected the one with a minimal deviation of regression (predictive) values of the X, Y and Z coordinates and the measured coordinate values for corresponding points at the ‘‘control’’ fibula. For additional points between Curves 3 and 4 the same statistical procedure was used, but some different models were tested because the values of the coordinates X, Y and Z depends from both parameters d3 and d4, according to correlation matrix. These models are: (1) linear, represented by the expressions (a) Xi, Yi, Zi = a + b di, i = 3, 4; (b) Xi, Yi, Zi = a + b d3 + c d4; (2) non-linear, represented by the expressions (a) Xi, Yi, Zi = a + b di2, i = 3, 4Xi, Yi, Zi = a + b d3 + c d42, Xi, Yi, Zi = a + b d4 + c d32; (b) Xi, Yi, Zi = a + b ln(di), Xi, Yi, Zi = a ln(di), i = 3, 4; (c) Xi, Yi, Zi = a exp(b di), i = 3, 4. The method of election of the appropriate model was the same as for the points of intersection, described above. Calculated predictive values of the X, Y and Z coordinates are used for constructing sets of points. Splines were constructed from matching points and presented in Fig. 4. By the loft function in CATIA Shape module surface predictive model was created. These sets of points were also transformed into the point cloud, from which the mesh could be created. These two models are further used for comparison and verification through the analysis of deviations and distances between initial polygonal model and obtained models.
Fig. 4 – Splines through regressive points.
3.
Results
Regression models of the points coordinates (Fig. 3, points A-P-1 and A-P-2 at Curves 1, 2, 3 and 4) are presented by equations in Tables 1–4. Different regression models are obtained for point coordinates X and Y at the Curves 1–4, while the regression equation for Z coordinate remains the same for all points at selected curve, because of the fact that proximal fibular end is cut by the planes perpendicular to the mechanical axis and Z axis which is parallel to the mechanical axis. Regression models of the points coordinates X, Y and Z between Curves 3 and 4 are presented in Tables 5 and 6 (as examples). For the Y coordinate for the points at the A-P curve the values for coefficients of regression are not calculated, because all the values are 0. Another statistical analysis of the deviations between the predictive coordinate values and actual (real) coordinate values was conducted. The results obtained by means of descriptive statistics are given in Tables 7–10 for the X and Y coordinates (for Curves 1–4), where Valid N represents the number of predictive (calculated) values. Analysis for Z coordinate was not conducted because there is a unique value at each curve. Based on the values of means and standard errors presented in Tables 7–10 we have concluded that all deviations between predicted and real values of coordinates can be treated as very good. That means our predicted values are very close to the real values. Similar analysis of the deviations was conducted for deviations of all X, Y and Z coordinates for additional points. The results are presented in Table 11, where Valid N represents the number of predictive (calculated) values. Some exceptions should be noticed in this case, especially for Y coordinate. The deviations could be reduced and regression models improved by increasing the number of samples. Another statistical analysis of the deviations between the predictive coordinate values and actual (real) coordinate
177
biocybernetics and biomedical engineering 36 (2016) 172–181
Table 1 – Regression equations for X, Y and Z coordinates of selected points at the parameter d1. Coord.
Points at Curve 1 A-P-1 X = 1.39372 exp(0.13076 d1) Y = 0.722038 + 0.157627 d1 0.008025 d1 2
X Y Z
X = 1.714685 ln(d1) Y = 0.495697 + 0.105480 d1
Discussion
where the first set of planes intersect polygonal model of the upper end of fibula (Curve 1 – parameter d1, etc.). Some other regularity can be seen from the obtained equation regression for: (a) X coordinates for points at Curve 1 – 3 exponential models, Y coordinates – 2 quadratic and 2 exponential models, all others are logarithmic;
Table 5 – Regression equations for X, Y and Z coordinates of points between the parameters d3 and d4 at the A-P curve. Coord.
Different dependences in the regression equations can be treated as a result of the variability of the upper fibular end (regarding to the sex, different types of the FACF and different geometry of the samples). But, some regularity can be noticed in regression models. The values of X, Y and Z coordinates for the points of intersections at Curves 1, 2, 3 and 4 directly depend from parameters which were measured at the places
X Z
Coord.
Table 2 – Regression equations for X, Y and Z coordinates of selected points at the parameter d2. Points at Curve 2
Coord. A-P-1
A-P-2
X = 3.89353 ln(d2) X = 1.410168 ln(d2) Y = 0.000436 ln(d2) Y = 0.004226 ln(d2) Z = 3.072500 ln(d2)
Table 3 – Regression equations for X, Y and Z coordinates of selected points at the parameter d3. Coord.
Points at Curve 3 A-P-2
X Y Z
A-P-2
X = 5.8614177 ln(d3) X = 28.1574 + 9.8685 ln(d3) Y = 0.006120 ln(d3) Y = 0.000389 ln(d3) Z = 4.069104 ln(d3)
Coord.
X Z
Coord.
X Y Z
Coord.
X Y Z
Coord.
Coord.
A-P-1 X Y Z
A-P-2
X = 7.75035 exp(0.03694 d4) X = 0.554902 ln(d4) Y = 0.002218 ln(d4) Y = 0.000945 ln(d4) Z = 10.334974 ln(d4)
A-P-1-1
A-P-2-1
X = 6.04468 ln(d3) Z = 6.134285 ln(d4)
X = 23.2087 + 9.2082 ln(d4) Z = 6.264310 ln(d4)
Points at the distance of 10 mm A-P-2-2
X = 9.77308-0.02950 d4 Z = 7.791366 ln(d4)
2
X = 0.882985 ln(d4) Z = 7.970058 ln(d4)
Points at the distance of 15 mm A-P-1-3
A-P-2-3
X = 6.71042 ln(d4) Z = 8.232533 ln(d3)
X = 34.7029 + 12.4725 ln(d4) Z = 9.438726 ln(d4)
Table 6 – Regression equations for X, Y and Z coordinates of points between the parameters d3 and d4 at the 158 curve.
Table 4 – Regression equations for X, Y and Z coordinates of selected points at the parameter d4. Points at curve 4
Points at the distance of 5 mm
A-P-1-2 X Z
X Y Z
0.005354 d1 2
Z = 1.270169 ln(d1)
values was conducted. The results obtained by means of descriptive statistics are given in Tables 7–10 for the X and Y coordinates (for Curves 1–4), where valid N represents the number of predictive (calculated) values. Analysis for Z coordinate was not conducted because there is a unique value at each curve (Tables 12–15). Based on the values of means and standard errors presented in Tables 7–10 we have concluded that all deviations between predicted and real values of coordinates can be treated as very good. That means our predicted values are very close to the real values.
4.
A-P-2
X Y Z
Points at the distance of 5 mm 15-1-1
15-2-1
X = 5.47192 ln(d3) Y = 1.79856 exp(0.05366 d4) Z = 6.099127 ln(d4)
X = 1.749891 ln(d3) Y = 0.268283 ln(d3) Z = 6.574582 ln(d4)
Points at the distance of 10 mm 15-1-2
15-2-2
X = 5.49791 ln(d3) Y = 1.65101 ln(d4) Z = 7.719516 ln(d4)
X = 1.256039 ln(d3) Y = 0.220073 ln(d3) Z = 7.930495 ln(d4)
Points at the distance of 15 mm 15-1-3
15-2-3
X = 5.28828 ln(d3) Y = 1.41699 ln(d3) Z = 9.246085 ln(d4)
X = 48.8380 + 15.6115 ln(d3) Y = 31.9905 12.5535 ln(d3) Z = 8.965555 ln(d4)
178
biocybernetics and biomedical engineering 36 (2016) 172–181
Table 7 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 1 (mm). Valid N X Y
24 24
Mean 0.363746 0.214581
Median 0.360879 0.155869
Minimum 1.29951 1.04009
Maximum
Variance
Std. dev.
Standard error
0.399901 1.210011
0.270316 0.305367
0.519919 0.552600
0.106128 0.112799
Table 8 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 2 (mm). Valid N X Y
24 23
Mean 0.804095 0.030057
Median 0.770134 0.405123
Minimum 2.64936 2.05593
Maximum
Variance
Std. dev.
Standard error
1.022455 3.304122
0.811303 2.396157
0.900724 1.547953
0.183859 0.322770
Table 9 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 3 (mm). Valid N X Y
24 24
Mean 0.443053 0.189732
Median 0.277745 0.185959
Minimum 2.44994 2.45246
Maximum
Variance
Std. dev.
Standard error
1.056277 3.300190
0.849222 2.141073
0.921533 1.463241
0.188107 0.298683
Table 10 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 4 (mm). Valid N X Y
24 24
Mean 0.258754 0.328384
Median 0.002746 0.092647
Minimum 0.67487 2.65817
Maximum
Variance
Std. dev.
Standard error
2.062212 1.126076
0.636094 0.631663
0.797555 0.794773
0.162800 0.162232
Table 11 – Statistical analysis for deviation between the predictive and the real values of the coordinates X, Y and Z for additional points between Curves 3 and 4 (mm). Valid N X Y Z
66 66 72
Mean 0.158863 0.198208 0.044425
Median 0.163819 0.068180 0.037413
Minimum 1.55211 2.44699 2.04981
Maximum
Variance
Std. dev.
Standard error
2.915473 2.576518 1.359783
0.463297 0.873151 0.254633
0.680660 0.934426 0.504612
0.083783 0.115020 0.059469
Table 12 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 1 (mm). Valid N X Y
24 24
Mean 0.363746 0.214581
Median 0.360879 0.155869
Minimum 1.29951 1.04009
Maximum
Variance
Std. dev.
Standard error
0.399901 1.210011
0.270316 0.305367
0.519919 0.552600
0.106128 0.112799
Table 13 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 2 (mm). Valid N X Y
24 23
Mean 0.804095 0.030057
Median 0.770134 0.405123
Minimum 2.64936 2.05593
Maximum
Variance
Std. dev.
Standard error
1.022455 3.304122
0.811303 2.396157
0.900724 1.547953
0.183859 0.322770
Table 14 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 3 (mm). Valid N X Y
24 24
Mean 0.443053 0.189732
Median 0.277745 0.185959
Minimum 2.44994 2.45246
Maximum
Variance
Std. dev.
Standard error
1.056277 3.300190
0.849222 2.141073
0.921533 1.463241
0.188107 0.298683
179
biocybernetics and biomedical engineering 36 (2016) 172–181
Table 15 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 4 (mm). Valid N X Y
24 24
Mean 0.258754 0.328384
Median
Minimum
0.002746 0.092647
0.67487 2.65817
(b) X and Y coordinates for points at Curve 2 – all models are logarithmic; (c) X coordinates for points at Curve 3 – 1 exponential model, 1 quadratic, while all others are logarithmic, the same for Y coordinates; (d) X coordinates for points at Curve 4 – 3 exponential models, Y coordinates – 2 exponential models, all others are logarithmic; (e) Models for Z coordinates on all 4 curves are logarithmic. The values of X, Y and Z coordinates for additional points depend from parameters d3 and d4, in such manner that the points which are nearer to Curve 3 show dependency from d3, while the points that are nearer to Curve 4 show dependency from d4. Certain variability can be observed in regard to regression models. The models for X coordinates are mostly logarithmic, with few exceptions (1 quadratic and 5 linear). For Y coordinates, the situation is slightly different – 1 exponential, 4 linear and all the rest are logarithmic, while all 72 models for Z coordinate are logarithmic. Based on the results presented above we can conclude that our method has sufficient accuracy and in this phase of development is suitable for application in practice. In the final stage of modeling, in order to additionally implement accuracy of our predictive models, 2 types of analysis in CATIA Shape module for both obtained models are carried out. First is analysis of deviations between the predictive surface models obtained by loft function and the upper extremity of the ‘‘control’’ bone (with accuracy of 0.01 mm). It has been conducted and presented in Fig. 5. The
Maximum
Variance
Std. dev.
Standard error
2.062212 1.126076
0.636094 0.631663
0.797555 0.794773
0.162800 0.162232
mean deviation is 0.789 mm, the standard deviation is 1.05 mm and 74.53% of analyzed points have deviations in 1.46 up to 1.8 mm. Second is analysis of range from deviations between the upper extremity of the control bone the predictive model obtained in the form of point cloud (with accuracy of 0.01 mm). It has been conducted and presented in Fig. 6. The mean deviation is 2.27 mm, the standard deviation is 0.951 mm, but only 33.33% of analyzed points have deviations in range from 0 up to 1.81 mm. Although the standard deviation is less in the second case, the main conclusion is that the model obtained by loft function
Fig. 6 – Verification of the second model by analysis of deviation.
Fig. 5 – Verification of the first model by analysis of deviation.
Fig. 7 – Verification of the first model by analysis of the distance.
180
biocybernetics and biomedical engineering 36 (2016) 172–181
Financial support This research was conducted as a part of the project III 41017 Virtual Human Osteoarticular System and its Application in Preclinical and Clinical Practice is sponsored and supported by the Ministry of Education, Science and Technological Development of the Republic of Serbia, for the period of 2011–2015.
Acknowledgements
Fig. 8 – Verification of the second model by analysis of the distance.
The paper is part of the project III 41017 Virtual Human Osteoarticular System and its Application in Preclinical and Clinical Practice sponsored by the Ministry of Education, Science and Technological Development of the Republic of Serbia, for the period of 2011–2015.
references is better. Analysis in the form of measuring distance in perpendicular direction, without modifying elements, between the initial model of upper fibular extremity and both predictive models is also performed. The results from this 3D comparison are presented in terms of color coded map (by use of min/max option) of distances are shown in Figs. 7 and 8. In the first case, most of the area of the model (67.17%) has been built within 0–1 mm, while in the second case value for this area is 62.04%, although the maximum value is slightly lower in the second case.
5.
Conclusions
The proposed predictive models of a bone can be used for a quick and efficient creation of a 3D model of human long bones. This is especially true for the model obtained from splines with the aid of loft function, because of the values of deviations and distances between initial model and obtained model. This model also better represents the geometry of the upper end of fibula bone. Such 3D models can be widely used in computer assisted preoperative planning in orthopedic surgery, as well as for designing and manufacturing of personalized fixators, implants and scaffolds. Based on the distances analysis results, we have concluded that accuracy of the 3D model is acceptable, but there is room for improvement, especially in the area of lateral part. It can be achieved by increasing the number of samples. This research provides a starting point for the development of the predictive model of fibula as a whole. In order to obtain prediction models faster using CATIA software, there is a need to write VBSCRIPT and CATSCRIPT routines for the input of the values of the measured parameters (d1–d4), for prediction values for the coordinate points, as well as for linking results in CATIA Shape module.
[1] Ritter L. 3D interactive segmentation – first applications for computer-aided craniofacial surgical planning. Technischen Universität München; 2005. [2] Grey H. Gray's Anatomy: The Anatomical Basic of Clinical Practice. 40th ed. Elsevier; 2008. [3] Radakovich M, Malone T. The superior tibiofibular joint: the forgotten joint. J Orthop Sports Phys Ther 1982;3:129–32. [4] Radu C. 3D modeling and static finite element analysis of human tibia. Adv Eng 2008;2:99–104. [5] Radu C, Roşca I. Some contributions to the design of osteosynthesis implants. Est J Eng 2009;15:121. http://dx.doi.org/10.3176/eng.2009.2.05 [6] Mojtaba K, Ke BG, LePing L. Load response to articular cartilage and ligaments valus loading – a fibril-reinforced model of the knee. Ann Meet Am Soc 2010;415. [7] Faustini MC, Neptune RR, Crawford RH. The quasi-static response of compliant prosthetic sockets for transtibial amputees using finite element methods. Med Eng Phys 2006;28:114–21. http://dx.doi.org/10.1016/j.medengphy.2005.04.019 [8] Ozkan A, Mutlu I, Cirpici IY, Buluc L, Muezzinoglu US, Kisioglu Y. Effects of fibula and talus on the tibial stress distribution. EnginSoft Int Conf. 2009. pp. 1–4. [9] Lethaus B, Kessler P, Boeckman R, Poort LJ, Tolba R. Reconstruction of a maxillary defect with a fibula graft and titanium mesh using CAD/CAM techniques. Head Face Med 2010;6:16. http://dx.doi.org/10.1186/1746-160X-6-16 [10] Trajanović M, Tufegdžić M, Arsić S, Veselinović M, Vitković N. Reverse engineering of human fibula. 11th Int. Sci. Conf., Novi Sad, Serbia: Faculty of Technical Sciences, Novi Sad. 2012. pp. 527–30. [11] Majstorovic V, Trajanovic M, Vitkovic N, Stojkovic M. Reverse engineering of human bones by using method of anatomical features. CIRP Ann Manuf Technol 2013;62:167–70. http://dx.doi.org/10.1016/j.cirp.2013.03.081 [12] Stojkovic M, Milovanovic J, Vitkovic N, Trajanovic M, Arsic S. Analysis of femoral trochanters morphology based on geometrical model. J Sci Ind Res (India) 2012;71:210–6. [13] Vitković N, Milovanović J, Korunović N, Trajanović M, Stojković M, Mišić D, et al. Software system for creation of
biocybernetics and biomedical engineering 36 (2016) 172–181
human femur customized polygonal models. Comput Sci Inf Syst 2013;10:1473–97. http://dx.doi.org/10.2298/CSIS121004058V [14] Marciniec A, Miechowicz S. Stereolitography – the choice for medical modelling. Acta Bioeng Biomech 2004;6:13–24. [15] Zhang S, Zhang K, Wang Y, Feng W, Wang B, Yu B. Using three-dimensional computational modeling to compare the
181
geometrical fitness of two kinds of proximal femoral intramedullary nail for Chinese femur. Sci World J 2013;2013:6. http://dx.doi.org/10.1155/2013/978485 [16] Tufegdžić M, Trajanović M, Vitković N, Arsić S. Reverse engineering of the human fibula by the anatomical features method. FACTA Univ Ser Mech Eng 2013;11:133–9.
Available online at www.sciencedirect.com
ScienceDirect journal homepage: www.elsevier.com/locate/bbe
Original Research Article
Predictive geometrical model of the upper extremity of human fibula Milica Tufegdzic a,*, Stojanka Arsic b, Miroslav Trajanovic c a
Department for Mechanical Engineering, Machine-electrotechnical School, Krusevac, Serbia Department of Anatomy, Faculty of Medicine Nis, University of Nis, Nis, Serbia c Department for Production, IT and Management, Faculty of Mechanical Engineering, University of Nis, Nis, Serbia b
article info
abstract
Article history:
Computer assisted preoperative planning in orthopedic surgery, as well as designing and
Received 27 June 2015
manufacturing of personalized fixators, implants and scaffolds requires a good three-
Received in revised form
dimensional model of bone(s) of the treated patients. Existing methods that convert the
30 November 2015
Computer Tomography (CT) images into the polygonal three-dimensional models are time-
Accepted 11 December 2015
consuming and inefficient. Therefore, we propose a predictive model that allows quick
Available online 29 December 2015
creation of three-dimensional (3D) surface model of a particular bone by measuring the relevant parameters from an X-ray or CT image.
Keywords: Geometrical model
In this paper, we present the process of creating a predictive geometrical model using the case of proximal end of fibula as an example. The predictive model is built by defining the
Parameters
referential geometric entities that correspond to anatomical features, based on which
Fibula
appropriate points, axes, planes and curves are created. Using the method of linear and
Prediction
nonlinear regression with four different parameters, which can be measured from X-ray
Regression model
images or anterior-posterior projection of fibula at CT scans, the equations for X, Y and Z coordinates of the selected 168 points are obtained and their predictive values are calculated. These values are used for creating 3D surface model with the aim of two different methods: using loft function and converting these coordinates into point cloud. These models were compared and verified through analysis of deviations and distances between initial model and predictive models. The resulting 3D model has satisfactory accuracy, and the process of its building is much shorter. # 2015 Nałęcz Institute of Biocybernetics and Biomedical Engineering. Published by Elsevier Sp. z o.o. All rights reserved.
* Corresponding author at: Department for Mechanical Engineering, Machine-electrotechnical School, Cirila i Metodija 26, 37000 Krusevac, Serbia. E-mail addresses: [email protected] (M. Tufegdzic), [email protected] (S. Arsic), [email protected] (M. Trajanovic). http://dx.doi.org/10.1016/j.bbe.2015.12.003 0208-5216/# 2015 Nałęcz Institute of Biocybernetics and Biomedical Engineering. Published by Elsevier Sp. z o.o. All rights reserved.
biocybernetics and biomedical engineering 36 (2016) 172–181
1.
Introduction
Three-dimensional (3D) models with accurate geometric representation of the long bones are being increasingly used for various aspects of clinical practice and research. They are an indispensable basis for preoperative planning in orthopedic surgery, designing and manufacturing of personalized fixators, implants and scaffolds. There are many programs that allow obtaining the polygonal 3D models of bones based on Computer Tomography (CT) images in DICOM format. Independently from the program, building a 3D surface model of a bone from CT scans must be conducted through several steps: data management, image enhancement, segmentation and 3D reconstruction [1]. Defining the boundaries of the bone is necessary for healing of the model, but it often cannot be done accurately, the process is time consuming and requires a great deal of patience. It is well known that in clinical practice the patient cannot undergo a CT scan because of previous exposure to very high levels of ionizing radiation. In such cases it is impossible to build a 3D model of a bone. Our idea is to create a method for obtaining a predictive model of long bones from X-ray or CT images, using Computer Aided Design (CAD) program and statistical tools. In this paper, we are presenting the method using the upper extremity of human fibula as an example. According to the anatomical description, the human fibula is a long, slender bone, located at the lateral aspect of the leg. It has three main parts: a proximal head, a narrow neck, a long shaft and a distal lateral malleolus [2]. The upper extremity of the fibula is composed of two parts, head and neck. The head of the fibula has a most prominent part called the apex and a flattened articular surface, which forms the proximal tibiofibular joint with a corresponding fibular articular surface on the lateral condyle of the tibia. Articular capsule is reinforced by the anterior and posterior tibio-fibular ligament of the fibular head. Main parts of human fibula are presented in Fig. 1 (down below). It is important to point out, as an additional biomechanical and functional description of the fibular upper end, that there are two types of superior tibio-fibular joints, horizontal and oblique, due to the anatomical variations [3]. They provide a compensatory motion in internal and external rotational movements of the tibia. The primary functions of the superior tibio-fibular joint are: (1) the dissipation of torsional stresses which occur at the ankle joint, (2) wasting of the lateral tibial flexing, and (3) tensile weight bearing [3]. The previous work on the 3D modeling of the fibula or the parts of the fibula has resulted in generating a 3D model of fibula bone ankle joint (including the lower end of the fibula) by using software for processing the medical images and the CAD software. The process of the 3D modeling of human fibula using medical image software is described in [4]. IGES surface of the human fibula using MIMICS software and a 3D surface model in Solid Works, in order to test a newly developed and adapted osteosynthesis implant, used for transsindesmotic fibula's fracture is obtained in [5]. In another paper a finite element model of the knee joint, including distal femur, tibia, fibula, menisci, articular cartilages and ligaments is obtained from MR images of a healthy male [6]. The models of different
173
bones such as distal femur, patella and upper extremities of fibula and tibia are obtained from the CT data. These CT slices were converted into the point clouds (using Mimics) and a CAD surface model was obtained in Geomagic. These CAD models were further processed for the purpose of Finite Element Analysis [7]. In order to study the influence of the fibula and talus on the distribution of stress on the tibia, a 3D solid model of the fibula was created using CT images for generating 3D solid models of tibia, fibula and talus, with the help of Mimics in the form of point clouds [8]. These point clouds were exported in Geomagic where the failure fixing and smoothing on the surfaces was performed. The 3D models were saved in STL format and later assembled using Mimics FEA module. A case of reconstruction of a maxillary defect and orbital floor with a micro vascular fibula graft is shown in [9]. An individualized titanium mesh, as a rapidly prototyped template for the missing part of fibula, is designed from a CT scan using Materialize. 3D surface model of the fibula was obtained in the process of Reverse Engineering (RE) [10], but this process was timeconsuming, as all processes mentioned above.
2.
Materials and methods
In the process of generating a predictive model of the proximal end of human fibula, a method of anatomical features (MAF) was used. This method was previously described in [11] and has been implemented on femur and tibia bones. In the case of fibula bone, first few steps of this method were applied in our research. But, because of the fact that linear regression was not suitable for fibula bone, the method was expanded by introducing different nonlinear regression models. These models were tested by statistical tools and for the level of significance a value less then 0.05 was adopted. In order to determine the most suitable model, from those which were statically significance, we have chosen one bone from the sample for which the value of the length of mechanical axis is the closest to the mean value of the sample. Regression values of coordinates for corresponding points on this bone were calculated for each model and compared with initial values. So, for the models were adopted the ones for which the difference from regression values was the least. Finally, we have created predictive geometrical models from the regression coordinates, with the aim to compare and verify our results using the deviation and distance analysis. So, our research was conducted through the following steps: CT scanning and digitalization of data, definition and determination of anatomical landmarks on the upper and lower end of fibula, creation of the referential geometric entities (RGEs) (points, axes, planes, surfaces) according to the anatomical and morphological characteristics of the bone [12,13], selection and measurement of the parameters on the upper end of the fibula, identifying and measuring the coordinates of the points for which the regression models are formed, mathematical models determination of the regression equations,
174
biocybernetics and biomedical engineering 36 (2016) 172–181
calculating the regression values for X, Y and Z coordinates, calculating the difference between regression and initial values for ‘‘control’’ bone and adoption of the best regression model, creation of the models using loft function and converting point coordinates into point cloud, and verification and comparison between the initial and obtained models, by using the deviation and distance analysis.
2.1.
CT scanning and digitalization of data
The two available scanning methods for quantitative 3D imaging are CT and MRI. Both CT and MRI provide accurate information for quantifying of the anatomical structures in a 3D environment. Both imaging methods have considerable application in clinical praxis, but CT is mainly used for the bone imaging and MRI for the soft tissue imaging. To create a predictive model of the upper extremity of the human fibula, a series of traverse CT images was taken from 11 patients (6 men and 5 women), 43–70 years old (58 years on average). These images were obtained on a Toshiba MSCT scanner Aquillion 64 (120 kV, 150 mA, thickness 1 mm, inplane resolution 0.781 0.781 (pixel size), acquisition matrix 512 512; field of view (FOV) 400 mm 400 mm). All the CT data were saved in DICOM format. After digital processing, these 2D data are converted to 3D STL (STereoLitography) models [14], considering that 3D modeling is more precise than 2D methods because it is less influenced by body position [15]. These models were imported in the CAD software and then further processed by eliminating the model errors (isolated triangles, non-manifold vertices, edges, etc.) and removing all unnecessary entities, healing, mesh creation, optimization and mesh smoothing. Model obtained in this manner was used for determining the RGEs in a CAD program.
2.2.
A-P (anterior–posterior) plane – defined by mechanical axis and ACF. ACF, AML, mechanical axis and the main parts of human fibula are presented in Fig. 1. The next stage was determining and measuring the 4 parameters in A-P plane (d1, d2, d3 and d4). These parameters are shown in Fig. 2: d1 – distance from the highest point of FACF in the A-P plane to the most distal point on the upper end of the fibula, measured in the direction perpendicular to the mechanical axis,
Methods
For the purpose of our research, RGEs at all fibula bones were identified. It was done according to the definitions described in [10,16]: (1) at the upper extremity (lat. epiphysis proximalis s. extremitas proximalis): ACF (lat. apex capitis fibulae) – apex of the fibular head, FACF (lat. facies articularis capitis fibulae) – surface for articulation with the lateral condyle of tibia, (2) at the lower extremity (epiphysis distalis s. extremitas distalis): FAML (lat. facies articularis malleoli lateralis) – surface for articulation with the talus faced medially, AML (lat. apex malleoli lateralis) – top of the lateral malleolus, as the most distal point on the lateral malleolus. At the basis of the previously defined RGEs the following entities were defined: the mechanical axis – which connects the centers of the articular surfaces at the upper (FACF) and the lower extremity (FAML) of fibula, and
Fig. 1 – Main parts and RGEs at the right human fibula, medial aspect.
biocybernetics and biomedical engineering 36 (2016) 172–181
Fig. 2 – Parameters on the upper extremity of the right human fibula.
d2 – distance from the center of FACF and the most distal point on the upper end of the fibula, measured in the direction perpendicular to the mechanical axis, d3 – distance from the lowest point on FACF in the A-P plane to the most distal point on the upper end of the fibula, measured in the direction perpendicular to the mechanical axis, d4 – the distance at the level (the intersection) where the head goes into the neck of the fibula, measured in the direction perpendicular to the mechanical axis. With the aim to get coordinates of points and their values, we have used two sets of curves, obtained with two sets of planes. These planes were utilized for intersections of polygonal model of the upper end of fibula. Measured values were applied for obtaining predictive values of X, Y and Z coordinates. The first set of planes is perpendicular to mechanical axis and set at places of the selected parameters (d1, d2, d3 and d4). So, the polygonal models of the upper end of fibula are intersected with these planes. As result, we got a corresponding set of curved lines (Curve 1, 2, 3 and 4, Fig. 3). The second set of planes was constructed at certain angles in relation to the A-P plane (158, 308, 458, 608, 758, 908, . . .), where a line parallel to the mechanical axis which passes through the ACF was used as the axis of rotation, according to our previous research [10]. The sections of the polygonal model of the upper end of fibula with these planes (including A-P plane), have resulted in corresponding set of 12 rotational curves, with ACF as common point. These two sets of curves always have intersections in two points and we have got 96 points. Because of the fact that the distance between the parameters d3 and d4 was great and we could not get satisfactory accuracy of the model, we have added 6 points on every curve from the second
175
Fig. 3 – Curves and points of intersections at the upper end of fibula.
set of planes. These points were chosen at the same distance (5, 10 and 15 mm) from the points of intersections of the curve 3 (Fig. 3) and the curves from the second sets of planes. At Fig. 3 are presented points of intersections between curves obtained at the section of polygonal model of the upper end with the A-P plane and the planes which are rotated for the angles of 158, 308, 458, 608, as an examples (respectively identified as A-P curve, 158 curve, 308 curve, 458 curve, 608 curve), and Curves 1, 2, 3 and 4. Points of intersection of rotational curves are marked as A-P-1, A-P-2 on Curves 1, 2, 3 and 4. The same methodology is applied for labeling corresponding points at the other rotational curves (15-1, 15-2, etc.), while additional points between Curves 3 and 4 are labeled as A-P-1-1, A-P-1-2, A-P-1-3, A-P-2-1, A-P-2-2, A-P-2-3, 15-1-1, and so on (for the simplicity of Fig. 3 only these points are marked). The coordinates of points of intersections in the newly formed coordinate system of the right orientation were measured. Because of the fact that ACF is the common point for the second set of curves, it was adopted for the origin, the X-axis is perpendicular to the A-P plane, Y-axis ‘‘lies’’ in the A-P plane (perpendicular to the mechanical axis), while the Z-axis is parallel to the mechanical axis. So, the sets of X, Y and Z coordinates were obtained and their values were measured. These sets of coordinates represent specific locations on the upper end of the fibula. The values of these coordinates for 10 fibula bones are the input data for the formation of the predictive model, while the 11th bone is separated and used as a ‘‘control’’ fibula for determination of the most suitable regression model. The main criterion for the selection of ‘‘control’’ fibula was mechanical axis length. It was selected in such manner that its length has the smallest deviation from the mean value of mechanical axis length of the sample.
176
biocybernetics and biomedical engineering 36 (2016) 172–181
Values of measured coordinates (Xi, Yi, Zi) were the input data for the statistical analysis, treated as dependent variables, while the values of the parameters d1, d2, d3 and d4 were treated as independent variables. The relationships between these two sets of variables were established by means of correlation coefficients and scatter plots were used for choosing appropriate regression models. Based on this analysis, 6 regression models (1 linear and 5 non-linear) were chosen and tested for every predicted value of the coordinates: (1) Linear model represented by the expression Xi, Yi, Zi = a + b di, (2) Non-linear models: (a) quadratic models, represented by the expressions Xi, Yi, Zi = a + b di2 and Xi, Yi, Zi = a + b di + c di2; (b) 2 logarithmic models represented by the expressions Xi, Yi, Zi = a + b ln(di) and Xi, Yi, Zi = a ln(di); (c) exponential model represented by the expressions Xi, Yi, Zi = a exp(b di); i = 1, 2, 3, 4; a, b – coefficients of linear regression, while di stands for the values of the parameters d1–d4. Coefficients in regression equations were calculated using statistical software and taking into account that the value of statistical significance p is p < 0.05. Models that have p-values greater than 0.05 for any of the coefficients in the regression equations were discarded. From the models in which the p-value was less than 0.05 was elected the one with a minimal deviation of regression (predictive) values of the X, Y and Z coordinates and the measured coordinate values for corresponding points at the ‘‘control’’ fibula. For additional points between Curves 3 and 4 the same statistical procedure was used, but some different models were tested because the values of the coordinates X, Y and Z depends from both parameters d3 and d4, according to correlation matrix. These models are: (1) linear, represented by the expressions (a) Xi, Yi, Zi = a + b di, i = 3, 4; (b) Xi, Yi, Zi = a + b d3 + c d4; (2) non-linear, represented by the expressions (a) Xi, Yi, Zi = a + b di2, i = 3, 4Xi, Yi, Zi = a + b d3 + c d42, Xi, Yi, Zi = a + b d4 + c d32; (b) Xi, Yi, Zi = a + b ln(di), Xi, Yi, Zi = a ln(di), i = 3, 4; (c) Xi, Yi, Zi = a exp(b di), i = 3, 4. The method of election of the appropriate model was the same as for the points of intersection, described above. Calculated predictive values of the X, Y and Z coordinates are used for constructing sets of points. Splines were constructed from matching points and presented in Fig. 4. By the loft function in CATIA Shape module surface predictive model was created. These sets of points were also transformed into the point cloud, from which the mesh could be created. These two models are further used for comparison and verification through the analysis of deviations and distances between initial polygonal model and obtained models.
Fig. 4 – Splines through regressive points.
3.
Results
Regression models of the points coordinates (Fig. 3, points A-P-1 and A-P-2 at Curves 1, 2, 3 and 4) are presented by equations in Tables 1–4. Different regression models are obtained for point coordinates X and Y at the Curves 1–4, while the regression equation for Z coordinate remains the same for all points at selected curve, because of the fact that proximal fibular end is cut by the planes perpendicular to the mechanical axis and Z axis which is parallel to the mechanical axis. Regression models of the points coordinates X, Y and Z between Curves 3 and 4 are presented in Tables 5 and 6 (as examples). For the Y coordinate for the points at the A-P curve the values for coefficients of regression are not calculated, because all the values are 0. Another statistical analysis of the deviations between the predictive coordinate values and actual (real) coordinate values was conducted. The results obtained by means of descriptive statistics are given in Tables 7–10 for the X and Y coordinates (for Curves 1–4), where Valid N represents the number of predictive (calculated) values. Analysis for Z coordinate was not conducted because there is a unique value at each curve. Based on the values of means and standard errors presented in Tables 7–10 we have concluded that all deviations between predicted and real values of coordinates can be treated as very good. That means our predicted values are very close to the real values. Similar analysis of the deviations was conducted for deviations of all X, Y and Z coordinates for additional points. The results are presented in Table 11, where Valid N represents the number of predictive (calculated) values. Some exceptions should be noticed in this case, especially for Y coordinate. The deviations could be reduced and regression models improved by increasing the number of samples. Another statistical analysis of the deviations between the predictive coordinate values and actual (real) coordinate
177
biocybernetics and biomedical engineering 36 (2016) 172–181
Table 1 – Regression equations for X, Y and Z coordinates of selected points at the parameter d1. Coord.
Points at Curve 1 A-P-1 X = 1.39372 exp(0.13076 d1) Y = 0.722038 + 0.157627 d1 0.008025 d1 2
X Y Z
X = 1.714685 ln(d1) Y = 0.495697 + 0.105480 d1
Discussion
where the first set of planes intersect polygonal model of the upper end of fibula (Curve 1 – parameter d1, etc.). Some other regularity can be seen from the obtained equation regression for: (a) X coordinates for points at Curve 1 – 3 exponential models, Y coordinates – 2 quadratic and 2 exponential models, all others are logarithmic;
Table 5 – Regression equations for X, Y and Z coordinates of points between the parameters d3 and d4 at the A-P curve. Coord.
Different dependences in the regression equations can be treated as a result of the variability of the upper fibular end (regarding to the sex, different types of the FACF and different geometry of the samples). But, some regularity can be noticed in regression models. The values of X, Y and Z coordinates for the points of intersections at Curves 1, 2, 3 and 4 directly depend from parameters which were measured at the places
X Z
Coord.
Table 2 – Regression equations for X, Y and Z coordinates of selected points at the parameter d2. Points at Curve 2
Coord. A-P-1
A-P-2
X = 3.89353 ln(d2) X = 1.410168 ln(d2) Y = 0.000436 ln(d2) Y = 0.004226 ln(d2) Z = 3.072500 ln(d2)
Table 3 – Regression equations for X, Y and Z coordinates of selected points at the parameter d3. Coord.
Points at Curve 3 A-P-2
X Y Z
A-P-2
X = 5.8614177 ln(d3) X = 28.1574 + 9.8685 ln(d3) Y = 0.006120 ln(d3) Y = 0.000389 ln(d3) Z = 4.069104 ln(d3)
Coord.
X Z
Coord.
X Y Z
Coord.
X Y Z
Coord.
Coord.
A-P-1 X Y Z
A-P-2
X = 7.75035 exp(0.03694 d4) X = 0.554902 ln(d4) Y = 0.002218 ln(d4) Y = 0.000945 ln(d4) Z = 10.334974 ln(d4)
A-P-1-1
A-P-2-1
X = 6.04468 ln(d3) Z = 6.134285 ln(d4)
X = 23.2087 + 9.2082 ln(d4) Z = 6.264310 ln(d4)
Points at the distance of 10 mm A-P-2-2
X = 9.77308-0.02950 d4 Z = 7.791366 ln(d4)
2
X = 0.882985 ln(d4) Z = 7.970058 ln(d4)
Points at the distance of 15 mm A-P-1-3
A-P-2-3
X = 6.71042 ln(d4) Z = 8.232533 ln(d3)
X = 34.7029 + 12.4725 ln(d4) Z = 9.438726 ln(d4)
Table 6 – Regression equations for X, Y and Z coordinates of points between the parameters d3 and d4 at the 158 curve.
Table 4 – Regression equations for X, Y and Z coordinates of selected points at the parameter d4. Points at curve 4
Points at the distance of 5 mm
A-P-1-2 X Z
X Y Z
0.005354 d1 2
Z = 1.270169 ln(d1)
values was conducted. The results obtained by means of descriptive statistics are given in Tables 7–10 for the X and Y coordinates (for Curves 1–4), where valid N represents the number of predictive (calculated) values. Analysis for Z coordinate was not conducted because there is a unique value at each curve (Tables 12–15). Based on the values of means and standard errors presented in Tables 7–10 we have concluded that all deviations between predicted and real values of coordinates can be treated as very good. That means our predicted values are very close to the real values.
4.
A-P-2
X Y Z
Points at the distance of 5 mm 15-1-1
15-2-1
X = 5.47192 ln(d3) Y = 1.79856 exp(0.05366 d4) Z = 6.099127 ln(d4)
X = 1.749891 ln(d3) Y = 0.268283 ln(d3) Z = 6.574582 ln(d4)
Points at the distance of 10 mm 15-1-2
15-2-2
X = 5.49791 ln(d3) Y = 1.65101 ln(d4) Z = 7.719516 ln(d4)
X = 1.256039 ln(d3) Y = 0.220073 ln(d3) Z = 7.930495 ln(d4)
Points at the distance of 15 mm 15-1-3
15-2-3
X = 5.28828 ln(d3) Y = 1.41699 ln(d3) Z = 9.246085 ln(d4)
X = 48.8380 + 15.6115 ln(d3) Y = 31.9905 12.5535 ln(d3) Z = 8.965555 ln(d4)
178
biocybernetics and biomedical engineering 36 (2016) 172–181
Table 7 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 1 (mm). Valid N X Y
24 24
Mean 0.363746 0.214581
Median 0.360879 0.155869
Minimum 1.29951 1.04009
Maximum
Variance
Std. dev.
Standard error
0.399901 1.210011
0.270316 0.305367
0.519919 0.552600
0.106128 0.112799
Table 8 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 2 (mm). Valid N X Y
24 23
Mean 0.804095 0.030057
Median 0.770134 0.405123
Minimum 2.64936 2.05593
Maximum
Variance
Std. dev.
Standard error
1.022455 3.304122
0.811303 2.396157
0.900724 1.547953
0.183859 0.322770
Table 9 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 3 (mm). Valid N X Y
24 24
Mean 0.443053 0.189732
Median 0.277745 0.185959
Minimum 2.44994 2.45246
Maximum
Variance
Std. dev.
Standard error
1.056277 3.300190
0.849222 2.141073
0.921533 1.463241
0.188107 0.298683
Table 10 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 4 (mm). Valid N X Y
24 24
Mean 0.258754 0.328384
Median 0.002746 0.092647
Minimum 0.67487 2.65817
Maximum
Variance
Std. dev.
Standard error
2.062212 1.126076
0.636094 0.631663
0.797555 0.794773
0.162800 0.162232
Table 11 – Statistical analysis for deviation between the predictive and the real values of the coordinates X, Y and Z for additional points between Curves 3 and 4 (mm). Valid N X Y Z
66 66 72
Mean 0.158863 0.198208 0.044425
Median 0.163819 0.068180 0.037413
Minimum 1.55211 2.44699 2.04981
Maximum
Variance
Std. dev.
Standard error
2.915473 2.576518 1.359783
0.463297 0.873151 0.254633
0.680660 0.934426 0.504612
0.083783 0.115020 0.059469
Table 12 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 1 (mm). Valid N X Y
24 24
Mean 0.363746 0.214581
Median 0.360879 0.155869
Minimum 1.29951 1.04009
Maximum
Variance
Std. dev.
Standard error
0.399901 1.210011
0.270316 0.305367
0.519919 0.552600
0.106128 0.112799
Table 13 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 2 (mm). Valid N X Y
24 23
Mean 0.804095 0.030057
Median 0.770134 0.405123
Minimum 2.64936 2.05593
Maximum
Variance
Std. dev.
Standard error
1.022455 3.304122
0.811303 2.396157
0.900724 1.547953
0.183859 0.322770
Table 14 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 3 (mm). Valid N X Y
24 24
Mean 0.443053 0.189732
Median 0.277745 0.185959
Minimum 2.44994 2.45246
Maximum
Variance
Std. dev.
Standard error
1.056277 3.300190
0.849222 2.141073
0.921533 1.463241
0.188107 0.298683
179
biocybernetics and biomedical engineering 36 (2016) 172–181
Table 15 – Statistical analysis for deviation between the predictive and the real values of the coordinates X and Y at Curve 4 (mm). Valid N X Y
24 24
Mean 0.258754 0.328384
Median
Minimum
0.002746 0.092647
0.67487 2.65817
(b) X and Y coordinates for points at Curve 2 – all models are logarithmic; (c) X coordinates for points at Curve 3 – 1 exponential model, 1 quadratic, while all others are logarithmic, the same for Y coordinates; (d) X coordinates for points at Curve 4 – 3 exponential models, Y coordinates – 2 exponential models, all others are logarithmic; (e) Models for Z coordinates on all 4 curves are logarithmic. The values of X, Y and Z coordinates for additional points depend from parameters d3 and d4, in such manner that the points which are nearer to Curve 3 show dependency from d3, while the points that are nearer to Curve 4 show dependency from d4. Certain variability can be observed in regard to regression models. The models for X coordinates are mostly logarithmic, with few exceptions (1 quadratic and 5 linear). For Y coordinates, the situation is slightly different – 1 exponential, 4 linear and all the rest are logarithmic, while all 72 models for Z coordinate are logarithmic. Based on the results presented above we can conclude that our method has sufficient accuracy and in this phase of development is suitable for application in practice. In the final stage of modeling, in order to additionally implement accuracy of our predictive models, 2 types of analysis in CATIA Shape module for both obtained models are carried out. First is analysis of deviations between the predictive surface models obtained by loft function and the upper extremity of the ‘‘control’’ bone (with accuracy of 0.01 mm). It has been conducted and presented in Fig. 5. The
Maximum
Variance
Std. dev.
Standard error
2.062212 1.126076
0.636094 0.631663
0.797555 0.794773
0.162800 0.162232
mean deviation is 0.789 mm, the standard deviation is 1.05 mm and 74.53% of analyzed points have deviations in 1.46 up to 1.8 mm. Second is analysis of range from deviations between the upper extremity of the control bone the predictive model obtained in the form of point cloud (with accuracy of 0.01 mm). It has been conducted and presented in Fig. 6. The mean deviation is 2.27 mm, the standard deviation is 0.951 mm, but only 33.33% of analyzed points have deviations in range from 0 up to 1.81 mm. Although the standard deviation is less in the second case, the main conclusion is that the model obtained by loft function
Fig. 6 – Verification of the second model by analysis of deviation.
Fig. 5 – Verification of the first model by analysis of deviation.
Fig. 7 – Verification of the first model by analysis of the distance.
180
biocybernetics and biomedical engineering 36 (2016) 172–181
Financial support This research was conducted as a part of the project III 41017 Virtual Human Osteoarticular System and its Application in Preclinical and Clinical Practice is sponsored and supported by the Ministry of Education, Science and Technological Development of the Republic of Serbia, for the period of 2011–2015.
Acknowledgements
Fig. 8 – Verification of the second model by analysis of the distance.
The paper is part of the project III 41017 Virtual Human Osteoarticular System and its Application in Preclinical and Clinical Practice sponsored by the Ministry of Education, Science and Technological Development of the Republic of Serbia, for the period of 2011–2015.
references is better. Analysis in the form of measuring distance in perpendicular direction, without modifying elements, between the initial model of upper fibular extremity and both predictive models is also performed. The results from this 3D comparison are presented in terms of color coded map (by use of min/max option) of distances are shown in Figs. 7 and 8. In the first case, most of the area of the model (67.17%) has been built within 0–1 mm, while in the second case value for this area is 62.04%, although the maximum value is slightly lower in the second case.
5.
Conclusions
The proposed predictive models of a bone can be used for a quick and efficient creation of a 3D model of human long bones. This is especially true for the model obtained from splines with the aid of loft function, because of the values of deviations and distances between initial model and obtained model. This model also better represents the geometry of the upper end of fibula bone. Such 3D models can be widely used in computer assisted preoperative planning in orthopedic surgery, as well as for designing and manufacturing of personalized fixators, implants and scaffolds. Based on the distances analysis results, we have concluded that accuracy of the 3D model is acceptable, but there is room for improvement, especially in the area of lateral part. It can be achieved by increasing the number of samples. This research provides a starting point for the development of the predictive model of fibula as a whole. In order to obtain prediction models faster using CATIA software, there is a need to write VBSCRIPT and CATSCRIPT routines for the input of the values of the measured parameters (d1–d4), for prediction values for the coordinate points, as well as for linking results in CATIA Shape module.
[1] Ritter L. 3D interactive segmentation – first applications for computer-aided craniofacial surgical planning. Technischen Universität München; 2005. [2] Grey H. Gray's Anatomy: The Anatomical Basic of Clinical Practice. 40th ed. Elsevier; 2008. [3] Radakovich M, Malone T. The superior tibiofibular joint: the forgotten joint. J Orthop Sports Phys Ther 1982;3:129–32. [4] Radu C. 3D modeling and static finite element analysis of human tibia. Adv Eng 2008;2:99–104. [5] Radu C, Roşca I. Some contributions to the design of osteosynthesis implants. Est J Eng 2009;15:121. http://dx.doi.org/10.3176/eng.2009.2.05 [6] Mojtaba K, Ke BG, LePing L. Load response to articular cartilage and ligaments valus loading – a fibril-reinforced model of the knee. Ann Meet Am Soc 2010;415. [7] Faustini MC, Neptune RR, Crawford RH. The quasi-static response of compliant prosthetic sockets for transtibial amputees using finite element methods. Med Eng Phys 2006;28:114–21. http://dx.doi.org/10.1016/j.medengphy.2005.04.019 [8] Ozkan A, Mutlu I, Cirpici IY, Buluc L, Muezzinoglu US, Kisioglu Y. Effects of fibula and talus on the tibial stress distribution. EnginSoft Int Conf. 2009. pp. 1–4. [9] Lethaus B, Kessler P, Boeckman R, Poort LJ, Tolba R. Reconstruction of a maxillary defect with a fibula graft and titanium mesh using CAD/CAM techniques. Head Face Med 2010;6:16. http://dx.doi.org/10.1186/1746-160X-6-16 [10] Trajanović M, Tufegdžić M, Arsić S, Veselinović M, Vitković N. Reverse engineering of human fibula. 11th Int. Sci. Conf., Novi Sad, Serbia: Faculty of Technical Sciences, Novi Sad. 2012. pp. 527–30. [11] Majstorovic V, Trajanovic M, Vitkovic N, Stojkovic M. Reverse engineering of human bones by using method of anatomical features. CIRP Ann Manuf Technol 2013;62:167–70. http://dx.doi.org/10.1016/j.cirp.2013.03.081 [12] Stojkovic M, Milovanovic J, Vitkovic N, Trajanovic M, Arsic S. Analysis of femoral trochanters morphology based on geometrical model. J Sci Ind Res (India) 2012;71:210–6. [13] Vitković N, Milovanović J, Korunović N, Trajanović M, Stojković M, Mišić D, et al. Software system for creation of
biocybernetics and biomedical engineering 36 (2016) 172–181
human femur customized polygonal models. Comput Sci Inf Syst 2013;10:1473–97. http://dx.doi.org/10.2298/CSIS121004058V [14] Marciniec A, Miechowicz S. Stereolitography – the choice for medical modelling. Acta Bioeng Biomech 2004;6:13–24. [15] Zhang S, Zhang K, Wang Y, Feng W, Wang B, Yu B. Using three-dimensional computational modeling to compare the
181
geometrical fitness of two kinds of proximal femoral intramedullary nail for Chinese femur. Sci World J 2013;2013:6. http://dx.doi.org/10.1155/2013/978485 [16] Tufegdžić M, Trajanović M, Vitković N, Arsić S. Reverse engineering of the human fibula by the anatomical features method. FACTA Univ Ser Mech Eng 2013;11:133–9.