=Paper=
{{Paper
|id=None
|storemode=property
|title= A Dynamic Programming Approach for Heart Segmentation in Chest Radiographs
|pdfUrl=https://ceur-ws.org/Vol-841/submission_29.pdf
|volume=Vol-841
|dblpUrl=https://dblp.org/rec/conf/maics/RahejaKF12
}}
== A Dynamic Programming Approach for Heart Segmentation in Chest Radiographs==
A Dynamic Programming Approach for Heart Segmentation
in Chest Radiographs
Aarti Raheja, Jason Knapp, Chunsheng Fang
{araheja,jknapp,vfang}@riverainmedical.com
Riverain Technologies
3020 South Tech Boulevard
Miamisburg, OH 45342-4860
Abstract hybrid voting scheme was used to combine the three
Chest radiographs are the most routinely acquired exams, methods. Shape models, such as the ASM, have the
which makes their use for diagnosis cost effective. In this drawback that their fitting routine can get caught in local
paper we present a dynamic programming approach for optima [5]. This effect can become quite pronounced when
automated heart segmentation on posterior-anterior (PA) applied to images that differ significantly from those used
chest radiographs. The goal of the proposed algorithm is to to build the model. This point is particularly important in
provide an accurate and reproducible method for heart our application as abnormal hearts are precisely what we’re
segmentation, which can then be used to detect certain trying to detect. For this reason, we opted for a different
cardiac abnormalities. Our method has several advantages
approach.
over previous methods, and provides superior performance
to previously published results. One important use of heart segmentation is the
measurement of the cardiothoracic ratio. The CTR is an
important measurement that can imply cardiomegaly
Introduction (abnormally large heart) [1]. The CTR is defined as the
maximum transverse diameter of the heart, divided by the
Heart segmentation in chest radiographs is a maximum internal diameter of the thoracic cage [6]
challenging task. One major difficulty in segmenting the Research in to methods for automatic CTR extraction has a
heart is the low contrast found in the mediastinum and long history [6]. Later in the paper we show how the CTR
diaphragmatic regions. These areas are difficult to can be used for assessing the quality of a heart
visualize even by radiologists. Other aspects that make the segmentation. Although the CTR can be computed without
problem challenging include: the significant variation in segmenting the heart, segmentation is useful as it can help
heart size across patients, the presence of disease in the radiologists validate the result. Figure 1 illustrates the idea.
lungs, and poor breadth holds by patients (leading to lower We use an algorithm based on dynamic programming
contrast on the heart boundary). Despite the challenges, (DP) to segment the heart. DP, an important algorithm in
development of an automated method for heart Artificial Intelligence [7], is used in applications such as
segmentation could provide significant clinical value [1]. finding the shortest path within a graph. DP decomposes a
Several methods have been proposed [1][2][3] for complicated problem into simpler sub problems; and, based
segmenting the heart. Nakamori et al [1] discuss a method on Bellman’s “Principle of Optimality”, the optimal
to segment the heart by detecting points along the heart solution to the original problem can be obtained by
boundary, which are then fitted using a Fourier shape combining the solutions to each sub problem.
model. This method was used in [2] to automatically In the proposed algorithm we formulate the DP sub
compute the cardiothoracic ratio (CTR) in 400 problem in an innovative way. The cost matrix is generated
radiographs. Out of the 400 radiographs, 20% required using image information assigning minimum cost to the
manual intervention. It was also shown in [3] that the heart pixels having heart edge characteristics. The cost matrix is
boundaries outlined by four experienced radiologists had a generated in the polar domain since the heart shape is
high degree of variability, which is an important result mostly circular. By using this method we allow the shape
when considering how to assess automatic methods. to vary in regions where enough information is present, but
Van Ginneken et al [4] discuss several approaches to force the shape to be circular in regions of uncertainty.
heart segmentation: active appearance model (AAM), In the next sections we describe our algorithm based
active shape model (ASM) and pixel classification. The on dynamic programming in detail, followed by a
individual methods performed comparably well, though presentation of extensive experimental results and a
significantly better performance was obtained when a conclusion.
method as discussed in [10]. The average of these two
locations, as shown in Figure 3, is used as the end row
value to define an approximate bounding box around the
heart region.
The top row of the bounding box, as shown in Figure 3,
is selected as the location where the heart and the left lung
first meet. The bounding box column locations, as shown
in Figure 4, are the locations along each lung mask that are
at a maximum distance from the central column.
Figure 1: Chest Radiograph with heart outline showing the
maximum internal diameter of the thorax (ID) and maximum
transverse diameter of the heart that is the sum of maximum
right heart width (MR) and maximum left heart width (ML).
Materials and Methods
We used the 247 chest radiographs from the JSRT
database to test the method. The JSRT database is
available publicly and consists of screen-film images
digitized with a 0.175mm pixel size and 2048×2048 image
size [8]. The heart annotations for this dataset [9] are Figure 3: Landmark points computed using curvature
available and were used to evaluate our method. information on the lung masks
In Figure 2 a flowchart of the method is shown.
This bounding box is used to define a center and a radius
around the approximate heart region. The center is selected
as the midpoint of the bounding box and the radius is
selected as half of the distance between the end column
locations.
Figure 4: The heart region determined using the lung mask
Figure 2: Algorithm Flowchart Polar Transform
Region of Interest around the Heart The border of the heart is roughly circular. For this reason
we apply a polar transform defined in equations (1)-(3) to
the approximate heart region.
We first obtain the ribcage mask and the segmented lung
masks from the chest radiographs. This is done using a
( ) ( ) (1)
method developed by Riverain Technologies. The lung
masks are then used to detect locations where the air,
heart, and diaphragm intersect as shown in Figure 3. These √ (2)
locations are computed based on a curvature detection
( ) (3)
enforcing a smoother result. The total cumulative cost
where ( ) is the image in the Cartesian coordinate matrix is defined as follows:
system and ( ) is the image in the polar coordinate
system. ( ) ( ) (5)
The polar transform is applied to the image as shown in
Figure 5(a) using the center and radius as defined in the ( )
previous section. To ensure all of the heart is included, the ( ) ( ) ( )
radius is multiplied by a factor α. In this paper we selected
a α value of 1.5. The polar domain image as shown in (6)
Figure 5(b) is used to compute a cost matrix for the
purpose of dynamic programming. where T represents the transition cost. The value “s” is the
offset between pixels when going from one column to the
next. The value of this offset is not allowed to be larger
than a specified value, “k”, depending upon the desired
path smoothness. The value of k for our experiments was
set to 3 pixels.
Pixels outside the lung mask, or those having cost
values above a maximum acceptable threshold, are set to
the maximum cost value as shown in Figure 6. This causes
a straight line to be the optimum path for these regions
Figure 5.(a): Cartesian system image with center and (circular arc in Cartesian domain).
radius marked for conversion in the polar domain
Figure 6.(a): Original Cost Image
Figure 5.(b): Polar co-ordinate image expressed in terms
of radial distance.
Figure 6.(b): Cost Image with non-air pixels suppressed
Dynamic Programming
In dynamic programming, the most important part is Figure 6.(c): Cost image with pixels having cost values above
constructing the cost matrix. Each pixel in the cost matrix a maximum acceptable threshold value set to the maximum cost
is assigned a local cost, where we use low cost values for value
pixels that have characteristics typical of the heart
boundary. The local cost is defined as a linear combination Once we obtain the optimal DP solution path, the heart
of individual cost images: segmentation is obtained by transforming the path to the
Cartesian domain, as shown in Figure 7.
(4)
where is the cost based on the gradient magnitude, (a)
is the weight assigned to , is the cost
based on a smoothed gray scale image and is the
weight assigned to . The gradient is calculated by
computing the derivative along each column (derivative in
the radial direction). The gray scale cost term is defined by
first computing a nominal value for the heart-lung border.
This is done for each column within a smoothed image.
These nominal values are then used to measure each
pixel’s deviation from the expected border value. Each (b) (c)
local cost term is scaled to the unit interval prior to Figure 7: (a) The optimum path obtained from dynamic
combining. programing solution is converted into (b) Cartesian co-ordinate
Given the local cost matrix, the next step is to compute system to obtain the heart segmentation with (c) some post
the cumulative cost. The cumulative cost accounts for processing.
both the local and transitional costs. The transitional term
weights the cost of going from one pixel to the next. The Some morphological post processing is applied to make
transitional cost we use increases with pixel distance, thus the heart shape smooth and convex, see Figure 7 for an
example.
Experiments Figure 10. The reason this can occur is that the source of
low overlap is generally from the mediastinal and sub
We carried out two experiments to validate the proposed diaphragmatic regions, which do not influence the
method. First, the algorithm output is compared to the transverse diameter of the heart.
manual outlines to evaluate the accuracy of the heart Figure 11 shows the only case with a low overlap score
segmentation. In a second experiment, we compared CTR that was not due to the mediastinal or sub diaphragmatic
values extracted from the algorithm against those extracted regions. The difficulty here is the fusion of the left lung
from manual outlines. The specific aim of this experiment and colon. This leads to an inaccurate estimate of the left-
was to evaluate if a reliable CTR estimate can be obtained lower landmark intersection location, which results in
even with a low overlap score. significant under segmentation. Fortunately, such an
The overlap score, Ω, between the manually outlined occurrence is rare and is left as an area for future
heart boundary and the output of our method is defined in improvement.
equation (7).
20
CTR Relative Percentage Difference
18
(7)
16
Example 1
14
where TP is the true positive area, FP is the false positive
12
area, and FN is the false negative area.
10
Figure 8 illustrates a summary of the overlap scores
obtained by our method. 8
6
Example 2
Overlap Score Distribution on 247 images 4
30 2
0
25 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Overlap Score
Number of Images
20
Figure 9: Scatterplot comparing the overlap score with the
15 relative CTR measure
10
5
0
0 0.2 0.4 0.6 0.8 1
Overlap Score
Figure 8: Histogram of overlap scores on 247 chest X ray
images. Most of the overlap scores concentrate around 0.87,
indicating the high accuracy of our method. Figure 10: Example 1 with Figure 11: Example 2
the lowest overlap score of having low overlap score,
The CTR values are computed by detecting the internal 0.63 but a more accurate CTR
diameter (ID) of the thorax and the transverse diameter of value
the heart (TD = MR+ML, Figure 1).
Some typical output segmentations are presented in
Figure 12. As can be seen, our proposed method captures
(8) the actual heart contour fairly accurately in most of the
cases.
The ID value was derived from the ribcage mask. The
TD values were computed using the heart mask derived Discussion
from the algorithm output and the manual outlines.
A relative difference between the CTR values was An average overlap score of 0.867 ± 0.046 was obtained
computed using the above TD and ID values. Figure 9 from the 247 JSRT images. We find that our method
shows a scatterplot comparing the overlap score with the produces outputs that are close to the human observer,
relative CTR measure. From this plot we can deduce that a while comparing favorably to the other methods discussed
good CTR estimate can be obtained even with a low in the survey paper [4]. The overlap scores in Table 1 are
overlap score. An example of such a case is shown in for the three hybrid methods discussed in [4]. These hybrid
methods make use of multiple methods making them References
computationally intensive. In addition, these methods are
supervised approaches whose outputs might not extend to [1] N. Nakamori, K. Doi, V. Sabeti, and H. MacMohan,
more atypical cases. By comparison, our method is far less “Image feature analysis and computer-aided diagnosis in
complex and has the advantage of making very few digital radiography: Automated analysis of sizes of heart
assumptions about the shape of the heart. and lung in chest images,” Med Phys., vol.17: pp. 342-350,
1990.
63 [2]N. Nakamori, K. Doi, H. MacMohan, Y.Sasaki, and S.
3 Montner, “Effect of heart-size parameters computed from
digital chest radiographs on detection of cardiomegaly:
Potential usefulness for computer-aided diagnosis,”
Investigat. Radiol., vol. 26: pp. 546-550, 1991
[3] R. Kruger, J. Townes, D. Hall, S. Dwyer, S.
65 80 93
Lodwick. “Automated radiographic diagnosis via feature
extraction and classification of cardiac size and shape
descriptors,” IEEE transaction on Biomedical Engineering.,
vol. BME-19:pp. 174-186, 1972
[4] B. Van Ginneken, M. Stegmann, M. Loog.
“Segmentation of anatomical structures in chest
radiographs using supervised methods: a comparative study
on a public database”, 2004.
69
78
[5] T. F. Cootes, C. J. Taylor, D. Cooper, and J. Graham.
“Active shape models – their training and application,”
Computer Vision and Image Understanding., vol. 61(1):pp.
38–59, 1995.
[6] H. Becker, W. Nettleton, P. Meyers, J. Sweeney, Jr CM
Nice. “Digital computer determination of a medical
Figure 12: Heart Segmentation. Blue represents user annotation
and red represents output of the current method. diagnostic index directly from chest X-ray images,” IEEE
Transaction on Biomedical Engineering., vol. BME-11:pp.
67-72, 1964.
Conclusion [7] Stuart Russell, Artificial Intelligence: A Modern
We presented an algorithm for segmenting the heart region Approach .
using dynamic programming. The proposed algorithm [8] J. Shiraishi, S. Katsuragawa, J. Ikezoe, T. Matsumoto,
provided an accurate and reproducible method for heart T. Kobayashi, K. Komats, M. Matsui, H. Fujita, Y. Kodera,
segmentation. The presented method makes few and K. Doi. “Development of a digital image database for
assumptions about the heart shape, has a simple chest radiographs with and without a lung nodule: Receiver
implementation, and provides superior performance to operating characteristic analysis of radiologists' detection
previously published results. of pulmonary nodules”. AJR., vol. 174:pp. 71-74, 2000.
Future work will involve the collection of more data, [9] Image Sciences Institute Research Databases.
which is needed for further evaluation and the http://www.isi.uu.nl/Research/Databases/.
development of strategies for handling outlier cases. Also, [10] S. Muhammad, M Asif and M. R. Asin. “A new
additional image features for improving the local cost term approach to corner detection,” Computer imaging and
will be explored. vision, vol. 32:pp. 528-533, 2006.
Table 1: Overlap score results compared to a human observer and various methods discussed in [4].
Heart µ±σ min Q1 median Q3 max
Human Observer 0.878±0.054 0.571 0.843 0.888 0.916 0.965
Dynamic Programming 0.867±0.0460 0.636 0.846 0.875 0.898 0.944
Hybrid Voting 0.860±0.056 0.651 0.833 0.870 0.900 0.959
Hybrid ASM/PC 0.836±0.082 0.430 0.804 0.855 0.889 0.948
Hybrid AAM/PC 0.827±0.084 0.499 0.791 0.846 0.888 0.957