Reprinted from

REGULAR PAPER

Automated Identification of the Heart Wall Throughout the Entire Cardiac Cycle
Using Optimal Cardiac Phase for Extracted Features
Hiroki Takahashi, Hideyuki Hasegawa, and Hiroshi Kanai

Jpn. J. Appl. Phys. 50 (2011) 07HF16# 2011 The Japan Society of Applied Physics

Person-to-person distribution (up to 10 persons) by the author only. Not permitted for publication for institutional repositories or on personal Web sites.

REGULAR PAPER

Japanese Journal of Applied Physics 50 (2011) 07HF16
DOI: 10.1143/JJAP.50.07HF16

Automated Identification of the Heart Wall Throughout the Entire Cardiac Cycle
Using Optimal Cardiac Phase for Extracted Features
Hiroki Takahashi1 , Hideyuki Hasegawa1;2 , and Hiroshi Kanai1;2 
12Graduate School of Biomedical Engineering, Tohoku University, Sendai 980-8579, Japan
Graduate School of Engineering, Tohoku University, Sendai 980-8579, Japan

Received December 13, 2010; accepted February 7, 2011; published online July 20, 2011
In most methods for evaluation of cardiac function based on echocardiography, the heart wall is currently identified manually by an operator.
However, this task is very time-consuming and suffers from inter- and intraobserver variability. The present paper proposes a method that uses
multiple features of ultrasonic echo signals for automated identification of the heart wall region throughout an entire cardiac cycle. In addition, theoptimal cardiac phase to select a frame of interest, i.e., the frame for the initiation of tracking, was determined. The heart wall region at the frame ofinterest in this cardiac phase was identified by the expectation–maximization (EM) algorithm, and heart wall regions in the following frames wereidentified by tracking each point classified in the initial frame as the heart wall region using the phased tracking method. The results for two
subjects indicate the feasibility of the proposed method in the longitudinal axis view of the heart. # 2011 The Japan Society of Applied Physics

1. Introduction

In conventional echocardiography, structural change andmacroscopic motion of the heart are diagnosed by measuringcross-sectional images, and quantitative evaluation of intracardiac blood ﬂow and cardiac wall motion can also be
performed by ultrasound Doppler measurements.1,2) Moreover, various valuable methods for evaluation of cardiac
function based on echocardiography, such as evaluation of
two-dimensional cardiac wall motion and strain rate by the
speckle tracking method,3–6) measurement of the propagationof vibration caused by the closure of valves,7) and measurement of the transient of myocardial contraction and relaxation
for about 10 ms,8,9) have been developed. In most of theabove-mentioned methods, the heart wall, which is the object
to be analyzed, is currently identiﬁed manually by anoperator. However, this task is very time-consuming andsuﬀers from inter- and intraobserver variability.10–12) For
elimination of operator dependence and facilitation of analysis, automated identiﬁcation of the heart wall is essential.
Various studies on segmenting regions in ultrasonicimages have been conducted. Papadogiorgaki et al. proposed techniques for automated detection of lumen andmedia–adventitia borders for intravascular ultrasonography(IVUS).13) Chang et al. classiﬁed the thyroid gland anddirectly estimated its volume in ultrasonic images.14) In thecited studies, objects in an ultrasonic image are segmented
using features extracted from the magnitude of the received
ultrasonic signal. Furthermore, Torp et al. proposed amethod for automatic detection of several candidate points
in the apical view of the left ventricle.15) This method also
extracts features from the gray-scale value, the velocity, andthe depth to detect the candidate points. The object ofinterest, i.e., the heart wall, can be identiﬁed automatically
based on features and can be identiﬁed manually using
observable features such as brightness and macroscopic
motion in an ultrasonic image. Therefore, the extraction ofvalid features from acquired ultrasonic images (or signals) is
important for the identiﬁcation of objects of interest.
Nillesen et al. introduced a method in which echogenicity
was smoothed in each local area to distinguish the heart wallfrom lumen.16) However, the method using only echogeni

E-mail address: kanai@ecei.tohoku.ac.jp

city would misclassify a region with a low echogenicity
inside the heart wall as lumen. To overcome this limitation,
Kinugawa et al. proposed a method of identiﬁcation of theheart wall using coherence between echoes.17) However, theseparability of this method is degraded by a stationary
component owing to echoes from external tissue such as ribs
because this component increases the coherence of echoesfrom the heart wall and lumen. The present paper introduces
multiple features to reduce the inﬂuences of such undesirable echo characteristics, which would degrade the accuracyof the identiﬁcation of the heart wall.
Segmentation of the ventricular endocardial surface ineach frame in an echocardiographic image sequence during
one cardiac cycle is made possible by tracking the positionof the endocardial surface determined at a certain frame12)
using one-dimensional (1D) or two-dimensional (2D)tracking methods. In the present study, the heart wall region
at a frame of interest (the frame where tracking is initiated)
was identiﬁed automatically by classiﬁcation using multiple
features extracted from echo signals, which are measured at
a high frame rate. Furthermore, to select the frame ofinterest, we identiﬁed the optimal cardiac phase forclassiﬁcation by multiple features.
2. Principles
2.1 Feature extraction

Figure 1(a) illustrates a plane in the heart, which wasscanned by an ultrasound beam. With respect to each scanline position j, as shown in Fig. 1(b), ultrasonic RF echo
signals scattered by the heart wall [e.g., interventricular
septum (IVS) and left ventricular posterior wall (LVPW)]
and lumen [e.g., right ventricle (RV) and left ventricle (LV)]are received with a time interval corresponding to theinverse of a frame rate.
Just like in conventional methods for segmentation of
echocardiographic images, the magnitude of the echo signal
was used as the ﬁrst feature in the present study. The
envelope-detected signal Envðx; nÞ, at a depth of x in the n-thframe, is extracted from the received ultrasonic RF echo
signal sðx; nÞ as follows:Envðx; nÞ ¼ log10 ½LPF½ðsðx; nÞÞ2 ;ð1Þwhere LPF½ shows low-pass ﬁltering with a cutoﬀ
frequency of 0.3 MHz (the transmit frequency of 3.75 MHz).
07HF16-1# 2011 The Japan Society of Applied Physics

Person-to-person distribution (up to 10 persons) by the author only. Not permitted for publication for institutional repositories or on personal Web sites.

Jpn. J. Appl. Phys. 50 (2011) 07HF16

H. Takahashi et al.

sector probe
interventricular
septumright ventricle

left ventricle
left ventricularposterior wallNjj

1 0 scan line No.

(a)(the distance from the probe)

direction of depth

direction of frame
n

n −1

n +1

n+ k

frame No.

ii+1

sampling Δxintervalmeasured RF signals at j−th scan line(b)Fig. 1. (Color online) (a) Geometry of measurement by transthoracic
echocardiography. (b) Illustration of variables of RF signals.

Nillesen et al. reduced the overlap between gray levels,
corresponding to the magnitude of the backscattered signals,
in the heart wall and lumen from 35.4 to 19.6% bysmoothing in each local kernel for segmentation of the heart
muscle in echocardiographic images.16) However, using
such methods, which depend only on the amplitude of thejð f ; x; nÞj2 ¼

RF echo, the heart wall is not completely distinguishedfrom the lumen because there is a region with a low echo
level inside the heart wall. Hete and Shung found thatthe magnitude of the ultrasonic backscattered signals waschanged by the orientation of muscle collagen ﬁbers.18)
Wickline et al. suggested that the magnitude of ultrasonicbackscatter is determined by local diﬀerences in acoustic
impedance, which are aﬀected by myocardial contraction
and relaxation.19) Thus, additional second and third features
were adopted in this study to evaluate diﬀerences between
echo signals of the heart wall and lumen regardless of theamplitude of backscattered signals as described below. In
the ultrasonic RF signal measured at a high frame rate, thecardiac muscles in the heart wall can be tracked accuratelyby the phased tracking method.20) As shown in Fig. 2, thetemporal changes of waveforms of RF echo signals in aregion of interest (ROI) placed in the heart wall region are
fairly constant and can be tracked because the myocardium
in the ROI stays in the focal area of the ultrasonic beamin two consecutive frames. On the other hand, blood cells
in the cardiac lumen are diﬃcult to track because theyslip oﬀ from the focal area of the ultrasonic beam owing
to blood ﬂow. Therefore, the temporal change in the
waveform of an ultrasonic echo in an ROI assigned in thelumen is signiﬁcant. The displacement ux ðn þ kÞ between
the n-th and (n þ k)-th frames of tissue in an ROI assigned
around depth x in the n-th frame is estimated by the phased^ þ kÞ; n þ kÞ
tracking method, and complex spectra, Yð f ; xðn^ þ k þ 1Þ; n þ k þ 1Þ, of RF signals in theand Yð f ; xðn^ þ kÞ ½¼ x þ ux ðn þ kÞ and
ROIs at the positions, xðn^ þ k þ 1Þ ½¼ x þ ux ðn þ k þ 1Þ, in the (n þ k)-thxðnframe and (n þ k þ 1)-th frame are obtained by applying
the discrete Fourier transform (DFT) to RF signals in^ þ kÞ; n þ kÞg,the ROIs. Using the complex spectra fYð f ; xðnthe characteristics of phase changes of echoes in the
ROI between two consecutive frames, which depends onthe motion of an object, are evaluated by the magnitudesquared coherence (MSC) function17,21) jð f ; x; nÞj2 asfollows:^ þ kÞ; n þ kÞYð f ; xðn^ þ k þ 1Þ; n þ k þ 1Þj2
jEk ½Y  ð f ; xðn;^ þ kÞ; n þ kÞj2 Ek ½jYð f ; xðn^ þ k þ 1Þ; n þ k þ 1Þj2 
Ek ½jYð f ; xðn

where Ek ½ and  denote the time averaging (for 50 frames)
and complex conjugate, respectively. The numerator ofeq. (2) denotes the absolute value of the mean of thecomplex cross spectra. Therefore, the consistency of thetemporal changes in phases of complex spectra between twoframes signiﬁcantly aﬀects this averaging operation shownin the numerator of eq. (2). Thus, eq. (2) evaluates the
temporal variance in changes in the complex spectra of RF
signals in the ROI between two consecutive frames.17)
Kinugawa et al. reported that the variance in the phasechanges in the lumen shows relatively large values
compared with that of the heart wall in the high frequency
range.17) In this study, the MSC at a frequency of 4.4 MHz,
which was higher than the transmit frequency of 3.75 MHz,was used as the second feature.
As the third feature, the temporal mean of the changes
in phases of RF signals between two consecutive frames

ð2Þ

was adopted because it cannot be evaluated by the MSC.
The phase shift corresponds to the axial velocity of tissue
and is obtained using the complex correlation technique.22)
Hatle and Angelsen showed that the peak blood velocity
during systole of the left ventricle ranges from 0.7 to
1.1 m/s, while that during diastole of the left ventricle
ranges from 0.6 to 1.3 m/s using the ultrasonic Doppler
method.23) The absolute value of the peak of myocardial
velocities (IVS and LVPW), which were measured by
Kapusta et al. using tissue Doppler imaging, was up to
0.2 m/s.24) In this study, the phase shift due to the axial
motion of tissue was extracted by estimation based oncorrelation of analytic signals fya ðx; nÞg of RF signals
fsðx; nÞg, where analytic signals fya ðx; nÞg were obtainedby applying the quadrature demodulation to RF signals.
The complex correlation function, Cðx; nÞ, is calculated asfollows:

07HF16-2# 2011 The Japan Society of Applied Physics

Person-to-person distribution (up to 10 persons) by the author only. Not permitted for publication for institutional repositories or on personal Web sites.

Jpn. J. Appl. Phys. 50 (2011) 07HF16

H. Takahashi et al.

Imtt+ Δ tt+ 2Δtultrasonic beamreceived
signalhigh consistencyReheart wall
traceable

DFT

spectra in heart wallregion ofinterest

DFT

Im

untraceable

lumen
blood celltdispersion

t+ Δ ta few [ms]t+ Δ ttt+ 2ΔtRe

time [s]t+ 2Δtspectra in lumen

Cðx; nÞ
M
t 1
X

(a)

ya ðx; n þ kÞ  ya ðx; n þ k þ 1Þ

10 mm

Fig. 2. (Color online) Illustration of changes of waveforms of echo signals in ROI.

RV
IVS

k¼M

t
ﬃvﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ ;¼ vﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
u Mt 1u Mt 1u X
u X
tjya ðx; n þ kÞj2 t
jya ðx; n þ k þ 1Þj2

k¼Mt

M mode

LV
LVPW

k¼Mt

ð3Þ
where 2Mt determines the number of frames used for
calculation of the correlation function and is set at 10 in
the present study. The absolute value of the phase shift,jðx; nÞj, which is obtained by spatially averaging theestimator of eq. (3), was used as the third feature as
follows:

!
M


d 1
X


ð4Þ
jðx; nÞj ¼ <
Cðx þ mx; nÞ ;


m¼M

ECG
PCG

(b)

2.2 Enhancement of accuracy of features by filtering

The MSC of RF signals, the second feature, is increased
owing to the components of echoes from external tissue

time

shift frequency [Hz]dwhere x and 2Md are the interval of the sampled signal
in the axial direction and the number of sampled points
used for spatial averaging (2Md  x corresponds to apulse duration deﬁned by the width at 20 dB of theenvelope of an ultrasonic pulse, about 1.44 mm).
By replacing x with the sampled point number i ð i  xÞ
in the direction of depth and using the scan line number j, theenvelope-detected signal Envðx; nÞ, the MSC jð f ; x; nÞj2 at
4.4 MHz, and the absolute value of the phase shift jðx; nÞjat the frame of interest are redeﬁned as Envði; jÞ, jði; jÞj2 ,and jði; jÞj, respectively.
Kinugawa et al. proposed a method for the identiﬁcationof the heart wall using only MSC as the feature.17) However,the MSC of echoes from the heart wall and lumen are
sensitively increased by a stationary component owing toechoes from the ribs. In this study, therefore, the multiple
features (e.g., the amplitude and the phase change ofultrasonic RF signals) were extracted from ultrasoundechoes to overcome such problems, e.g., a low echogenicity
inside the heart wall and a stationary component owing toechoes mainly from the ribs.

130 msR−wave

400
350
300
250
200
150
100
50
0

power [dB]

0 0.2 0.4 0.6 0.8 1 1.2 1.4

40
30
20
10
0
−10
−20
−30

time from R−wave [s]Fig. 3. (Color online) (a) M-mode image (RV: right ventricle, IVS:
interventricular septum, LV: left ventricle, LVPW: left ventricular posteriorwall). Yellow lines show tracking positions. (b) Spatial averaged power
spectrum distribution obtained by applying DFT in the direction of frame(50 frames) to the yellow lines shown in (a).

such as the ribs17) contained in the received RF signals. In
addition, there is a considerable diﬀerence between the
acoustic impedances of a rib and soft biological tissue.
Therefore, the amplitude of the echo signal, the ﬁrst feature,backscattered particularly from a shallow region in theheart, also becomes higher because it contains echoes fromthe ribs. Thus, reduction of the stationary components fromthe external tissue is necessary to accurately extract thefeatures of the heart wall and lumen. In this study, to reduce
this component, a third-order Butterworth high-pass ﬁlter25)
was used for the moving target indicator (MTI) ﬁltering26) of
RF signals before calculating the ﬁrst and second features,Envði; jÞ and jði; jÞj2 . Figure 3 shows the spatially averaged
power spectrum distribution obtained by applying DFT with

07HF16-3# 2011 The Japan Society of Applied Physics

Person-to-person distribution (up to 10 persons) by the author only. Not permitted for publication for institutional repositories or on personal Web sites.

Jpn. J. Appl. Phys. 50 (2011) 07HF16

H. Takahashi et al.

RF data (corresponding to US image sequence)

Steps for classification
Assigning a point in each region
for initial parameters of EM algorithm

Selecting a frame of interest

Feature extraction
RF signals

MTI filtering

MTI filtering

for the component of external tissue

Envelope

for the clutter component

MSC

Phase shift

Extracting the initial parameters

Classification of the region
by using EM algorithm

Tracking points classified into heart wallusing 1−D or 2−D tracking method

Identification of heart wall at the other frames
Fig. 4. Schematic procedure of identiﬁcation of the heart wall proposed in this study.

respect to the direction of the frame to RF echoes from theposition of the RV in each frame, which was acquired at aframe rate of 751 Hz and contained many of the echo
components from the ribs. The displacement of each discrete
point to be analyzed by DFT was tracked using the phased
tracking method. As shown in Fig. 3(b), RF echoes in the RV
mainly have components of approximately 0 Hz (stationary).
It is necessary to reduce the components of the heart wall by
MTI ﬁltering as little as possible. Therefore, in this study, thecutoﬀ frequency fc1 for the ﬁrst and second features of thehigh-pass ﬁlter was set to 10 Hz (approximately 0 Hz).
On the other hand, for extracting the phase shift of echoes
caused by blood ﬂow in the lumen using the third featurejði; jÞj, clutter components, which are signals scatteredfrom slowly moving tissues except for blood cells, such as
ribs and cardiac muscle, need to be reduced in the received
RF signal. Therefore, the cutoﬀ frequency, fc2 , of the MTI
ﬁlter for the third feature has to be changed by referring to
the velocity of cardiac muscle for extracting only the echoesfrom blood particles. In this study, the optimal cardiac phase
was examined for various values of fc2 .
2.3 Classification method

Kinugawa et al. identiﬁed the heart wall on the basis of Bayes

decision rule using probability distributions of the MSC ofechoes from the heart wall and lumen, which had been
manually segmented.17) However, the parameters, used in themethod, of probability distributions could not consider features of the other points except regions manually segmented.
In this study, an expectation–maximization (EM) algorithm27)
for a mixture of Gaussian distributions was used to classify
feature vectors, which were obtained from RF signals in aframe of interest, into either the heart wall or the lumen. Using
the EM algorithm, parameters of the probability distribution,
which better ﬁts features of all points, are estimated. The
probability distribution Pðdi;j Þ of a feature vector di;j , which is
obtained from either of 2 classes, can be represented using a
mixture of two Gaussian distributions as follows:Pðdi;j Þ ¼ 0 gðdi;j j0 ; 0 Þ þ 1 gðdi;j j1 ; 1 Þ;

ð5Þ

where {h }, {h }, {h }, and gðÞ denote mixing coeﬃcients,
mean vectors, covariance matrices [lumen (h ¼ 0) and heart
wall (h ¼ 1)], and Gaussian distribution, respectively. These
parameters {h }, {h }, and {h }, which maximize thelikelihood for a mixture of Gaussian distributions shown in
eq. (5), can be estimated using an EM algorithm, one of theiterative calculation methods, based on extracted feature
vectors. In this study, the initial parameters for the EM

07HF16-4# 2011 The Japan Society of Applied Physics

Person-to-person distribution (up to 10 persons) by the author only. Not permitted for publication for institutional repositories or on personal Web sites.

Jpn. J. Appl. Phys. 50 (2011) 07HF16

10 mm

H. Takahashi et al.
20 ms

10 mm

R−waveRF aquisition area
RV

RV
IVS

IVS
LV

R−wave
1st heart sound

2nd heart sound

ECG
PCG

LV
10 mm

LVPW

LVPW

(a)

region for evaluation of heart wall(b)

region for evaluation of lumen(a)Fig. 5. (Color online) Cross-sectional images of the left ventricle inlongitudinal axis view obtained by conventional echocadiography for
subjects (a) A and (b) B. The regions surrounded by a red line show
acquisition areas of RF data.

20 ms
R−wave

1st heart soundR−wave

2nd heart sound

ECG
PCG

algorithm were obtained using feature vectors extracted fromthe area around the point manually assigned in each region ofa frame of interest. Using the estimated parameters, theposterior probability, ði; jÞ, of the heart wall at each discrete
point ði; jÞ is estimated as follows:1 gðdi;j j1 ; 1 Þ:ði; jÞ ¼
0 gðdi;j j0 ; 0 Þ þ 1 gðdi;j j1 ; 1 Þ

10 mmregion for evaluation of heart wallregion for evaluation of lumen

ð6Þ

In this study, discrete points were classiﬁed into the heart wall(h ¼ 1) when ði; jÞ exceeded 0.5.
Figure 4 shows the above-mentioned procedure proposedin this study for automated identiﬁcation of the heart wall
during an entire cardiac cycle. In the proposed method,
selection of a frame of interest, i.e., the frame for thedetermination of regions to be tracked, is important because
it aﬀects the results of the identiﬁcation of the heart wall in
all following frames. Kinugawa et al. showed the results ofidentiﬁcation of the heart wall based on the MSC in several
cardiac phases.17) However, the separability of the feature,
MSC, in each cardiac phase was not evaluated. To
accurately identify the heart wall in a frame of interest, inthe present study, we investigated the optimal cardiac phasefor classiﬁcation as described below.

(b)Fig. 6. (Color online) Regions manually assigned for evaluation of
separability between heart wall and lumen for subjects (a) A and (b) B.
Heart walls (IVS and LVPW) and lumen (LV) are shown as the areas
surrounded by red and blue lines, respectively.

sampling frequency of the measured RF signals was
15 MHz. Using the three-dimensional feature vector variable
di;j ¼ ðEnvði; jÞ; jði; jÞj2 ; jði; jÞjÞ extracted at each discrete
point in the acquisition areas, prior probabilities {ph }, mean
vectors {h }, and covariance matrices {h } are calculated
for classes {!h } [lumen (h ¼ 0) and the heart wall (h ¼ 1)]
as follows:

3. Evaluation of Separability in Each Cardiac PhaseNhðh 2 f0; 1gÞ;N
h ¼ Ei;j ½di;j j!h  ðh 2 f0; 1gÞ;3.1 High frame rate acquisition of RF signals

h ¼ Ei;j ½ðdi;j  h Þðdi;j  h Þ j!h TTo accurately diﬀerentiate the heart wall from the lumen
using features that are needed to estimate phase shift, it isessential that RF echo signals be measured at a highframe rate for the prevention of aliasing in the phased
tracking method. To do that, the number of transmits,which is required to obtain the same number of scan lines asthat in conventional sector scanning, was reduced bycreating many focused receiving beams in one unfocusedtransmitting beam such as a plane wave [parallel beamforming (PBF)].28,29)
3.2 Evaluation of separability

Figures 5(a) and 5(b) show acquisition areas of RF echoes in
the longitudinal axis view of the hearts of two healthy 23year-old males (subjects A and B) using a 3.75 MHz sectortype probe. Using PBF with plane wave transmission, thehigh frame rate (1010 Hz) measurement of RF echo data wasrealized with 80 scan lines with intervals of 0.375 . The

ð7Þ

ph ¼

¼

1
X

ph Ei;j ½di;j j!h ;

ð8Þ
ðh 2 f0; 1gÞ;

ð9Þ
ð10Þ

h¼0

where Ei;j ½ denotes the spatial averaging in the axial andlateral directions, and Nh , N, and  are the number offeature vectors in each class, the total number of vectors in
the area for evaluation, and the mean vector in all classes,
respectively. As shown in Figs. 6(a) and 6(b), the class !h ,to which each feature vector di;j belongs, was denoted bymanually assigning areas in the heart wall and lumen, i.e.,
the IVS, LVPW, and LV, at intervals of 20 ms from the Rwave to the next R-wave. To evaluate the separability of thefeature vector, the inter-class covariance matrix SB ineq. (11), which corresponds to the distance between mean
vectors of classes, and the intra-class covariance matrix SW
in eq. (12), which corresponds to the mean of the variancesof feature vectors within classes, were deﬁned using theseparameters as follows:

07HF16-5# 2011 The Japan Society of Applied Physics

Person-to-person distribution (up to 10 persons) by the author only. Not permitted for publication for institutional repositories or on personal Web sites.

Jpn. J. Appl. Phys. 50 (2011) 07HF16

R−waveΙΙ ΙΙΙ

Ι

H. Takahashi et al.
Ι

IV

(a)

1

normalized
magnitude [dB]

MSC
IVS
LVPW
LV

0.4
0
−20

normalized
magnitude [dB]

MSC

0.6
0.2

velocity [mm/s]

12
10
8
6
4
2

(b)

0.8

strain rate [(m/s)/m]

cutoff frequency [Hz]

300
250
200
150
100
50
0

(c)
−40
−60
−80
−150
−100
−50
0
50
100
150
10

(d)

(e)

5
0
−5
−10

0

0.2

0.4

0.6

0.8

1

1.2

time from R−wave [s]Fig. 7. (Color online) Results of subject A: (a) Criterion of separability
J. Spatial means of (b) the MSC and (c) the envelope in each region
manually assigned. (d) Velocities and (e) strain rate in IVS and LVPW on
the beam shown by the yellow line in Fig. 6(a).

SB ¼
SW ¼

IV

ECG
PCG

criterion of
separability

strain rate [(m/s)/m] velocity [mm/s]

cutoff frequency [Hz]

ECG
PCG

−15

ΙΙ ΙΙΙ

R−wave

1
X
h¼0
1
X

ph ðh  Þðh  ÞT ;ð11Þ

ph h ;

ð12Þ

h¼0

The separability of the feature vector is increased inproportion to the distance between mean vectors of classes,
whereas it is decreased in proportion to the variances of thevectors within classes.30) Therefore, the criterion of separability J is deﬁned by normalizing the inter-class covariancematrix SB by the intra-class covariance matrix SW as
follows:31)
J ¼ tr½S1
W SB ;

ð13Þ

where tr½ represents the sum of diagonal components.
Figure 7(a) shows the criterion of separability, J, ofsubject A estimated using diﬀerent values of fc2 of MTI
ﬁltering for the third feature. The spatial means of theenvelope Envði; jÞ and the MSC jði; jÞj2 in each frame,
which correspond to the components of h in eq. (8), ofsubject A are shown in Figs. 7(b) and 7(c). In addition,

criterion of
separability
12
(a)

300
250
200
150
100
50
0
1
0.8
0.6
0.4
0.2
0
−20

10
8
6
4
2

(b)
IVS
LVPW
LV

(c)

−40
−60
−80
−150
−100
−50
0
50
100
150
10
5
0
−5
−10
−15

(d)

(e)

0

0.2

0.4

0.6

0.8

1

time from R−wave [s]Fig. 8. (Color online) Results of subject B: (a) Criterion of separability
J. Spatial means of (b) the MSC and (c) the envelope in each region
manually assigned. (d) Velocities and (e) strain rate in IVS and LVPW on
the beam shown by the yellow line in Fig. 6(b).

Figs. 7(d) and 7(e) show velocities and strain rates in theaxial direction of IVS and LV estimated by the phased
tracking method8,20) along the scan line shown by the yellow
line in Fig. 6(a). Similarly, these values of subject B are
shown in Figs. 8(a)–8(e). Periods I, II, III, and IV denote theejection phase, isovolumic relaxation phase, rapid ﬁllingphase, and slow ﬁlling phase in Figs. 7 and 8.
3.3 Discussion of determination of optimal cardiac phase
to select a frame of interest

As shown in both Figs. 7(a) and 8(a), results of subject A
and subject B, the maximal values of J are found in thetransition period, which is approximately 100 ms from therapid ﬁlling phase (III) to the slow ﬁlling phase (IV). As
shown in Figs. 7(b) and 8(b), this period shows greater
separability owing to the large diﬀerence in the MSC
between the heart wall region and the lumen region. The
MSC in the lumen in this period is small because the blood
velocity is high just after the rapid ﬁlling phase. Moreover,
the MSC in the heart wall in this period is large because thestrain rate of the cardiac muscle, which is a deformation oftissue changing the spacing of scatters,32) is low, as shown in
Fig. 7(e). The decrease in the strain rate of cardiac muscle in

07HF16-6# 2011 The Japan Society of Applied Physics

Person-to-person distribution (up to 10 persons) by the author only. Not permitted for publication for institutional repositories or on personal Web sites.

H. Takahashi et al.

(b)

1

1

magnitude−squared coherence

(a)

0.8
0.6
0.4
0.2
0

−80
−60
−40
normalized magnitude [dB]

0.8
0.6
0.4

0

(d)

1

1

decrease of the MSC
in lumen
0.8
0.6
0.4
0.2
0

decrease of the MSC
in heart wall
−80

−60

−40

increase of the MSC
in lumen

0.2(c)magnitude−squared coherence

magnitude−squared coherence

magnitude−squared coherence

Jpn. J. Appl. Phys. 50 (2011) 07HF16

−80
−60
−40
normalized magnitude [dB]

0.8
0.6
0.4

increase of the MSC
in heart wall

0.2
0normalized magnitude [dB]

−80

−60

−40

normalized magnitude [dB]Fig. 9. (Color online) Distributions on two-dimensional feature space constructed by the ﬁrst and second features, i.e., the envelope and the MSC, at eachframe in (a) late ejection phase, (b) isovolumic phase, (c) rapid ﬁlling phase, and (d) the transition period from the rapid ﬁlling phase to the slow ﬁlling phasefor subject A.

this period is found because the rapid extension of LV,33)
which corresponds to the increase in volume of LV, nearlyﬁnishes at the end of the rapid ﬁlling phase. Therefore, thetransition period from the rapid ﬁlling phase to the slow
ﬁlling phase is the valid cardiac phase to select a frame ofinterest in this study. As described above, the rapid ﬁllingphase is not suitable for the selection of a frame of interest
because the MSC in the heart wall is lower owing to anincrease in the strain rate of cardiac muscle accompanied by
rapid extension of LV. From the late ejection phase to theisovolumic relaxation phase, the separability is low owing toan increase in the MSC in the lumen. Scatterers (blood cells)
are likely to remain in an ultrasonic beam between two
consecutive frames because the blood ﬂow velocity is low inthese cardiac phases. Tanaka et al. showed that themagnitude of the velocity vector of blood ﬂow is about
20 cm/s and increased slightly to 50–60 cm/s at the aortic
oriﬁce in late systole.34) The width of the emitted plane waveof PBF is 1.6 mm,6) and blood cells may completely slip oﬀ
from the beam during a period of about 3 to 8 ms in thesecardiac phases under this condition. Therefore, these cardiac
phases, when the blood ﬂow velocity is decreased, areconsidered to show low separability due to a high temporal
resolution of about 1 ms. The above-mentioned changes in
the MSC appear as distributions of feature vectors in thefeature space shown in Fig. 9.
For evaluation of the eﬀect of the change in values of fc2 ,the criterion of separability, J3 , of the third feature, jði; jÞj,is obtained by normalizing the distance between means ofclasses by the variances of the third feature within classes, aswell as by eq. (13), as follows:

1
X

J3 ¼

ph ð3h  3 Þ2

h¼0
1
X

;

ð14Þ

2

ph Ei;j ½ðjði; jÞj  3h Þ j!h 

h¼0

where 3h and 3 are the third component of h in eq. (8)and  in eq. (10), respectively. Figures 10(a) and 10(b)show the criterion of separability, J3 , of respective subjects
estimated using diﬀerent values of fc2 of MTI ﬁltering, as inFigs. 7(a) and 8(a). The separability J3 in each cardiac phaseis changed by the value of cutoﬀ frequency fc2 . In this study,
the cutoﬀ frequency fc2 was set to 200 Hz in the optimal
cardiac phase, the transition period from the rapid ﬁllingphase to the slow ﬁlling phase, because the separability, J3 ,of the third feature in this period was found to be high in therange of the cutoﬀ frequency fc2 from 200 to 300 Hz.
4. Identification of Heart Wall Throughout Entire
Cardiac Cycle

Figures 11(a) and 12(a) show B-mode and region-identiﬁed
images in frames of interest selected in the transition period
from the rapid ﬁlling phase to the slow ﬁlling phase forsubjects A and B. The initial parameters of the EM
algorithm were calculated by extracting feature vectors
from the manually assigned regions, which are surroundedby the red line (heart wall) and the blue line (lumen) in
Figs. 11(a) and 12(a). As shown in Figs. 11(a) and 12(a), theheart wall region was accurately identiﬁed in the frame ofinterest in the optimal cardiac phase. The heart wall in eachframe throughout the entire cardiac cycle was identiﬁed by

07HF16-7# 2011 The Japan Society of Applied Physics

Person-to-person distribution (up to 10 persons) by the author only. Not permitted for publication for institutional repositories or on personal Web sites.

Jpn. J. Appl. Phys. 50 (2011) 07HF16

H. Takahashi et al.
ΙΙΙ ΙΙΙ

IV

ECG
PCG

300
250
200
150
100
50
0

300
250
200
150
100
50
0

0

ΙΙ ΙΙΙ

IV

R−waveECG
PCG
cutoff frequency [Hz]cutoff frequency [Hz]R−waveΙ

0.2 0.4 0.6 0.8 1 1.2time from R−wave [s]separability
4
3
2
1
0

0

0.2

0.4 0.6 0.8
1time from R−wave [s]

(a)

(b)

Fig. 10. (Color online) Criterion of separability J3 of the third feature obtained with diﬀerent values of cutoﬀ frequency fc2 of MTI ﬁltering for subjects(a) A and (b) B.

assigned area in heart wall
for initial parameters of EM algorithm
assigned area in lumen
for initial parameters of EM algorithm

ECG
PCG

(a)

R−wave

1st heart sound

2nd heart soundearlyejection phaselateejection phaseisovolumic
relaxation phase

0.2 s

10 mm

(b)

heart walllumen

isovolumic
contraction phaseFig. 11. (Color online) (a) B-mode and region-identiﬁed images in the frame of interest in the optimal cardiac phase for subject A. (b) Region-identiﬁed
images obtained by tracking the points, which were classiﬁed as heart wall in the frame of interest, in the frame in isovolumic contraction phase (just afterR-wave), in the early ejection phase, in the late ejection phase, and in the isovolumic relaxation phase.

tracking the points identiﬁed as the heart wall region in theframe of interest in Figs. 11(a) and 12(a) using the phased
tracking method, as illustrated in Fig. 4. Figures 11(b) and12(b) show region-identiﬁed images in the frames in
diﬀerent cardiac phases, i.e., the period just after the Rwave, early ejection phase, late ejection phase, andisovolumic relaxation phase, using the phased tracking
method. The obtained results show that the method proposedin this study can identify the heart wall region during anentire cardiac cycle.
A serrated heart wall region is shown in Figs. 11(b) and12(b). It was caused by diﬃculty in tracking points in thelumen around the boundary between the heart wall andlumen misclassiﬁed as the heart wall in the frame of interest.
Therefore, a method for recovering the shape14) of the heart
wall in a frame of interest (or the following frames) should

be further investigated. Moreover, further investigation isrequired for improving the separability of the proposed
method, such as application of the 2D tracking method and amore eﬀective classiﬁcation method (e.g., changing theprobability distribution for the EM algorithm), and improvement of MTI ﬁltering.
5. Conclusions

In this study, we proposed a method for automatedidentiﬁcation of the heart wall region throughout an entire
cardiac cycle and determined the optimal cardiac phase to
select a frame of interest, i.e., a frame for the initiation of
tracking. To determine the optimal cardiac phase, theseparability of classes was quantitatively evaluated in eachcardiac phase for two subjects. The heart wall was accurately
identiﬁed in a frame of interest in the optimal cardiac phase

07HF16-8# 2011 The Japan Society of Applied Physics

Person-to-person distribution (up to 10 persons) by the author only. Not permitted for publication for institutional repositories or on personal Web sites.

Jpn. J. Appl. Phys. 50 (2011) 07HF16

H. Takahashi et al.

assigned area in heart wall
for initial parameters of EM algorithm
assigned area in lumen
for initial parameters of EM algorithm

R−wave

1st heart sound

2nd heart sound

ECG
PCG
0.2 s

(a)
10 mm

(b)

heart walllumen

isovolumic
contraction phaseearly
ejection phaselateejection phaseisovolumic
relaxation phaseFig. 12. (Color online) (a) B-mode and region-identiﬁed images in the frame of interest in the optimal cardiac phase for subject B. (b) Region-identiﬁed
images obtained by tracking the points, which were classiﬁed as heart wall in the frame of interest, in the frame in isovolumic contraction phase (just afterR-wave), in the early ejection phase, in the late ejection phase, and in the isovolumic relaxation phase.

by applying the EM algorithm, i.e., one of the classiﬁcation
methods, to extracted feature vectors, and heart wall regionsin the following frames were also identiﬁed by tracking eachpoint classiﬁed as the heart wall region using the phased
tracking method. The results indicate the feasibility of the
proposed method in the longitudinal axis view. To apply the
proposed method to ultrasonic images in various cardiac
views, the optimal cardiac phase should be further
investigated in the short axis and apical views of the heart.

 Brodin,
15) A. H. Torp, S. I. Rabben, A. Stylen, H. Ihlen, K. Anderson, L. A.
and B. Olstad: Proc. 2004 IEEE Ultrasonics Symp., 2004, p. 474.
16) M. M. Nillesen, R. G. P. Lopata, I. H. Gerrits, L. Kapusta, H. J. Huisman,

J. M. Thussen, and C. L. D. Korte: Ultrasoud Med. Biol. 33 (2007) 1453.

17) T. Kinugawa, H. Hasegawa, and H. Kanai: Jpn. J. Appl. Phys. 47 (2008)

4155.
18) B. Hete and K. K. Shung: IEEE Trans. Ultrason. Ferroelectr. Freq. Control

40 (1993) 354.
19) S. A. Wickline, L. J. Thomas III, J. G. Miller, B. E. Sobel, and J. E. Perez:

J. Clin. Invest. 76 (1985) 2151.
20) H. Kanai, M. Sato, Y. Koiwa, and N. Chubachi: IEEE Trans. Ultrason.

Ferroelectr. Freq. Control 43 (1996) 791.

21) H. Takahashi, H. Hasegawa, and H. Kanai: Cho-onpa Igaku 36 (2009) 679[in Japanese].
22) C. Kasai, K. Namekawa, A. Koyano, and R. Omoto: IEEE Trans. Sonics

Ultrason. 32 (1985) 458.

1) G. R. Sutherland, G. D. Salvo, P. Claus, J. D’hooge, and B. Bijnens: J. Am.
2)
3)
4)

5)

6)
7)
8)
9)
10)
11)
12)
13)
14)

Soc. Echocardiogr. 17 (2004) 788.
W. N. McDicken, G. R. Sutherland, C. M. Moran, and L. N. Gordon:
Ultrasound Med. Biol. 18 (1992) 651.
L. N. Bohs, B. J. Geinman, M. E. Anderson, S. C. Gebhart, and G. E.
Trahey: Ultrasonics 38 (2000) 369.
K. Kaluzynski, X. Chen, S. Y. Emelianov, A. R. Skovoroda, and M.
O’Donnell: IEEE Trans. Ultrason. Ferroelectr. Freq. Control 48 (2001)
1111.
J. D’hooge, E. Konofagou, F. Jamal, A. Heimdal, L. Barrios, B. Bijnens,
J. Theoen, F. V. de Werf, G. R. Sutherland, and P. Suetens: IEEE Trans.
Ultrason. Ferroelectr. Freq. Control 49 (2002) 281.
Y. Honjo, H. Hasegawa, and H. Kanai: Jpn. J. Appl. Phys. 49 (2010)
07HF14.
H. Kanai: IEEE Trans. Ultrason. Ferroelectr. Freq. Control 52 (2005) 1931.
H. Yoshiara, H. Hasegawa, H. Kanai, and M. Tanaka: Jpn. J. Appl. Phys.
46 (2007) 4889.
H. Kanai: Ultrasoud Med. Biol. 35 (2009) 936.
V. Chalana and Y. Kim: IEEE Trans. Med. Imaging 16 (1997) 642.
E. Brandt, L. Wigström, and B. Wranne: Proc. MICCI’99, 1999, p. 410.
Q. Duan, E. D. Angelini, S. L. Herz, C. M. Ingrassia, K. D. Costa, J. W.
Holmes, S. Homma, and A. F. Laine: Ultrasound Med. Biol. 35 (2009) 256.
M. Papadogiorgaki, V. Mezaris, Y. S. Chatzizisis, G. D. Giannoglou, andI. Kompatsiaris: Ultrasound Med. Biol. 34 (2008) 1482.
C. Y. Chang, Y. F. Lei, C. H. Tseng, and S. R. Shih: IEEE Trans. Biomed.
Eng. 57 (2010) 1348.

23) L. Hatle and B. Angelsen: Doppler Ultrasound in Cardiology (Lea &Febiger, St. Louis, MO, 1982) 2nd ed.
24) L. Kapusta, J. M. Thussen, M. H. M. Cuypers, P. G. M. Peer, and O.

Daniëls: Ultrasoud Med. Biol. 26 (2000) 229.
25) T. Higuchi and M. Kawamata: Digital Shingoshori — MATLAB Taio —

26)
27)
28)
29)
30)
31)
32)
33)
34)

07HF16-9(Digital Signal Processing: With MATLAB Exercises) (Shokodo, Tokyo,
2000) [in Japanese].
B. R. Mahafza: Introduction to Rader Analysis (CRC Press, Boca Raton,
FL, 1998).
A. P. Dempster, N. M. Laird, and D. B. Rubin: J. R. Stat. Soc., Ser. B 39
(1977) 1.
M. Tanter, J. Bercoﬀ, L. Sandrin, and M. Fink: IEEE Trans. Ultrason.
Ferroelectr. Freq. Control 49 (2002) 1363.
H. Hasegawa and H. Kanai: IEEE Trans. Ultrason. Ferroelectr. Freq.
Control 55 (2008) 2626.
C. M. Bishop: Pattern Recognition and Machine Learning (Springer,
Berlin, 2006).
K. Fukunaga: Introduction to Statistical Pattern Recognition (Academic
Press, San Diego, CA, 1990) 2nd ed.
H. Hasegawa and H. Kanai: IEEE Trans. Ultrason. Ferroelectr. Freq.
Control 55 (2008) 1921.
A. M. Katz: Physiology of the Heart (Lippincott Williams & Wilkins,
Philadelphia, PA, 2001).
M. Tanaka, T. Sakamoto, S. Sugawara, H. Nakajima, Y. Katahira, S.
Ohtsuki, and H. Kanai: J. Cardiol. 52 (2008) 86.

# 2011 The Japan Society of Applied Physics