Automated Vertebral Segmentation from CT Images for Computation of Lumbolumbar Angle

Abstract

Evaluation of lumbolumbar (LL) angle of the vertebral column from medical images is a common measure in case of patients suffering from lower back pain (LBP). There are five lumbar vertebrae (L1–L5). The angle formed between the slope of L1 and slope of L5 is the LL angle. This paper presents a computational technique for automated measurement of the LL angle from computerized tomography (CT) images using digital image processing techniques. The computation of LL angle comprises of four major steps: (i) image denoising (ii) segmentation of vertebral column (iii) morphological operations, and (iv) statistical estimation. To the best of our knowledge, the proposed technique for automated LL angle from CT images is the first of its kind. The segmentation based on s-t cut gives a promising result. We have performed our experiments over 30 CT images and have achieved satisfactory results in automated detection with an accuracy of 80 %.

1 Introduction

Lumbar vertebrae of human vertebral column carry major weight of the body. Lower back pain (LBP) [1] is a terminology that is commonly heard in orthopedics. It affects about 40 % [1] of people at some stage of their lives. Therefore, lumbolumbar (LL) angle [2] is a common measure that is followed to realize the abnormality of lumbar region. It also helps in monitoring pre- and post-treatment condition of the patient. CT [3] and magnetic resonance imaging (MRI) [3] are recommended by doctors for diagnosis and for having a better realization about the disease.

Computer-aided diagnosis (CAD) makes the treatment procedure much fast, easy, and cost-effective. Many research works are intended in developing a reliable and CAD for vertebral abnormalities using CT and MR images. The edges of the vertebrae are quite prominent in CT and MR images, and they have a much better contrast. But still there lies a challenge for CT and MRI in totally automating the diagnosis technique without any user intervention. The pixel intensities of a particular region may be dissimilar; again pixel intensities of different regions may be quite similar. To make the computer-aided computation reliable, it is very essential to have a robust segmentation technique.

There are a good number of algorithms that have been proposed for vertebral segmentation in recent years. In [4], the authors proposed a method that automatically segments the spinal cord and canal from 3D CT images. For segmenting the spinal canal, an extended region-growing technique is suggested whereas for spinal cord segmentation active contours are applied. The algorithm needs one seed point in the central location of the spine where the seed point selection is for determination of starting slice and for providing the positional hint for segmentation. A vertebral segmentation technique from 2D CT images was proposed by Graf et al. [5]. The method possesses five steps: noise removal, region extraction and weighting, candidate generation, candidate selection, and dynamic refinement algorithm. In [6], the authors presented a research work on segmentation of vertebra from MR images based on multiple feature boundary classification and mesh inflation. But this segmentation technique was semi-automated, and a start point was needed for vertebral initialization. In recent years, graph cut [7] segmentation has gained popularity for its reliable results. Many research works can be found based on graph cut segmentation for evaluation of different diseases. In 2008, a graph cut approach for fully automated brain tumor segmentation in 3-D MRI was proposed by Wels et al. [8]. Cui [9] proposed a work on fully automatic segmentation of white matter lesions from MRI. Ababneh et al. [10] proposed a work based on graph cut from knee bones MRI for arthritis research. It is a content-based system that automatically segments knee bones and does not require user interaction. Evaluation of spine deformity in orthopedics carries a good importance in medical treatment. Research works are going on for computerized diagnosis of various spine deformities from medical images. Bhole et al. [11] proposed a technique for automatic detection of lumbar vertebrae and disk structure from MR images. Ghosh et al. [12] proposed an automated lumbar vertebral segmentation from clinical CT image that helps in wedge compression fracture diagnosis. The technique of computerized diagnosis was based on a collection of image processing techniques. Egger et al. [13] proposed a local template-based s-t cut segmentation. This research efficiently extracts the vertebrae from the MRI. But, the major drawback of the technique lies in placing landmark (by user) on every vertebra for locating the individual template positions. This makes the total process time-consuming and laborious. The intensities that resemble the foreground and background are selected manually. In our research work, we have improved this template-based segmentation technique.

In this paper, we focus on CT images and have performed an automated computer-aided LL angle computation where the segmentation of the computer-aided technique is based on s-t cut. Our specific contributions are as follows:

• For segmentation, instead of considering all pixels of the image we have considered a single template whose dimension is similar to the image. The template distributes sampled nodes whose number is much less than the number of pixels of the image.

• The selection of source for s-t cut is done adaptively. The source intensity is computed from the histogram of the input image.

• The total procedure of LL angle evaluation is performed without any user intervention, and the research work aims in developing a totally automated CAD with high reliability.

The remaining sections of the paper are organized as follows. Section 2 describes the image processing techniques used in our proposed scheme of LL evaluation. Section 3 includes results and experimental data. Conclusion and future scope of research work are discussed in Sect. 4.

2 Proposed Method

2.1 Image Denoising

The medical CT images are sensitive to noise, which results in false intensity value. The contamination of image by noise may disturb the post-processing steps. So, a good denoising technique is necessary to make image noise free. State-of-the-art technique for image denoising is bilateral filtering [14]. It smoothes image as well as preserve edges. The performance of the filter is based on two kernels, the domain and the range kernel. The domain kernel is meant for noise removal in the homogeneous region, and the range kernel is designed for cleaning noise in discontinuous regions of the image. A mask of size 3 × 3 was selected for the denoising technique. A standard deviation of 3 and 0.1 was considered for the domain and range kernel of this preprocessing image enhancement. The segmentation is performed on this enhanced image. Denoised images after bilateral filtering are shown in Fig. 4a, c, and e.

2.2 Vertebral Segmentation

2.2.1 Graph Cut

Segmentation is partitioning the region of interest. We present a segmentation based on graph cut [7] that separates the vertebral column from the lateral view CT image. A template is considered that possess 7,200 points, which are sampled on the whole image. The CT image is rectangular in shape, and normally the height of the image is almost twice the width. It is observed that considering the row number as 120 and column number as 60 for the template gives a good segmentation result. This helps in maintaining an equidistant position between the nodes of the graph. The consideration of template for segmentation instead of all pixels of image allows reduction of computational complexity from N ² to T ² where N × N represents the dimension of the image and T ² is the number of nodes in the template. Here, we have considered T ² = 7,200, which is very small in comparison with number of pixels in an image.

The sampled points are the nodes n ∊ V of the graph. The graph G(V, E) is composed of set of nodes V and set of directed edges E. Every n of the graph possess grayscale intensity in the range [0 1]. The source s and sink t are the two virtual nodes of G. The basic principle of the proposed segmentation depends on intensities of s and t. There are E that connect the s and n of the graph and there are E between n and t. Intensity of s and t resembles the intensity of foreground (vertebral column) and background, respectively. Every E of the graph is assigned a positive weight w. n that has intensity value more similar to s will have more weight for the E that connects between s and n where as in that case n will have less weight for the E that connects the same n to the t. The reverse principle is applicable when the intensity value of n is more similar to the intensity of t. In our segmentation technique, the mathematical expression of weights (w) for s and t are as follows:

$$ w_{s} = \exp \left( { - \frac{{{\text{d}}(s,n)}}{h}} \right) $$

(1)

$$ w_{t} = \exp \left( { - \frac{{{\text{d}}(t,n)}}{h}} \right) $$

(2)

w _s and w _t are the weights for s and t, d(s, n) is the intensity difference between s and n, and d(t, n) is the intensity difference between t and n. h is the control parameter of the weight.

The s-t cut is grouping the n of the graph (template) into two subsets S and T such that n similar to s lies in S and the n similar to t lies in T. Figure 1 shows a simple example of the theory.

Fig. 1

a Template for s-t cut segmentation, where the color indicates the gray color intensity that the nodes possess. b Nodes having intensity similar to s have more weighted edge for s and vice versa. c Nodes of the template selected in the group of foreground are indicated in white and nodes selected in the group of background are indicated in black. d Figure represents only the border nodes of the foreground-segmented regions. The border nodes are marked white

2.2.2 Automated Selection of Source

The selection of s and t by manual landmarks is a tedious process. Even assigning predefined values to s and t may lead to erroneous results because the intensity of regions of CT varies from image to image. We have observed that selection of s and t values from the histogram make the procedure much adaptive. The wide range of intensity [0 1] of CT image can be grouped into five major subgroups. To have an idea, a CT of vertebral column and its histogram are shown in Fig. 2, where the five modes of the histogram represent five major subgroups of the image. The first three modes from left to right represent the background, and the remaining modes on the right are the representation of foreground (vertebral column).

Fig. 2

a CT image of vertebral column and b Corresponding histogram

The CT images may be of low, medium, and high contrast. For the automated s-t cut segmentation, it was observed that sink value of 0.15 for all types of contrast images generates a reliable segmentation. Whereas the formula for selection of source intensity can be expressed as follows:

$$ {\text{source}}_{\text{m}} = I_{\max} + (I_{\max} - t) \times 2.3 $$

(3)

$$ {\text{source}}_{\text{l}} = I_{\max} + (I_{\max} - t) \times 2 $$

(4)

$$ {\text{source}}_{\text{h}} = I_{\max} + (I_{\max} - t) \times 4 $$

(5)

where I _max is the maxima of the histogram above intensity 0.1 and source_m, source_l, and source_h are the estimated source intensity for medium-, low-, and high-contrast CT images, respectively.

2.3 Extraction of Vertebral Body and End Plates of Vertebrae

The segmentation of vertebral column was followed by morphological operations. Initially, the holes of the extracted regions were filled. And to identify the regions of the template separately, each of them were labeled uniquely. The computation of LL angle only requires the vertebral body end plate. So in the next step of the CAD, the interest was to keep the vertebral body and discard other parts of the vertebral column. The features, average intensity, and location of the regions were used to recognize the vertebrae. The vertebrae of the segmented vertebral column have darker intensity and lie to the left than other segmented foreground regions. The average intensities (of regions) those were more similar to the minimum average intensity of all regions were considered as vertebrae and regions with average intensity more near to maximum average intensity of all regions were discarded. Images of Fig. 5b and c give clear idea of the section. After selection of the vertebrae from the segmented regions, morphological operations such as removing isolated nodes, removing small regions, extraction of region boundary, and removal of vertical lines were performed to extract end plates of the vertebrae. All these processing were done with the nodes of the template, where the template was superimposed on the image.

2.4 Lumbolumbar Angle Computation

In manual measure of LL angle, pencil and scale are used to draw lines on the printed medical images. A line is drawn through the end plate corner points of the L1 vertebra, and similarly another line is drawn through the end plate corners of the L5 vertebra. The angle formed by the intersection of these two lines is the LL angle [2].

In our computer-aided process, the extracted end plates were relabeled so that end plates of every vertebra could be identified. As there are five lumbar vertebrae and every vertebra possess two end plates, it was observed that selection of second and tenth end plate (bottom to up direction of the template) helped us in getting the plates that represented the slope of extreme vertebrae (L1, L5) of the lumbar curvature. At this point, an intensity value of 1 was assigned to the nodes of one end plate and an intensity value of 2 was assigned to the nodes of second end plate. And the other nodes of the template were assigned value 0. Next, a scanning was performed on the template to determine the end nodes of each of the finally selected end plate. The slope of line obtained by joining the end nodes of an end plate can be considered similar to line obtained by joining the corner points of the end plates. A neighborhood window of size 3 × 3 was considered for this operation. Nodes for the first end plate were selected as end nodes where a sum of intensity of the window was 2 and end nodes of the second end plate were marked where the sum resulted as 4. Figure 3 shows the end node selection conditions for end plates one and two. It can be seen from Fig. 3b that sum of values of any condition is 2. If a node satisfies any one of these condition, it is end node of first line. Figure 3c shows the conditions for considering end nodes of second line where sum of values of every condition is 4. The selected end nodes of a particular end plate were joined with each other to get a line. The angle of the line was computed with respect to the horizontal axis. A sum of angles obtained from superior and inferior end plate gave the LL angle.

Fig. 3

a Shows the 3 × 3 neighborhood. b Illustrates conditions for regarding a node as end node of first line. c Shows the conditions of end node for second line. d Shows example of identification of a node as end node when it matches any one of the above conditions

3 Results and Discussion

MATLAB 2011b was used for implementing the CAD. The images that were used for the experiment varied in resolution. As a result, the time taken to perform the image denoising varied from one CT image to other. The post-processing steps were based on fixed number of template nodes. So, the computational time for later steps remained constant. After denoising, the overall time taken to perform the template-based segmentation and computation of LL angle from nodes was 7 s. A machine specification of Intel Core2 Duo CPU, 3 GHz, 3 GB RAM, Windows XP Professional Version 2002, and Service Pack 3 was used for this work. The time required to process the whole CAD on CT was about 14–25 s.

The vertebral body segmentation results from different patient’s CT image are shown in Fig. 4a, c, and e show the denoised images, and b, d, and f show the corresponding extracted vertebrae. Due to similar average intensity of regions (explained in Sect. 2.3), there remains a chance that sometimes some other parts of the vertebral column get extracted with vertebrae. Such an example is shown in Fig. 4d. As the main objects of interest are horizontal end plate lines from the square vertebrae, it is observed that post-morphological operations can easily remove these unwanted regions.

Fig. 4

a, c, and e are the bilateral filtered CT images. b, d, and f show corresponding extracted vertebrae. The white-colored nodes denote vertebrae, and the black are the background nodes. The rest of the computation for LL angle depends on these white foreground nodes

Figure 5 shows the results of main steps of our proposed technique. The final image of the figure shows the identified slope of extreme vertebrae (L1, L5) for the lumbar curvature.

Fig. 5

Steps of our technique is as follows : a denoised CT image and b shows the automated segmented vertebral column based on s-t cut. c Shows the selected vertebral body (vertebrae) from (b). The white nodes of d show the extracted end plate lines of the vertebrae from (c), and e shows the identified slopes of the superior and inferior lumbar vertebrae. The LL angle is evaluated by computing the angle formed by intersection of these two lines

Next, we have performed an experiment to judge the segmentation accuracy. A direct comparison is done between the ground truth developed by manual landmarks and the segmented vertebra obtained from our proposed s-t cut. The comparison is shown in Fig. 6. The white-colored nodes of Fig. 6b are the ground truth of vertebra boundary, and Fig. 6c is our segmentation result on the ground truth image, Fig. 6b. Figure 6c shows that the ground truth of vertebra boundary directly matches with the boundary of the automated segmented vertebra (white).

Fig. 6

a Shows the input image and b is ground truth. White points were marked manually on the boundary of middle vertebra of the image a. c Shows the superimposed automated segmented result on the ground truth where both the boundaries efficiently match each other

The LL angles (in degree) computed from our automated CAD are shown in Table 1. We have performed the experiment of computerized measurement over low-, medium-, and high-contrast CT images. In every case, manual measurement [2] of LL angle over printed CT image was recorded. For this, lines were drawn on the end plate of the upper lumbar vertebra and on the end plate of the lower lumbar vertebra. The interior angle formed by the intersection of these two lines is the LL angle. Manual LL angle measurement may lead to some variability if the lines that are running over the end plate of the vertebra are not drawn exactly through the extreme corners of the end plate. An intra-observer CAD was recorded where two trials were taken by the user for every CT image. Considering the manually measured LL angle as standard value, a variability of ±0.1° to ±5° has been noticed from our CAD technique. In medical an error of ±5° [15] can be tolerable for angle measurement. This demonstrates the clinical acceptability of our CAD technique.

Table 1 LL angle value from CAD

The constraints such as contrast and brightness play a major role in the success of any medical CAD. If the contrast and brightness are kept constant in the CT imaging, it will further help in reducing the variability of our CAD to a much greater extent.

4 Conclusion

In this paper, we have presented an automated LL angle evaluation process from clinical CT image. Initially, the source intensity is computed from histogram, and segmentation is performed using template-based s-t cut. A proper selection of s and t gives a very reliable segmentation result. The automated CAD of LL angle from CT image reduces user intervention. This research work is a good initiative for automated clinical application. Digital measurement also eases up the processing as it can be entirely done on the system, which captures the source image. In future, we will try to reduce the variability of this computer-aided evaluation and also we will try to extend our work of automated LL angle computation to digital X-ray image cases.