Myocardial Motion Analysis of Echocardiography Images using Optical Flow Radial Direction Distribution

Problem statement: Myocardial motion is important information for phys icians in diagnosing cardiac abnormalities. The motion vector of myocard ial can be computed using optical flow technique, which then can be further analyzed based on its mag nitude and angle. In practice, physicians are not concern about the angle of vector itself, but are m ore interested on whether a segment is moving to th e center or not. Approach: Therefore, in this study we propose a relative moti on direction with respect to the center of the cardiac cavity, called radial directi on, which is more useful for diagnosis. The radial rection is computed as the difference between the angle of optical flow at a point of interest and the angle b etween the point and the cavity center. Because of the dif ficulty in performing analysis based solely on indi vi ual vectors, it is helpful to visualize and extract the overall trend by representing motion vectors by th eir angular distribution. Results: This method has been tested on clinical echocardiog raphy sequences and has been shown to be successful in providing a radial d irection profile of every segment for each echocardiographic frame. A comparison between the n ormal angular distribution and the proposed radial direction profile was also presented. Conclusion: The proposed profile was shown to be successful in providing the pattern of segmental motion which is easier for physician to analyze the myocardial moti on compared with the normal angular distribution as we ll as more invariant to segment locations.


INTRODUCTION
In 2004, cardiovascular disease caused 17.1 million deaths, which corresponded to 29% of the global deaths for that year as reported by the World Health Organization (2010). From these data, 7.2 million died because of coronary heart disease, while the rest were affected by stroke. Coronary heart disease, which is also known as "heart attack," occurs when there is an interruption of blood supply to the heart muscle. It can be detected early by observing abnormality in the motion of the left ventricular. Realizing the importance of this motion, many researchers have focused on the estimation, quantification and interpretation of cardiac motion.
In medical practice, a physician observes the left ventricular motion by means of three main imaging modalities: Single Photon Emission Computed Tomography (SPECT), Magnetic Resonance Imaging (MRI) and ultrasound imaging. The use of computed tomography and magnetic resonance imaging on the Left Ventricular (LV) motion analysis has been successful because both modalities generate relatively noise-free images (Lin et al., 2006;Sundar et al., 2009). However, they have a number of drawbacks such as invasive, radiation and problems with device portability. In order to overcome these problems, an ultrasound scanning method was introduced; it is a promising tool that has proven to be useful in many clinical applications. Furthermore, many researchers have been attracted on the use of echo images to analyze myocardial motion using block matching method, 3 dimensional model and non-rigid registration (Boukerroui et al., 2003;Papademetris et al., 2001;Ledesma-Carbayo et al., 2005). Other popular methods for myocardial motion estimation have been based on the computation of optical flow (Veronesi et al., 2006;Duan et al., 2006;Riyadi et al., 2011).
Current research on the cardiac motion of echo images is mainly focused on the improvement of the accuracy of the motion estimation as there have been difficulties in the detection of motion using echocardiography images. On the other hand, quantification of the computed motion will be very useful for physicians to interpret the cardiac motion for making a clinical decision. Few publications report on quantification and interpretation of computed motion on echo images. Ledesma-Carbayo et al. (2005) computed the mean displacement vector and the mean strain tensor of each segment from computed myocardial motion. Both features showed significant differences between normal and abnormal cardiac motion. However, the presented features did not provide a complete segmental profile of cardiac motion in a temporal observation. In another proposed method, Duan et al. (2006) developed a four-dimensional optical flow frame-work to track the endocardial surface using real-time 3 dimensional ultrasound. Some features, such as flow magnitude, radial displacement, circumferential displacement, thickening, circumferential and longitudinal stretch and twist, were computed to provide dynamic measurement of cardiac motion. The extracted information was mostly used to visualize realtime 3 dimensional cardiac motion rather than to detect cardiac abnormality. In clinical practice, a cardiologist does not only need to observe the magnitude of myocardial motion, but also the direction of its motion to ensure an effective LV ejection fraction. Instead of the vector angle itself, in this article we propose a relative motion direction with respect to the center of cavity, called radial direction, which is more useful for diagnosis because this will enable the physician to infer whether a segment is moving to the center or not more clearly. The rest of this study describes the method of optical flow used, the computation of optical flow angle and radial direction and finally we present the comparison of direction profile of optical flow using its angle and proposed radial direction method.
Understanding LV motion: For imaging purpose, the American Heart Association has released four standard views for my ocardial segmentation and nomenclature. One of them is Parasternal Short Axis view (PSAX) as shown in Fig. 1. In this view, there are six segments of the myocardium which are associated with the distribution of three main coronary arteries. For normal cardiac function, the myocardial segments move inward to the cavity center with uniform displacement from end diastole to end-systole (Anderson, 2007).

MATERIALS AND METHODS
The extraction of myocardial motion property based on optical flow involves several steps. They are preprocessing stage, estimation of myocardial motion using optical flow, computation of angle and radial direction and generating the angular distribution profile.
We start with taking two consecutive images of echocardiography sequences. Since it is well-known that echocardiography images contain speckle noise, a quasi Gaussian-Discrete Cosine Transform filter adopted from (Riyadi et al., 2009) is then applied to improve the image quality. Subsequently, we manually select the cavity center and initial points on endocardium and then apply B-spline to generate other points along the endocardium.
After the creation of myocardium boundary, optical flows for every pixel along the endocardial boundary were computed. A statistical averaging approach is implemented to remove outlier data along the boundary and to generate a smoother boundary optical flow vectors. To assist with the interpretation of mycardial motion, we propose a radial direction distribution for each myocardial segment, instead of a single angular distribution of flow vector for all the segments. The radial direction distribution represents a direction pattern of segmental myocardial with respect to the cavity center in such a way that cardiologist can directly know the direction in which a segment moves with respect to the centre of the cavity. The rest of this section describes the method used to generate smoothendocardial optical flow and determine radial direction profile.
Generating smooth endocardial optical flow: We computed the optical flow between two consecutive frames using Horn-Schunck's method. This method was selected because of its simplicity, low computing time and suitability for estimation of small displacement. To compute optical flow, the intensity of the object is assumed constant. This can be represented as: The intensity constraint can be approximated using the 1st order Taylor series expansion since the δx, δy and δt are within small value, as shown below: c y I(x , y , t t) Simplifying Eq. 1 and 2, the basic constraint of the optical flow equation can be articulated as: x y t x y t I dx I dy I dt 0or I u I v I 0 Horn-Schunck's method introduces global smoothness error (E S ) and brightness constancy error (E C ) defined as (Horn and Schuck, 1981): In terms of motion estimation, pixels should not have abrupt changes in magnitude and direction of their vector flow compared with the adjacent pixels. Certain points in the neighborhood that have obviously different magnitudes or directions can be considered as errors due to miscomputation of optical flow. To smooth the variation of the vectors in a neighborhood, we apply statistical averaging to remove the outlier data along the boundary. The outlier data on the neighboring points are successively removed and replaced by the average of neighboring vectors. Figure 2a and b show an example of computed optical flow before and after the smoothing process, respectively.

Determining radial direction:
The two parameters of a motion vector are its magnitude and angle. The vector magnitude represents the displacement of a pixel for two consecutive frames, while the vector angle is its direction in terms of the axial (u) and the ordinal (v) component. The vector angle θ 1 is computed using: In mathematics, the two-arguments function atan2(y, x) is a variant of the arctangen function, which calculates an angle (in radians) between the positive xaxis and a point given by coordinate (x, y). Following the standard mathematical convention, atan2 function results angles that increase counter-clockwise. However, since the optical flow computation uses clockwise convention, the vector angle of optical flow θ 1 also follows the clockwise convention as shown in Fig. 3. The atan2 function produces results in the range (-π, π] that can be convert to [0, 2π) by adding 2π to negative value of θ 1 .
In clinical diagnoses, the physician usually does not use the vector direction to analyze the myocardial motion. The direction of the vector with respect to the center of the cavity is more useful information because the physician observes whether a segment moves to the center of the cavity or not. Therefore, in this study, we introduce a new motion property, which is called radial direction defined as the direction of a motion vector with respect to the center of the cavity as illustrated in Fig. 3. The radial direction α of a point of interest (x i , y i ) is the difference between the angle of optical flow at the point of interest (θ 1 ) and the angle between the point and the cavity center (θ 2 ), or defined as: where, the angle θ2 can be calculated by the following equation: The Eq. 7 produces radial direction values in the range (-2π, 2π). To make the value simpler and more informative, we normalize this value using the following rules: if 2 θ −π≤θ ≤ π   θ = θ + π − π≤ θ < −π   θ − π π >θ ≥ π  (9) therefore the radial direction angle α will be in the range (0, π). The 0 means that the optical flow vector direction is towards the cavity center whereas the π indicates the opposite direction.
To facilitate the interpretation of myocardial motion, we generated the distribution of radial direction for every segment for each frame. This distribution shows the pattern of direction in such a way that a physician will be able to quickly recognize the direction to which a segment moves.

RESULTS
Smoothed endocardial optical flow: Our proposed method was tested with a set of clinical echocardiography data, which were previously recorded by a physician. The data were acquired on a parasternal short axis standard view using a Siemens Acuson Sequoia scanner with a 25 frame rate. Manual interpretations on the recorded data have been also performed by two cardiologists to provide manual diagnosis on the myocardial motion. At the end, this manual diagnosis will be used to validate the result obtained by the proposed method. Figure 4 shows an example of four frames of ten consecutive sequences of PSAX view echocardiography images. We then computed the optical flow vectors of two consecutive images (frames 2 and 3) and applied the proposed method to these optical flow vectors.
In the Materials and Methods, we have seen that smooth and correct directions were obtained for the pixels located on the endocardium, whereas many incorrect optical flow vectors occurred inside the myocardium. For this reason, we should extract only the optical flow vectors on the endocardium to provide the myocardial motion profile precisely. Prior to extracting the endocardial features, we have to manually select a few points distributed along the boundary of the cavity and then generate other points along the boundary using the B-spline algorithm as shown in Fig. 5a and b, respectively. In the next stage, the vector flows on the boundary were obtained by averaging 5×5 neighboring vector flows along the boundary and then smoothed as shown in Fig. 5c.

Distribution of radial direction:
To illustrate the midanterior-septum (upper center segment) of Fig. 5c, the zoomed image of the motion vectors and their angle distribution are shown in Fig. 6a and c, respectively. This angle distribution shows that these vectors have dominant directions in the interval 70-90°. The standard vector angle convention follows a clockwise direction orientation where 0, 90, 180 and 270° indicate the right, down, left and up direction, respectively. According to this convention, a segment with the dominant directions in the interval 70-90° moves downward to the center. Similarly, Fig. 6d shows an angle distribution of the mid-posterior (lower center segment) in Fig. 6b. The dominant direction is in the interval 250-290° and it indicates that part of the segment moves upward to the center.

DISCUSSION
These two examples in Fig. 6a-c and b-d show that although the value of dominant directions is largely different, both segments move towards the center. In other words, because the location of mid-anteriorsegment is on the upper center of myocardial, physician should remember that it will move towards the center only if the angle direction is nearby 90°. For the midposterior which is on the lower center segment, it moves towards the center if the dominant direction is about 270°. To enable such an inference to be made, it is necessary for the physician to remember the angle convention and also the parts of the cardiac segment that are being referred to. Hence, it is not easy to correlate the angle distribution and the center-directed motion for a particular cardiac segment.
If we evaluate the distributions of motion vector based on radial direction Fig. 6e and f, they are easier for a physician to analyze myocardial motion instead of the angle of the vector flow itself. These figures show that both segments which move towards the center have a same dominant angle, i.e., in the interval 0° and 30°. Physicians do not need to remember the angle convention and also the location of the segment of interest. They just to know a simple convention that a 0° radial direction means the vector moves to the center of the cavity whereas a 180° radial direction means that the vector is moving to the opposite direction of the cavity center. We conclude that the radial direction is more invariant to the segment location because we have considered the vector location with respect to the center in the computation of the radial direction. Additionally, using this orientation, the pattern of segmental motion direction can be simply obtained by considering the dominant value of the radial direction distribution.

CONCLUSION
A new myocardial segmental motion property based on radial direction has been described.. The radial direction is computed by considering the angle of optical flow at a point of interest and the normal angle between the point and the cavity center. The individual radial direction along the boundary is then used to generate a radial direction profile for every myocardial segment of each frame. The result shows that the radial direction profile is easier for physician to analyze the myocardial motion instead of normal angular distribution as well as more invariant to segment locations.