Statistical properties of dynamical chaos

1.   Introduction

Cardiac Magnetic Resonance (MR) Imaging provides high spatial and temporal resolution images and has been proven to be a powerful imaging modality ^[1,2]. It can capture the anatomical and functional information of a heart within a cycle ^[3], and thus, be well suited for the clinical diagnosis and research of cardiovascular diseases (CVDs) for radiologists and physicians. Quantification of the Left Ventricle (LV) of the heart usually acts as index of cardiac health ^[4]. Therefore, segmentation of the LV from cardiac MR images is an essential prerequisite for the diagnosis. From accurate segmentation results, many clinical diagnosis indices such as wall thickness, myocardial strain, ejection fraction, and circumferential shortening of myocardial fibers can be evaluated ^[5].

Clinically, delineating the boundaries of the LV is often performed manually. However, manual segmentation often takes a long time, and is tedious, and subject to high inter-observer variability. Thus, there is a flurry of works devoted to segmenting the LV from MRI in the past decades, such as classification ^[6,7,8], watershed ^[9], graph cuts ^[10], active contour models ^{[11,12,13,14,15,16]}, level set ^{[16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31]}, active shape models ^{[32,33,34,35,36]}, and other methods ^[37,38]. Petitjean ^[39] and Peng ^[40] presented a comprehensive review respectively on this topic, interested readers can refer to ^[39,40] for details and references, and we merely present a brief survey here.

Active contour models, or the snake models, were one of the most popular tools for segmentation of the region of interesting (ROI), and found successful applications in extracting the LV from MRI. The basic principle remains the same as those presented in the seminal work on the snake model proposed by Kass et al. ^[41]. Wu et al. coined the circle shape energy to aid the snake model, which considered the LV as a circle, and extracted the endo- and epicardium of the LV by coining the shape similarity constraint ^[11]. However, the delineation merely relying on the shape constraint still results in mistakes when the discontinuities are found in the object boundaries. The discontinuities may be caused by the papillary muscle or the organs nearby such as the liver. Santarelli also used the snake model for this task ^[13], in which the Gradient Vector Flow ^[42] acts as an external force. However, there is no solution for artifacts and weak edges. Ranganath ^[14] modified the internal energy of the snake by taking into account the curvature of the snake contour to smooth the boundary of the LV. McInerney et al. coined an active surface model to extract the surface of the LV by representing the surface using finite element ^[15].

Geometric active contour, or level set, has been proven to be very powerful in a variety of segmentation tasks including the LV extraction. Li et al. coined a bi-layer level set ^[17] to extract the inner and outer boundaries of the LV, then the bi-ventricle ^[18]. Ismail et al. manually extract the LV in the first frame in a sequence of images, and supposed the intensity distribution of the myocardium in the first image is similar to the rest of the ones in the sequence ^[25]. Pluempitiwiriyawej introduced the ellipse-fitting scheme into the level set model and proposed an integrated model for cardiac MR image segmentation ^[28].

With the development of deep learning, deep learning-based methods have shown impressive performance in medical image processing tasks. For example, Zhang et al. proposed the BLU-Net ^[43] to segment LV, Wang et al. combined the DenseNet segmentation and convolution neural network(CNN) regression to quantify full LV metrics ^[44], Shen et al. proposed a pure dilated residual U-Net ^[45] to segment the femur and tibia bones. Besides, deep learning combined with traditional methods has also been proved to be useful. Ngo combined deep learning and level set ^[46] to segment the left ventricle of the heart from cardiac cine MR data automatically, in which the deep learning method is responsible for locating the LV, and the level set is responsible for refining the boundary of the LV. Avendi et al. combined deep learning algorithms with deformable models ^[47] to evaluate automatically the LV function from short-axis cardiac MRI datasets.

In this study, we propose an extension of our previous work in ^[11] based on the snake model. For a MR image sequence, it is often sampled from a cardiac cycle, and it is noticed that the shape of the LV in one image is similar to that of the LV in nearby ones in the sequence. In the aforementioned works, segmentation is carried out in an image-by-image manner, while in this work, we propose to segment all images in the sequence simultaneously and the shape of the LV in each image is supposed to be similar to that of the LV in nearby images in the sequence. We coined the novel shape similarity constraint, and it is called sequential shape similarity (SSS in short). Our approach integrates the prior knowledge from the image sequence into the snake model.

The organization of this paper is as follows. Section II presents a brief introduction of our previous work in ^[11]. In Section III, the proposed sequential shape similarity (SSS) is detailed. In Section IV, the experimental results are reported, and the conclusion marks are presented in the final Section.

2.   Background

We first made a brief review of our previous work ^[11], which include the snake model, GVC external force, Circle Shape Constraint (CSC) for endocardium segmentation, shape similarity constraint (SSC), and novel edge map for epicardium.

2.1.   GVC snakes

Active contour models or snakes, coined by Kass et al. in 1988 ^[41], are dynamic curves $c\left(s\right) = \left(x\left(s\right), y\left(s\right)\right), s\in \left[\mathrm{0, 1}\right]$ that evolve by minimizing the following energy functional,

where ${c}^{\text{'}}\left(s\right)$ and ${c}^{\text{'}\text{'}}\left(s\right)$ are the first and second derivatives of $c\left(s\right)$ with respect to arc length $\mathrm{s}$. $\alpha$ and $\beta$ are the weights controlling the smoothness and rigidity of the curve, respectively. ${E}_{ext}$ is derived from the image data and takes the smallest values at boundaries. The snake models gained great success in image processing.

However, the traditional snakes are limited by the external force based on image gradient. Therefore, the gradient vector flow (GVF) for snakes is proposed by Xu and Prince ^[42]. The GVF is effective, yet, with high computation costs. To overcome this shortcoming, we previously proposed a new external force called Gradient Vector Convolution(GVC) ^[48] based on convolution and fast Fourier transform (FFT). The GVC can be calculated as follows:

where ${K}_{con}$ is the convolution kernel, $\otimes$ is the convolution operator. The GVC model behaves similarly to the well-known GVF model on enlarging capture range and robustness to initialization. Moreover, using the FFT, the GVC can be implemented in real-time since it is expressed in convolution, and thus, the snake model can also benefit from this fact. In practice, we take ${K}_{con} = 1/{\left({r}_{h}\right)}^{n}$, where ${r}_{h} = \sqrt{{x}^{2}+{y}^{2}+h}$ $h\in {R}^{+}, n\in {R}^{+}$. ${K}_{con}$ always shows high performance in terms of extending and smoothing gradient vectors. Thereby, the GVC snake can preserve edges and dive into deep concavities.

2.2.   Circle shape and shape similarity constraints

For the segmentation of the endo- and epicardium of the LV, the difficulties originate primarily from the gaps or the weak boundaries caused by the artifacts and papillary muscles. In the seminal work in Eq (1), the internal energy, which consists of first and second derivatives, controls the continuity and smoothness of the contour, respectively. However, the first and second derivatives are merely local geometrical properties of the contour and if there were weak edges or gaps, the snake contour would fail to bridge the gaps and leak out. Moreover, if there is a local minimum caused by the external force, the snake contour would be trapped. To address these issues, the information about the overall shape of the endo- and epicardium should be incorporated into the snake energies. The shape information will guide the snake contour to delineate the target object correctly. The circle shape constraint (CSC), which roughly regards the shape of the endocardium as a circle, is adopted for the segmentation. This shape constraint tries to keep the contour being a circle during the movement. The circle shape constraint can be represented as follows,

where $R\left(s\right) = \sqrt{{(x\left(s\right)-{x}_{c})}^{2}+{(y\left(s\right)-{y}_{c})}^{2}}$, $\stackrel{-}{R} = {\int }_{0}^{1}R\left(s\right)ds$, ${x}_{c} = {\int }_{0}^{1}x\left(s\right)ds$, ${y}_{c} = {\int }_{0}^{1}y\left(s\right)ds$. $({x}_{c}, {y}_{c})$ is the center of the snake contour. The energy function in Eq (3) characterizes the deviation of the snake contour from a circle of radius $\stackrel{-}{R}$ located at $({x}_{c}, {y}_{c})$. Both radius $\stackrel{-}{R}$ and $({x}_{c}, {y}_{c})$ change along with the motion of the snake contour. If the snake contour is attracted by some nuisances such as the or papillary muscle, tit is still approximately a circle under this constraint, which is the global shape of the LV.

For segmentation of the epicardium, we proposed the Shape Similarity Constraint (SSC) that is based on the observation of shape similarity between the endocardium and the epicardium. The segmented endocardium is referred to as a prior shape constraint to extract the epicardium. The SSC constraint is formulated as,

The parameters in Eq (4) are meaningly similar to those in Eq (3), but $R$ and $\stackrel{-}{R}$ are for epicardium and $r$ and $\stackrel{-}{r}$ are for segmented endocardium.

2.3.   Modified edge map for epicardium segmentation

After the endocardium is segmented, the result is employed as the initial contour to segment the epicardium contour. However, the endocardium is merely a local minimum of the external energy derived from GVC (the red circle approximates the endocardium in the upper-right of Figure 1), and it would prevent the contour from converging towards the epicardium (the yellow circle approximates the epicardium in the upper-right of Figure 1. To address the issue, the original edge map is corrected by setting to zero the original edge map around and within the endocardium. This new edge map is used to calculate a new external force by GVC, and the endocardium is no longer the local minimum, and the new external force can directly put the snake contour forward to the epicardium.

3.   Proposed method: sequential shape similarity

Let $I$ be an MR image sequence containing $K$ images. Suppose that there are 20 samples in a cardiac cycle, and the number K is equal to 20. Generally speaking, the purpose of the LV segmentation is to track the endocardium ${R}_{k}^{endo}$ and epicardium ${R}_{k}^{epi}$ in each image ${I}_{k}$, $k = 1, \dots, K$. The aforementioned SSC can be described as follows, ${R}_{k}^{endo}$ and ${R}_{k}^{epi}$ in image ${I}_{k}$ are generally similar in shape and when segmenting ${R}_{k}^{epi}$, the segmented ${R}_{k}^{endo}$ is referred to as the prior shape constraint. For SSC, only the similarity in a single image is integrated into the snake model and some useful information among the images in a sequence is discarded. For example in Figure 2, there are gaps in the myocardium in Figure 2(b), see the location in the red circle, and it is hard to extract the myocardium in an image-by-image manner. However, in Figure 2(a), (c), the myocardium is intact, therefore, the shape of the myocardium in Figure 2(a), (c) can be employed to construct the constraint to segment the endo- and epicardium in Figure 2(b).

Based on the observation in Figure 2, we propose the sequential shape similarity (SSS in short). For an image ${I}_{k}$ (say the one in Figure 2(b)) and its neighbors ${I}_{k-1}$ (the one in Figure 2) and ${I}_{k+1}$ (the one in Figure 2(c)) in sequence $I$, the shape of the endocardium ${R}_{k}^{endo}$ is similar to ${R}_{k-1}^{endo}$ and ${R}_{k+1}^{endo}$. Thus when segmenting ${R}_{k}^{endo}$, besides the CSC in Eq (3), ${R}_{k-1}^{endo}$ and ${R}_{k+1}^{endo}$ can be also referred to as the prior shape constraint. Furthermore, the similarity range can be extended from ${R}_{k\pm 1}^{endo}$ to ${R}_{k\pm m}^{endo}$, $m\in [1, \cdots K/2]$. For ${R}_{k}^{endo}$, there are 2 m endocardium in the sequence to construct a more effective similarity constraint. For ${R}_{k}^{epi}$, besides the SSC in Eq (4), ${R}_{k-1}^{epi}$ and ${R}_{k+1}^{epi}$ can be also referred to as the prior shape constraint. The proposed constraint is called sequential shape similarity (SSS). For clarification, we neglect the upper indices (epi, endo) in the derivation of the SSS, and the SSS can be represented as follow

where $m\in [1, \cdots K/2]$. The parameters in Eq (5) are meaningly similar to those in Eq (4), but ${R}_{k+j}\left(s\right)$ and ${\stackrel{-}{R}}_{k+j}$ are for image ${I}_{k+j}$ and ${R}_{k}\left(s\right)$ and ${\stackrel{-}{R}}_{k}$ are for image ${I}_{k}$. In Eq (5), it is clear that ${R}_{k}\left(s\right)-{\stackrel{-}{R}}_{k}$ characterize the deviation of the snake contour from a circle of radius ${\stackrel{-}{R}}_{k}$ in ${I}_{k}$ and ${R}_{k+j}\left(s\right)-{\stackrel{-}{R}}_{k+j}$ measures the deviation of the snake contour from a circle with a radius ${\stackrel{-}{R}}_{k+j}$ in ${I}_{k+j}$. Minimizing ${E}_{SSS}^{k}$ makes the snake contour resembling the shape of LV in all the images in sequence $I$.

Suppose there are n discrete points on the snake contour, ${\stackrel{-}{R}}_{k}$ and ${R}_{k}$ in Eq (5) can be estimated by ${\stackrel{-}{R}}_{k} = \frac{1}{n}{\sum }_{i = 1}^{n}{R}_{i}^{k}$, ${R}_{i}^{k} = \sqrt{{\left({x}_{i}^{k}-{x}_{c}^{k}\right)}^{2}+{\left({y}_{i}^{k}-{y}_{c}^{k}\right)}^{2}}$, and ${x}_{c}^{k} = \frac{1}{n}{\sum }_{i = 1}^{n}{x}_{i}^{k}$, ${y}_{c}^{k} = \frac{1}{n}{\sum }_{i = 1}^{n}{y}_{i}^{k}$, $i = \mathrm{1, 2}, \cdots, n$, $\left({x}_{c}^{k}, {y}_{c}^{k}\right)$ is the centroid of the contour in ${I}_{k}.$ By calculus of variation, the discretized equation of Eq (5) reads,

here, the following function is employed to determine the weight ${\rho }_{j}$,

${\rho }_{j}$ is to weigh the degree of similarity of ${R}_{k}^{endo}$ and ${R}_{k+j}^{endo}$. For the image ${I}_{k+j}$, it is obvious that when $j = \pm 1$, ${I}_{k+j}$ is the neighbor of ${I}_{k}$, and there should be a larger ${\rho }_{j}$ to weigh the similarity; when $j$ increases or decreases, the ${I}_{k+j}$ is far from ${I}_{k}$ and the weight ${\rho }_{j}$ should be decreased accordingly. Figure 3 illustrates the weight for the image ${I}_{k+j}$ with varying$j$, where $a = 1, b = 2$.

The solution of Eq (5) is based on the shape similarity of the LV in the sequence, and the contours are changing during the snake deforming process. Thus, the deforming process of all of the contours should be carried on simultaneously, and all the images in the sequence are supposed to be segmented simultaneously. With the SSS and the strategies in Section II, we summary the segmentation algorithm in Algorithm 1.

4.   Experimental results

In order to verify the performance of the proposed strategy, the parameters of our model are $\alpha = 0.5, \beta = 0.5, \lambda = 0.4, \eta = 0.6, n = 2.6, h = 1, a = 1, b = 2, m = 1$. In our experiments, two datasets are employed for evaluation: MICCAI 2009 ^[43] and the dataset published by Shuo Li (Shuo Dataset in short) ^[43,44,49]. Cardiac cine MRI data of MICCAI 2009 is provided by the MICCAI 2009 cardiac MR LV segmentation challenge organizers ^[43], and there are 2900 2D short-axis cine MR images of 145 subjects in the Shuo Dataset, and interested readers can find a detailed description in ^[49].

Two metrics are employed to evaluate the segmentation performance, i.e., the Mean Absolute Distance (MAD) and the Hausdorff Distance (HD). Suppose the snake contour $S = \left\{{s}_{1}, {s}_{2}, \cdots, {s}_{p}\right\}$ has p points, and the ground truth $M = \left\{{m}_{1}, m, \cdots, {m}_{q}\right\}$ has q points, the MAD and HD are respectively defined as

where $d\left({s}_{i}, M\right) = \underset{{s}_{i}\in S, {m}_{j}\in M}{\mathit{min}}‖{s}_{i}-{m}_{j}‖$ is the distance from point ${s}_{i}$ to the closet point on contour $M$ and $d\left({m}_{j}, S\right) = \underset{{s}_{i}\in S, {m}_{j}\in M}{\mathit{min}}‖{m}_{j}-{s}_{i}‖$ is that from point ${m}_{j}$ to the closet point on contour $S$.

The MICCAI 2009 dataset is first employed for comparison, and the radial GVF method (RGVF) ^[16], level set method ^[23], the max-flow method (MFM) ^[25], and our previous work ^[11] are applied for comparison on this dataset, the detailed description of this dataset has been presented in ^[43]. The endocardia in Figure 5, especially in images #6–#10, are challenging for the snake model. However, the work in ^[11] regards the shape of the endocardium as a circle and employs the CSC as constraint, while the proposed strategy supposes that the shape of the endocardia in the consecutive image is similar and employs the SSS and CSC as constraints. The two methods both get pleasant results. But for epicardium segmentation, the artifacts and noise make some of the boundaries 'missing', and the snake moves to a false contour shown just with SSC in Figure 5(a), the results are a little far from pleasant. Figure 5(b) shows that the extra SSS effectively pushes the snake contour to the expected solution. Figure 6 shows the MAD errors corresponding to Figure 5. Overall, Integrating SSS to the GVC snake model can effectively reduce the effect of artifacts, and prevent the snake from weak edge leakage. One can learn from ^[11] the average MAD of RGVF, LSM, and MFM and our previous method for the endocardium segmentation are 9.06 pixels, 7.21 pixels, 4.79 pixels and 5.06 pixels, respectively, and those of the epicardium are 7.85 pixels, 5.70 pixels, 5.23 pixels and 5.18 pixels, respectively. Comparing them with the results shown in Figure 6, we can see the proposed strategy is able to achieve much better results. More representative segmentation results on MICCAI 2009 dataset are presented in Figure 7.

To further demonstrate the effectiveness of SSS, the Shuo dataset is also employed for test ^[49]. Besides, in this experiment, we adopt the method of component edge map ^[48] after extracting the endocardium. First, with the endocardial center as a reference, we coin two edge maps by selectively opting the gradient in direction of x and y according to the grey level characteristics of images. In the left half part, the septum links the blood pools of RV and LV while the right ventricle is brighter. Thus, ${f}_{l}^{x} = \left|\mathit{min}\left({I}_{x}, 0\right)\right|$ is utilized as an x-direction gradient to neglect the edge between the septum and the LV blood. For the right half, the lateral is surrounded by fat so we adopt ${f}_{r}^{x} = max\left({I}_{x}, 0\right)$ to erase the edge between the fat and the lung as well as that between lateral and the LV pool. Similarly, we construct the y-direction gradient as follows: ${f}_{u}^{y} = \left|min\left({I}_{y}, 0\right)\right|$, ${f}_{d}^{y} = max\left({I}_{y}, 0\right)$ where $\mathrm{u}, \mathrm{ }\mathrm{d}$ represent the upper and lower half of the center respectively. The new edge map can be obtained by setting the ${f}_{x}, {f}_{y}$ to zero around and within the endocardium. Then, this new edge map is employed to calculate a new external force by GVC. Figure 8 shows our results on a representative sequence of Shuo dataset, it can be seen that SSS combined with this new external force can get satisfactory results. The corresponding MAD indexes of endocardium and epicardium in this sequence were 1.15 pixels and 1.35 pixels respectively. We also calculated that the HD metrics of endocardium and epicardium in this sequence were 2.8210 mm and 3.5118 mm separately, while according to ^[43], the same indexes obtained from Graph cuts and Level set on some images in this dataset were 9.290 mm, 5.739 mm on endocardium and 4.200 mm, 5.484 mm on the epicardium. We can see that our method performs more excellently. Another set of experimental results of this dataset are shown in Figure 9.

5.   Conclusions

In this work, we proposed a new strategy for left ventricle segmentation in cardiac cine MRI based on GVC snake. In the proposed strategy, all the images in a sequence are segmented simultaneously and the shape of the LV in each image is supposed to be similar to that of the LV in nearby images in the sequence. We coined the novel shape similarity constraint, and it is called sequential shape similarity (SSS). For the endocardium, the SSS is combined with the GVC snake with a circle constraint. With the SSC and SSS, the epicardia are segmented successively. Experimental results on the MICCAI 2009 dataset and Shuo dataset shown that the proposed segmentation approach performs well on conquering papillary muscles, weak edges, artifacts, and noise in the LV segmentation.

Acknowledgments

This work was supported in part by the National Science Foundation Program of China (NSFC) (grant number: 61976241), and the International Science and technology cooperation plan project of Zhenjiang (grant number: GJ2021008).

Conflict of interest

The authors declare there is no conflict of interest.