Texture Based Bone-Radiograph Image Analysis for the Assessment of Osteoporosis

Texture Based Bone-Radiograph Image Analysis for the Assessment of Osteoporosis

Kavya R

Dept. of Bio Medical signal processing and Instrumentation

B.M.S College of Engineering Bangalore, India

Dr.Joshi Manisha Shivram Asst.Professor Dept. of Bio Medical signal processing and Instrumentation

1. College of Engineering Bangalore, India

   Abstract-In recent years, there is a growing interest in finding the efficient diagnostic methods for skeletal system disease. Osteoporosis is one such kind of a skeletal system disease mainly categorized by alterations in bone mass density and its structure that makes the bone prone to fracture. The assessment of osteoporotic disease from bone radiograph images presents a major challenge for pattern recognition. Textured images of osteoporotic and healthy subjects show high degree of resemblance, it increases the difficulties in discriminating such textures.This is an IEEE open challenge and the images used in this research work have been provided by the university Dorleans. The database consists of 87 bone radiograph images of control subjects and 87 subjects with bone pathology.

   In this paper a new method has been proposed for the early detection of osteoporosis disease through 2D texture analysis. The results obtained will help in discriminating the osteoporotic patient from that of normal subject. Feed forward neural network classifier has been used for the discrimination and obtained 95% sensitivity.

   Keywords: Image Preprocessing; Feature Extraction; Neural Network Implementation.

   1. INTRODUCTION

      Texture Characterization of Bone radio- graph images (TCB) is a task in the osteoporosis diagnosis.Osteoporosis is mainly categorized by decrease in the bone mass and structure to increases the hazard of fracture. It is difficult to identify the normal subjects fromthat of osteoporotic because they show high degree of resemblance in their texture. The database mainly consists of 87 images of patient with osteoporotic fracture and 87 images are normal of subjects. It occurs most commonly in women because of hormonal changes (menopause).Generally occursmainly in spine, hips and wrist. By knowing the quality of the bone we can prevent the fracture risk.

      Figure 1: Example of textures of a normal Subject (left) and a patient with osteoporotic fractures (right).

      The assessment of discrimination of the texture images become a major problem for pattern recognition, because the textured images of both patient with osteoporotic fracture and normal subject shows high degree of resemblance, it increases the difficulties in classifying such textures. The main goal is to develop a procedure which will facilitate the early detection of the disease through 2D texture analysis.

      The proposed method is less expensive because X-ray images have been used rather than the magnetic resonance imaging and CT tomography requires expensive equipment. This is a simple diagnostic technique helps in the detection of osteoporosis.

   2. RELATED WORK

      Taleb-ahmed, dubois, duquenoy[1]developed the algorithm used for the discrimination of the textured CT images. The ultimate objective of this method is the discrimination of control subjects from that of subjects with pathology. They explained several methods of texture analysis. These methods are based on the fractal geometry whose claim to the analysis of texture images is recent. Here they have calculated fractal dimension from area not from volume.

      In this method they have calculated the fractal dimension, theoretical dimension and fractal signature .For the planned clinical use, the goal is not to measure the exact value of the fractal dimension, but mainly to control its variations according to the bone texture alteration. In this case a maximal range of variation of this fractal dimension can give best sensitivity in discriminating healthy subjects from that of bone pathologic patients.

      The unique reason for using fractals for the analysis is that, parameter (fractal dimension) derived from fractal geometry varies with the alteration of trabecular bone structure especially in osteoporotic patient and also thinning of the trabecular bone. In this method in order to estimate the fractal dimension they have selected the ROI (Region of Interest) manually. This ROI region is covered by using circular element. After obtaining the ROI they estimated the fractal dimension by using morphological covering method and variation method. It involves calculating the surface covering area, to determine this area two parameters must define they are lower surface and upper surface, these both should make the covering of the original surface. The concept called erosion and dilation have been used to calculate the lower and upper surface. The differences in the erosion and dilation values are summed for each pixel. The results obtained from these methods showed that there is a correlation

      between the results obtained from the fractal dimension and BMD amount by using DEXA.

      Ahmed Salmi Boumnini Hassani, Rachid Jennane, Mohammed Rziza1, Eric Lespessailles [2] developed a texture analysis method for the trabecular bone X-ray images. The main goal of this project is to study the effect of preprocessing the data of bone radiograph images for the diagnosis of osteoporosis. In the preprocessing step they enhanced the image by using Retinex algorithm, in the second step,these enhanced images are analyzed by using anisotropic morlet wavelet.

      The exploitation of the fully anisotropic morlet enabled them solving the problem of orientation which is caused by the non-uniform changes. In the third step,the renyi entropy has been used for the description of anisotropic textures.

      The renyi entropy results from the generalization of the entropy of shannon. It is an efficient tool which has shown good performances. After extracted the texture features they used Neural Network classifier for obtaining the good discrimination. The results obtained showed that this method is 92% efficient in classifying the patient with osteoporotic fracture from that of normal subject but they have not validated this result.

      Andrzej Materka [3]made an attempt to apply the digital image analysis technique for the detection of bone mass and its structure. Here the distal forearm bones are investigated. They have included the calibration phantom to improve the image intensity. They extracted first order texture parameters and fractal dimensions are evaluated. These derived texture features are correlated with the bone mineral density by using DEXA (Dual Energy X-ray Absorptiometry).In the methodology they have used image preprocessing to remove the noise as well as to extract the region of interest (ROI).Resultsobtained showed that by measuring the changes in statistical texture parameters and fractaldimensions of X-ray images it is possible to monitor changes in calcium contents and internal structure of the bone. Texture analysis shows potential usefulness in the diagnosis of skeletal diseases. This initial research was carried out by using first-order texture features only.

      Florian Yger [4] proposed a new method for texture analysis based on covariance matrices and wavelet marginal. In their work they focused mainly on covariance matrices and wavelet marginal rather than the complicated features.

      Covariance matrices have been studied as image descriptor in wide variety of applications from license plate detection to pedestrian detection, but in this case samples taken are more compared to the parameters (features) and they are more sensitive to its outlies. In order to overcome from this issue, they have used minimum covariance determinant and aims in giving the lower determinant for their experiment, they have used two variants of features they are, gradient based and gabor based and these covariance matrices belongs to a non-euclidean space where distances are not computed on straight line rather they computed on curves.

      Wavelet marginal are mainly based on wavelet decomposition mainly used to extract the frequential information for translation invariants. Once they defined the wavelet decomposition of a 2D image they extended it to the marginal, in their experiment marginal are used as baseline.

      In the methodology they resized the image from 400400 to 128128, 6464,256256, then they applied both the methods and the same procedure has been carriedoutfor the validation are compared, finally results showed that marginal haar wavelet on 128128 yields better discrimination between the healthy subjects from that of patient withosteoporotic fracture.

      Namita Aggarwal, R. K. Agrawal [5]proposed a new method mainly used to characterize the texture, by using statistical and wavelet based parameters, to discriminate the patient with alzheimer disease from that of normal subject. In order to do this they have used T2 weighted MRI brain image and also feature extraction technique for the classification. Then they made comparison for features derived from both methods.

      For the wavelet based feature extraction method, they have used daubechies wavelet coefficient function and also haar wavelet. After extracting all the features they fed into the SVM (Support Vector Classifier), and less accuracy has been found. In order to improve the accuracy reduced the extracted features by using a method called PCA (principal component analysis).The selected features are more efficient in giving the accurate results than all the extracted features. In the first and second order statistical feature extraction derived first order features such as mean, variance, skewness, and kurtosis and also derived the second order features from GLCM (Gray Level C-occurrence Matrix).

      Finally classification has done by using SVM classifier and calculated different statistical parameters to compare the results and finally they came to a conclusion that texture features from first and second order features gives significant performance than the wavelet based features.

   3. MATERIALS

      The images used in this research have been provided by the University Dorleans. The data consist of 1742D radiographic images in TIF format. In that 87 images are of control subjects and 87 are of pathological images. These images are 16 bit and their size is of 400400.

   4. METHODOLOGY

      In the proposed method 174 2D radiographic images have been used. In that 87 are normal subject images and another 87 are images of patient with osteoporotic fractures.

      Figure2: Methodology

      1. Image Preprocessing

         Preprocessing is one of the key step in texture analysis, the raw images are not free from the noise, this may cause variation in the statistical measures derived from the 2D radiographic images. In order to make the image free from noise and to enhance the quality of the image, preprocessing is the initial procedure in the development of image analysis algorithm.

         In our proposed algorithm wiener filter has been used to remove the random noise and histogram equalization has been performed to enhance the contrast. Image preprocessing involves noise removal and enhancement but the suitability of which task to be performed first can be judged from the statistical features calculated for the images.

         Figure3: Image Preprocessing

         In order to decide the best suitable flow for the preprocessing two approaches are considered.Figure3 represents the considered methods in the preprocessing. The first approach is to perform filtering first and then enhancement and the later one considers enhancement first and the filtering.

         The objective analysis for both the flows is as shown in Figure4 and Figure5 .Based on the power signalto noise ratio(PSNR),flow1 has been considered in the development of the algorithm.Statistical measures for both the flows are as shown in Table1.

         Figure4: Preprocessed Images by using filtering then enhancement approach

         (a) Original Image,(b) Increased contrast by Imadjust ,(c) Filtered image by using wiener filter,(d) Enhanced contrast of the image by using histogram

         equalization.

         Figure5: Preprocessed images by filtering then enhancement approach (a) Original image,(b)Increased contrast by Imadjust (c) Enhanced contrast of the image by using histogram equalization,(d) Filtered image by using wiener filter.

         PSNR Values when filtering operation is performed after the enhancement	PSNR values when filtering operation is performed before the enhancement
         30.8967	31.05428
         27.034	27.7251
         22.789	23.1061
         28.6432	29.9855
         23.0733	25.5308
         24.3945	24.9094
         27.048	28.4189

         Table1: PSNR Values for both the approaches

      2. Feature Extraction

         1. First and Second order statistical features

            1. First order texture features: Mean,median,variance,standardeviation,skewness,kurtosis,the se are the first order features extracted.

            2. Second order statistical fetures

               The second order texture features are derived from GLCM(Gray Level Cooccurrence Matrix).It is a statistical approach mainly used to characterize the texture based on the spatial relationship between the pixels.

               First and Second order features used in this work are:

               Mean: It is defined as the average value of the samples. Its value is calculated by taking sum of all the observed results from the sample divided by the total number of measures and it is given by the below formula

               Where n represents the sample size and x represents the observed value.

               Mode: It is defined as the number with the highest frequency, which means, for a particular set of data, the number which occurs maximum number times will be the mode of that data set.

               Variance: It is the measure of the spread between the numbers in a dataset. It is calculated by subtracting the number by mean value, squaring their differences and dividing the sum of the squares by the total number of events in the set. It is given by the below formula

               Standard deviation: It is obtained by raising the power ½ to the variance value.

               Energy:By summing all the squared elements in the GLCM(Gray Level Co-occuurrence Matrix) gives the energy value of that matrix. And it is represented by the below mathematical form

               Correlation: It measures how one pixel correlated to its neighbor. Its range is from -1 to 1.

               1 indicates perfectly positively correlated.

               -1 indicates perfectly negatively correlated. NaN indicates for a constant image.

               K=E(x-µ)4 /4

               Whereµ represents the mean value of x, represents the standard deviation of x, and E(t) represents the predictable value of the quantity t.

         2. Hus Moment Invariant Features

            Moments are scalar measures used to describe a function and to capture its important features. These features are nothing but shape descriptors. From this we extracted seven invariant moments. The major challenge for the pattern analysis is the recognition of objects (pattern) irrespective of the position, size and orientation The seven moment invariant features are explained below:

            For the digitally sampled 2D RR image, gray function is given by f(x,y),(x,y=0,1..R-1)and its central moments rpqor

            µpqare computed as x1=m10/m10 and y1=m01/m00..And its mathematical formula is given by/p>

            µpq=(x-x1)p.(y-y1)q f(x,y).

            Central moment expression, when the scaling normalization has been applied as illustrated below

            pq/ 00

            pq= µ µ where = [(p+q)/2]+1 rpq=(x) p.(y)q f(x,y). p=0, 1, 2

            when the moment f(x,y) translated by an amount (a,b),it is defined as

            Its mathematical form is given by

            Contrast:It yields a degree of the intensity contrast between a pixel and its neighbor for the whole image.And its range is given by the beolw formula

            Range = [0 (size(GLCM,1)-1)^2]

            Its value is 0 for a constant image.Its mathematical form is given by

            Homogeinity:Ityields a value that measures the closeness of the distribution of elements in the GLCM to the GLCM diagonal.

            Range varies from 0 to 1.

            For diagonal GLCM its value is 1.

            Skewness: It yields a value that measures the unevenness of the data around the sample mean. If its value is negative then the data has moved more to the left than to the right. If it is positive then the data is moved more to the right than to the left. Its value is zero when it is normally distributed.

            Mathematical form of the skewness is given by the below formula

            Whereµ represents the mean value of x, represents the standard deviation of x, and E(t) represents the predictable value of the quantity t.

            Kurtosis: It yields a degree of spreading of the outlier prone. Its value for normal distribution is 3. Distributions thatare more outlier-prone than the normal distribution have the kurtosis value greater than 3.Its mathematical form is given by the below formula

            µpq=(x+a) p. (y+b) q f(x, y).

            Hus defined seven moment values in terms of central moments are given as follows:

            R1 = (20 +02),

            R2 = (20 02)2+ 4211,

            R3 = (30 312)2+ (321 03)2, R4 = (30 + 12)2+ (21 + 03)2,

            R5 = [(30 312)(30 + 12)[( 30 + 12)2 3(21 +

            03)2]

            + (321 03)(21 + 03)[3(30 + 12)2 (21 + 03)2], R6 = [(20 02)[( 30 + 12)2 (21 + 03)2

            + 411 (30 + 12)(21 + 03)]

            R7 = [(321 03)(30 + 12)[( 30 + 12)2 3(21 +

            03)2]

            (30 + 312)(21 + 03)[3(30 + 12)2 (21 + 03)2].

      3. Neural Network Implementation Results for the Extracted Features

         Implementation results obtained for all the29 features by assigning different number of hidden neurons in each hidden layer and the best performance result obtained are as illustrated in Table2 and in Figure6.

         Table2: Performance measures for all 29 features.

         Figure6:Best performance result for two hidden layer

      4. Feature Selection

      Features extracted from the textures of the images are all not able to give the good discrimination results and the high dimentionality of the data requires large amount of time as well as memory.In order to improve the classification accracy feature selection method has been proposed.In this all the extracted features are reduced to a minimum set features to improve the accurcy of discrimination.Two feature selection methods have been proposed in this work. They are as follow:

      1. Experimentation with Wrapper Method

         Selecting the best subset of features from the available features, wrapper method has been proposed to achieve the finest possible performance by choosing a specific learning algorithm for a precise training data set. We explore the relation between best feature subset selection and relevance as shown in Figure7.

         Set No	Features
         1	Median, Energy[0:1],F2(2nd order moment),correlation [-1:1],contrast[-1:0]
         2	F1(1st order moment),contrast [-1:-1],contrast[0:1],kurtosis
         3	Energy[-1:-1],homogeneity[- 1:1],Energy[-1:0],F4(4thorder moment0,correlation[0:1]
         4	F3(3rd order moment),variance, energy[-1:1],F6(6thorder moment)
         5	Standard deviation, correlation [-1:1],skewness, correlation[- 1:0],F1(1st order moment)
         6	F5(5thordermoment),mean, homogeneity[0:1],homogeneity [-1:-1],F7(7th order moment)

         Table3: Feature sets for wrapper approach

      2. Experimentation with Fischer Ratio Method

         It is generally used method for feature selection. It selects the features based on the scores under the fischer criterion. Higher the score for the featuremeans it gives more optimistic result than the feature which is having lower score.In particular,forthe selected k features the input data matrix X

         Rd*reduces to Z k*n. Results obtained are as shown in Table3.

         Then feature score is computed as

         1

         2

         FDR= (µ1- µ2)2/ 2+ 2

         The intersection of both the approaches gives the final feature set during the experimentation. It has been noticed that both wrapper method and fischer ratio method have given the same features which could be used in the classification.

         Figure7: Wrapper method approach for feature selection

         In this approach all the 29 features are divided into 6 sets they are shown in Table1.In order to find out the best possible feature set, trial and error method has been used. Mean square error have been calculated for each set, higher the mean square error less efficient in discrimination. Out of the 6 sets, last 3 sets shown in Table1,yields the best performance and again permutation and combination among these sets to obtain the best 10 features which are more efficient in obtaining the good discrimination as shown in Table3.

      1. Classification

      Feed forward neural network is a type of artificial neural network is a technique of artificial intelligence that has the skill to learn from involvements, cultivating its performance by adapting to the changes in the environment. The main advantages of this network are: the possibility of efficient manipulation of large amounts of data and its ability to generalize results.

      Figure8: Basic structure of Neural Network

      In this type of neural network the information moves in only one direction that is from input node to output node.

      The data has been divided into training set and testing set. Training set is larger than that of the testing set in order to train up the network properly.

   5. STATISTICAL MEASURES

      In the proposed method data has been divided into two sets, in such a way that 60% of the data have been used for the training and remaining 40% of the data for testing. Training set should be larger than that of the testing set in order to train up the network properly. From the result obtained following parameters were calculated.

      TP It stands for true positive: The number of subjects with Osteoporosis that are correctly identified.

      FP It stands for false positive: The number of Control subjects which are incorrectly identified.

      TN It stands for true negative: The number of Control subjects correctly identified.

      FN It stands for false negative: The number of subjects with Osteoporosis that are incorrectly identified.

      Sn – Sensitivity: defined as Sn = TP/TP+FN Sp – Specificity: defined as Sp = TN/FP+TN

   6. VALIDATION PROCEDURE

      Figure9:Best performance result for two hidden layers

      C. Classification Results Test Results

      Validity of an assessment is the degree to which it measures what it is supposed to measure. Test validation is an assessment to measure the degree of correlation between the test and a criterion. In the proposed work 58 images have been used for validation purpose. In that 32 images are used to train the network properly and remaining 26 have been used for testing the network. The gold standard images used for the validation have been provided by the University Dorleans.

   7. EXPERIMENTAL RESULTS

      Blind Results

      Tble5: Result obtained for test data

      1. Feature Extraction

         The final feature set has been selected through rigorous feature selection process. Using wrapper and fischer ratio method the final feature set of 10 most suitable features have been considered for final classification step.

      2. Neural Network Implementation Results

      Implementation results obtained for the best 10 features by assigning different number of hidden neurons in each hidden layer and the best performance result obtained as shown in Table4 and in Figure9.

      Table4: Mean square error for all the 10 features with different number of hidden layers and neurons.

      Table6: Result obtained for Blind Data

   8. CONCLUSION

In this work analysis and the characterization of the textured images for the discrimination of the healthy subject from that of the pathologic subject has done by using the statistical processing of the X-ray images, Extraction of first order features is not sufficient in obtaining the good discrimination. Therefore in combination with the first order features, Second order features derived from GLCM and Hus invariant moments were used for the better discrimination. They helped us in classifying the healthy subjects from that of patient with osteoporotic fracture.By using feed forward neural network classifier achieved 97% sensitivity; it shows the proposed algorithm is efficient in discriminating the healthy subject from that of subject with bone pathology.

REFERENCES

1. Taleb-ahmed, Dubois, Duquenoy, Analysis methods of CT-scan images for the characterization of the bone texture.Laboratories dAnalyse des Systems du Littoral Cote dOpale.

2. Entropy Ahmed Salmi EL Boumnini El Hassani, Mohammed El Hassouni2, RachidJennane3, Mohammed Rziza1, Eric Lespessailles4,Texture Analysis for Trabecular Bone X-Ray Images Using Anisotropic Morlet Wavelet and Renyi .

3. Andrzej Materka,Texture Analysis Of X-Ray Images For Detection Of Changes In Bone Mass And Structure. Institute of Electronics, Technical University of ód, Poland,

4. Florian Yger. Challenge IEEE-ISBI/TCB: Application of Covariance matrices and wavelet marginal. October 13, 2014. Tokyo Institute of Technology, October 13, 2014.

5. Florian Yger. Challenge IEEE-ISBI/TCB: Application of Covariance matrices and wavelet marginal. October 13, 2014. Tokyo Institute of Technology, October 13, 2014.

6. J.C.Buckland-Wright, J.A.Lynch, J. Rymer, I. Fogelman, "Fractal signature analysis of macroradiographs measures trabecular organization in lumbar vertebrae of postmenopausal women", Calcif. Tissue Int., Vol 54, pp106-112, 2003.

7. J.A.Lynch, D.J. Hawkes, J.C. Buckland-Wright, "Analysis of texture in macro radiographic of osteoarthritic knees using the fractal signature", Phys. Med. Bio, Vol 36, pp 709-722, 2009.

8. S.Naor, R. Hartley, D. Avnir, "Multiple resolution texture analysis andclassification", I.E.E.E Trans. Pattern. Anal. Mach. Intell, Vol 6, pp 661-674, 2001.

9. C.Tricot, "Courbe et dimension fractals", Springer Verlag, Ed. Science et Culture, 1993.

1. M Tuceryan and A K Jain, Texture analysis, in the handbook of pattern recognition and computer vision (2nd edition)," 1998.

2. C.L. Benhamou, S. Poupon, E. Lespessailles, S. Loiseau, R. Jennane, V. Siroux, W. Ohley, and L. Pothuaud, Fractal analysis of radiographic trabecular bone texture and bone mineral density: two complementary parameters related to os- teoporotic fractures.," J Bone Miner Res, vol. 16, no. 4.

3. R. Jennane, W. J. Ohley, S. Majumdar, and G. Lemineur,

   \Fractal analysis of bone x-ray tomographic microscopy projections.,"

4. Olmez T, Dokur Z 2003 Classification of heart sounds using an artificial neural network. Recogn.Lett.24: 617629