Automatic Computation of Left Ventricular Volume Changes Over a Cardiac Cycle from Echocardiography Images by Nonlinear Dimensionality Reduction

Cardiac Cycle

A cardiac cycle can be divided into seven phases (Fig. 1) [1]. First, “Atrial Contraction”, that is the result of the atrial musculature contraction which gives an additional thrust to inflow blood into the ventricles. Second, “Isovolumetric Contraction”, during which ventricular volume does not change because all valves are closed. In addition, the ventricular volume is maximum called the end diastolic volume. Third, “Rapid Ejection”, the first third of the ejection period in which 70 % of ventricular emptying occurs. Fourth, “Reduced Ejection”, it is the last two third of the ejection period in which the rate of ejection decreases so that only 30 % of ventricular emptying occurs. Fifth, “Isovolumetric Relaxation”, the ventricular volume remains constant during this phase because all valves are closed. Moreover, the ventricular volume is minimum called the end systolic volume. Sixth, “Rapid Filling”, in this phase, blood flows into the ventricles rapidly from the atrium. Seventh, “Reduced Filling”, during this phase, the change in ventricular volume is very little due to low blood entering to the ventricles and less compliant of ventricles.

LV Volume Computation from 2-D Echocardiography Images

The American Society of Echocardiography (ASE) recommends several methods for LV volume computation from 2-D echocardiography images based on geometric models of the LV, dimensions, and area measurements obtained from long axis apical and short axis views [2, 3, 21]. According to the ASE, the biplane method of discs (modified Simpson’s rule) based on paired long axis apical views (two and four chamber) is the best method to compute LV volume (Fig. 2). In this process [2, 3, 21], the inner wall of the LV from the two and four chamber views is separated. Then, the main axis of LV in two views (maximal distance from the mid-mitral annulus to the LV apex) is determined and next LV contour is divided into 20 equal sections so that each section is perpendicular to the main axis (Fig. 2). Finally, LV volume is calculated according to the summation of areas of 20 elliptical discs of equal height whose centers are all in the main axis. Generally, for LV volume calculation using modified Simpson’s rule, the two apical views are used. However, for routine use, one apical view is accepted and can be applied to either long axis apical views [3]. In this case, the area of the disc is assumed to be circular and equation for LV volume calculation is as follows:

where d_i is the diameter of the ith disc LV contour and L is the length of the LV main axis from the two or four chamber views.

Data Description

The study was performed according to the Helsinki Declaration and approved by the Regional Committee for Medical Research Ethics. Two consecutive cardiac cycles of 2-D images of ten healthy volunteers and six patients (ischemic heart disease) from apical four chamber view were acquired using an ultrasound machine (General Electric Vivid 3, Horten, Norway), including the electrocardiogram (ECG) recording. Each cardiac cycle was identified by selecting the two consecutive R-wave of ECG signal synchronized with the end of diastole. The heart rates ranging are from 45 to 95 beats per minute (bpm) and the frame rates of image recording are from 40 to 70 frames per second (fps). Moreover, the curve of the LV volume changes in each cardiac cycle is calculated by a highly experienced echocardiographer. To do this, the echocardiographer first manually detects the inner wall and main axis of the LV in each image of a cardiac cycle. Then, LV volume is calculated using single plane modified Simpson’s rule method (Eq. (1)).

NLDR Algorithms

NLDR algorithms embed data sets from high to lower dimensional parameterizations while the intrinsic geometry of the data sets is respected [17–20]. Images can be thought of as points in a high dimensional space with dimension size equal to the number of pixels in the image. However, when image sets vary due to a small number of factors (lighting, deformation, etc.), these images have a natural low dimensional parameterization in high dimensional space which can be extracted by NLDR algorithms [17–20].

NLDR algorithms have been used for visualization of image sets (rendered faces, handwritten digits, head pose, and writing style) [17, 19, 22]. Moreover, these algorithms have been widely used for medical image analysis such as visualization of cardiopulmonary MR images [23] and echocardiography images [24], classification in brain MR images [25], tracking of the LV in 3-D echocardiography images and heart in 3-D CT images [26], segmentation of breast MR images [27], assessment of regional and global wall motion abnormalities in echocardiography images [28], and detection of polyps in CT colonography images [29, 30].

Isometric feature mapping (Isomap) [17], Locally Linear Embeddings [18, 19], and Laplacian eigenmaps [20] are the most popular NLDR algorithms. In this paper, we have used Isomap algorithm. This algorithm attempts to extract low dimensional parameterizations for data sets in a high dimensional space in such a way that pairwise geodesic distances are preserved so that nearby and far points in high dimensional space map to nearby and far points in low dimensional space [17]. The most common image distance metric for computing the pairwise geodesic distances in Isomap algorithm is the square root of the sum of squared pixel intensities difference between homologous pixels [17]. However, when image sets vary due to the nonrigid deformation of a specific object, an alternative image distance metric based on measuring similarity in terms of nonrigid deformation in image sets is required [22, 31]. This leads to capturing the variation of image sets in the resulting embedding more accurately. Echocardiography images of one cardiac cycle vary due to the nonrigid deformation of the LV caused by the subject’s heartbeat. Therefore, in this study, we have used an image distance metric based on the parameters of the nonrigid transformation model produced by a fully automated image registration method. This nonrigid transformation model is based on an affine transformation for modelling the global LV motion and a B-spline free-form deformation transformation for modelling the local LV deformation [32, 33]. Image registration as an optimization problem finds the optimal parameters of the transformation model so that the mapped image at the time t₁, is as close as possible to the reference image taken at time t₂. In conclusion, the final parameters of nonrigid transformation model calculated by image registration method between image pairs are used to define pairwise image distance measure [28].

In summary, suppose that there are n images as the high dimensional data points, the improved version of Isomap algorithm is as follows:

Compute the distance between all pairs of images in the input space using the method described above. Then, create a neighborhood graph, in which each image is taken as a node and is connected to its k nearest neighbors (the value of k is specified by the user and usually selected 5 to 10). Thus, in this graph, the edge between two nodes exists if these two nodes are neighbors and the length of this edge is the distance between the two images. After that, the shortest path for each pair of nodes through neighborhood graph named geodesic distance is computed using Dijkstra’s or Floyd-Warshall shortest path algorithm [34, 35]. By calculation of geodesic distance between all pairs of nodes, the matrix of pairwise geodesic distance is constructed. The final step applies multidimensional scaling (MDS) [36] on the resulting geodesic distance matrix to construct the low dimensional data points.

Statistical Analysis

In this paper, the mean absolute distance (MAD) is used for measuring the difference between paired curve symbols as follows:

where L_C and L_H are the sets of symbols along the curve of LV volume changes obtained by the proposed method (or other methods) and the reference one, respectively. N is the number of symbols in each curve, and i is the number of sequence of data set. The overall performance measure for each of method is the averaged distance on the whole test set of sequences. We also record the standard deviation of the distance metric for across sequences.

Visualization of Echocardiography Images

By manually defining a rectangular region of interest (ROI) around the LV on ED image, the approximate region of the LV is extracted from echocardiography images (Fig. 3). This image is chosen because of its maximum LV volume at one cardiac cycle. Then, the coordinates of extracted rectangular ROI are applied on all images of a cardiac cycle automatically. Figure 4 demonstrates ten selected images over one cycle of heart in two subjects. These images are considered as points in a high dimensional space (dimension is equal to the number of pixels in the image). Thus, the dimension of first space is 52000 for these 260 × 200 pixel images. Although the input dimension is very high, these images which are different due to LV nonrigid deformation caused by the subject’s heartbeat have a natural low dimensional parameterization in high dimensional space. This low dimensional parameterization can be extracted by improved Isomap algorithm. Consequently, the echocardiography images of one cycle of heart have been embedded in a 2-D space by the improved Isomap algorithm with k = 8 neighbors (Fig. 5). In this figure, the symbols are the nonlinear projection of the consecutive images of one cycle of heart which are joined by a line based on the frame order and in this way we generate a curve in the 2-D space that is called manifold. In all cases, the cyclic nature of the heart motion leads to a cyclic closed manifold. It is noted that the horizontal and vertical axes do not have units.

Improved Isomap algorithm attempts to extract a low dimensional (2-D) data from LV echocardiography images (52000-D) so that geodesic distances between these images are preserved. Therefore, nearby and far images in high dimensional space map to nearby and far symbols in 2-D space. Consequently, in Fig. 5, manifold symbols visualizes the relationship between these images based on distance between them in high dimensional space or other words based on the nonrigid deformation of the LV which results in LV volume changes. As demonstrated in these manifolds (Fig. 5), changes in LV volume in seven phases of the cardiac cycle are observed. The highly experienced echocardiographer by investigating ECG signal have determined manifold symbols corresponding to the seven phases of the cardiac cycle for two above cases and demonstrated with different symbols and colors in Fig. 5. As illustrated in these figures (especially for the normal case), in three phases of the cardiac cycle (second, fifth, and seventh phases), where there are no remarkable changes in LV volume, the distances between symbols are very small. In addition, during the third and fourth phases of the cardiac cycle, in which there are considerable changes in LV volume, symbols are far apart and the distances between them are high. Besides, there is significant distance between consecutive symbols in the manifold during the first and the sixth phases because of considerable LV volume changes. Moreover, because of the highest volume difference (highest LV deformation) between ED and ES images compared to any other two images, the distance between these two symbols is approximately maximum compared to the distance of any other two symbols (Fig. 5). Finally, with increased LV wall motion in healthy case compared with the patient case (increasing the difference between LV echocardiography images), the distances between manifold symbols have increased and consequently the areas of resultant manifold have increased. Consequently, this new visualization demonstrates hemodynamic information of a cardiac cycle and allows a concise evaluation of various stages of ventricular systolic and diastolic function.

Automatic Computation of LV Volume Changes

As described, manifold symbols visualize the relationship between LV echocardiography images of one cycle of heart based on LV volume changes. Therefore, the curve of LV volume changes can be extracted from the resultant manifold. To do this, after automatically determining the symbol corresponding to ES image using the method described in [24], the Euclidean distance between ES symbol and all other symbols in the resultant manifold is calculated. The curve obtained from this distance calculation based on frame number corresponds to LV volume changes of one cycle of the heart (Fig. 6). The LV volume changes curve obtained from the proposed method for a healthy case is similar to LV volume changes curve observed in Fig. 1 for a normal heart. It is noteworthy that the obtained curve of LV volume changes is not scaled based on milliliter. If we want to determine the exact volume at each symbol of this curve, we must scale this obtained curve. To do this, the inner wall (endocardium) and main axis of the LV from the ED and ES images (extracted using the method described in [24]) are detected by the highly experienced echocardiographer. The detected endocardium contour and main axis of the LV are then used to compute LV volume according to the single plane modified Simpson’s rule (Eq. (1)) in each of these two images. Then, calculated ED and ES volumes are used to scale the derived curve of LV volume changes so that ED volume corresponds to maximum Euclidean distance and ES volume corresponds to minimum Euclidean distance (minimum Euclidean distance is zero). In general, the Euclidean distance of each symbol to ES symbol in the resultant manifold is scaled as follows:

Thus, the curve obtained from the distance of each symbol to ES symbol in the resultant manifold (curve of LV volume changes in Fig. 6) is scaled based on milliliter. Figure 7 demonstrates the obtained curve of LV volume changes after scaling (solid line). This figure also demonstrates the obtained curve of LV volume changes in one complete cardiac cycle for two above cases using the highly experienced echocardiographer (dotted line), standard Isomap method (dashed line), and an automatic echocardiographic image segmentation method using level set [5, 6] (dash-dot line).

For a systematic quantitative evaluation, the curves of LV volume changes obtained by the proposed method and two other methods mentioned above are compared to those obtained by the reference experienced echocardiographer using MAD in 32 echocardiography image sequences. Figure 8 and Table 1 show the results of quantitative evaluation using this metric. Result showed that the proposed method achieved higher performance in term of MAD (lower average and standard deviation) than the automatic echocardiography image segmentation and the standard Isomap methods. Also, we have found that for all analyzed sequences there is agreement between reference experienced echocardiographer and automatic proposed assessment.