=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== https://ceur-ws.org/Vol-841/submission_29.pdf
           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