2450

IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control , vol. 56, no. 11, November 2009

Phase-Sensitive Lateral Motion Estimator for
Measurement of Artery-Wall Displacement—
Phantom Study
Hideyuki Hasegawa and Hiroshi Kanai
Abstract—Artery-wall motion due to the pulsation of theheart is often measured to evaluate mechanical properties of
the arterial wall. Such motion is thought to occur only in thearterial radial direction because the main source of the motion
is an increase of blood pressure. However, it has recently beenreported that the artery also moves in the longitudinal direction. Therefore, a 2-D motion estimator is required even when
the artery is scanned in the longitudinal direction because thearterial wall moves both in the radial (axial) and longitudinal (lateral) directions. Methods based on 2-D correlation of
RF echoes are often used to estimate the lateral displacement
together with axial displacement. However, these methods require much interpolation of the RF echo or correlation function to achieve a sufficient resolution in the estimation of displacement. To overcome this problem, Jensen et al. modulated
the ultrasonic field in the lateral direction at a designed spatial
frequency to use the lateral phase for the estimation of lateralmotion. This method, namely, the lateral modulation method,
generates complex signals whose phases change depending on
the lateral motion. Therefore, the lateral displacement can be
estimated with a good resolution without interpolation, although special beamformers are required. The present paper
describes a method that can be applied to ultrasonic echoes
obtained by a conventional beamformer to estimate lateral displacements using the phases of lateral fluctuations of ultrasonicechoes. In the proposed method, complex signals were generated by the Hilbert transform, and the phase shift was estimated by correlation-based estimators. The proposed method
was validated using a cylindrical phantom mimicking an artery.
The error in the lateral displacement estimated by the proposed method was 13.5% of the true displacement of 0.5 mm
with a kernel size used for calculating the correlation functionof 0.6 mm in the lateral direction, which was slightly smaller
than the width at −20 dB of the maximum lateral ultrasonicfield (about 0.8 mm).

I. IntroductionM

ethods for estimating the artery-wall motion due
to the pulsation of the heart have been developed todetect atherosclerotic changes of the artery-wall mechanical properties because they are well known to be alteredby atherosclerosis [1], [2]. In these methods, the radial
motion of the arterial wall, such as changes in diameter[3]–[6] and radial strain [7]–[9], are measured because it
is considered that the source of the artery-wall motion isthe change in internal pressure (blood pressure) and that
there is no longitudinal motion. However, Cinthio et al.
showed that the artery also moves in the longitudinal di-

Manuscript received April 3, 2009; accepted August 2, 2009.
H. Hasegawa and H. Kanai are with the Graduate School of Biomedical Engineering/Graduate School of Engineering, Tohoku University,
Sendai, Japan (e-mail: hasegawa@ecei.tohoku.ac.jp).
Digital Object Identifier 10.1109/TUFFC.2009.1332

rection [10]. Therefore, a 2-D motion estimator is required
to estimate both the axial (radial) and lateral (longitudinal) displacements, even when an artery is scanned in thelongitudinal direction.
de Korte et al. introduced a motion estimator based on
the correlation between RF echoes to measure the 2-D displacement of the arterial wall in the cross-sectional scan[11]. Such a motion estimator based on correlation between RF echoes had previously been developed and thoroughly investigated in the field of tissue elasticity imaging[12]–[15], results indicating that it accurately estimatesthe 2-D displacement. However, much interpolation is required to realize an accurate estimation.
Jensen et al. introduced a method, namely, the lateralmodulation method, in which the ultrasonic field is modulated in the lateral direction at a designed spatial frequency to realize a lateral displacement estimation using
the lateral phase induced by the modulated field [16]. Thismethod generates complex signals whose phases change
depending on the lateral displacement and, therefore, thelateral displacement can be estimated with a good resolution because the phase change can be directly converted
into the lateral displacement. However, this method requires special beamformers that are not available in conventional equipment.
To overcome this problem, Chen et al. recently proposed a method to estimate lateral displacements using
the lateral phases of ultrasonic echoes obtained by conventional beamformers [17]. The accurate estimation of lateral displacements would be very useful, particularly whenit could be done based on conventional beamformers. The
present paper describes a method, which also uses thelateral phases of echoes obtained by conventional beamformers, for estimation of lateral displacements of arterial
walls. In the proposed method, complex signals are generated by the Hilbert transform, and the phase shift due
to the lateral motion is estimated by a correlation-based
estimator. The proposed method was validated using acylindrical phantom mimicking an artery and comparedwith the lateral modulation method.

II. Materials and Methods
A. Difficulties Encountered in Estimation of Lateral
Motion Using Phases of Ultrasonic Echoes
This section describes the fundamental theory of ultrasonic fields to show why it is difficult to use the phase0885–3010/$25.00 © 2009 IEEE

hasegawa and kanai: phase-sensitive lateral motion estimator–phantom study

2451

oscillate at the same spatial frequency, fx0, but with aphase difference of 90 degrees, are produced to create acomplex signal g(x; n), whose phase changes depending on
the lateral displacement of the object. Based on the relationship in (1), g(x; n) can be considered to be the complex
version of s(x; n). In this case, the 2 point spread functions
hr(x) and hi(x) are approximately expressed by cos(2πfx0x)
and −sin(2πfx0x), respectively. Therefore, complex signal
g(x; n) can be obtained based on (1) as follows:g(x; n) ==Fig. 1. Geometry for measurement.

information of ultrasonic echoes for estimation of lateralmotion.
When the ultrasonic field is focused at a depth of interest z using a linear array probe, a point spread function(PSF) h(x) is created as illustrated in Fig. 1, where only
its profile in the lateral direction x at depth z is considered.
By defining the spatial distribution of the amplitude reflection coefficient of an object in the nth frame as r(x; n),the amplitude s(x; n) of an echo at depth z obtained by an
ultrasonic beam focused at the point of interest (x, z) isexpressed as follows:s(x; n) ==¥h(x) × r(x - x; n)dxò-¥

h(x) × r(x - x; n)dx

(1)

= h(x ) * r(x; n),where * denotes convolution.
Let ux(n) be the lateral displacement of the object between the nth and (n + 1)th frames. Reflection coefficient r(x; n + 1) in the (n + 1)th frame is expressed byr(x; n + 1) = r(x − ux(n); n), where it can be assumed that
there is only the lateral motion when the axial motion iscompensated by an axial motion estimator (in this study,
the method proposed in [18] was used). In addition, it
was assumed that there is no distortion in r(x; n) between
the nth and (n + 1)th frames. Under such conditions, theamplitude s(x; n + 1) is given bys(x; n + 1) = h(x ) * r(x; n + 1)
= h(x ) * r(x - u x(n); n).

(2)

In the lateral modulation method proposed by Jensen
et al. [16], 2 point spread functions, hr(x) and hi(x), which¥(3)ò-¥r(x - x; n) × exp(-j2p f x0x)dx.

As can be seen in (3), the complex signal g(x; n) obtainedby the lateral modulation method is the Fourier coefficientR(fx0; n) of reflection coefficient r(x; n) at spatial frequency
fx0. Therefore, g(x; n + 1) in the (n + 1)th frame can be
expressed by g(x; n + 1) = g(x; n) × e -2p f x 0u x (n). Under suchcondition, the lateral displacement ux(n) can be estimatedby the phase shift −2πfx0ux(n) from g(x; n) to g(x; n + 1)
using the conventional correlation technique because thespatial modulation frequency fx0 can be appropriately obtained by designing the point spread functions hr(x) andhi(x).
In the present study, the Hilbert transform was applied
to RF echoes obtained by conventional beamforming to
use the lateral phase. Complex spectrum S(fx; n) of s(x; n)of (1) is expressed by H(fx; n) and R(fx; n) of the point
spread function h(x) and reflection coefficient r(x; n) at
spatial frequency fx as S(fx; n) = H(fx) · R(fx; n). Complex
spectra H(fx) and R(fx; n) are described as follows:

ò-¥
¥

¥

ò-¥r(x - x; n) × (h r(x) + jh i(x))dx¥

H (f x ) =ò-¥h(x) × e -j2pf xdx,

R(f x; n) =ò-¥r(x; n) × e -j2pf xdx.

¥

(4)

x

x

(5)

Analytic signal y(x; n) of s(x; n) is obtained by the inverse Fourier transform of S(fx; n) in the range of positive
spatial frequencies as follows; see Appendix for derivationof (6):¥y(x; n) =ò0 H (f x ) × R(f x; n) × e j2pf x df x=ò-¥H ¢(f x ) × R¢(f x; n) × e j2pf x df x

¥

x

x

(6)

= h ¢(x ) * r ¢(x; n),whereìï H (f x ) (f x ³ 0),
H ¢(f x ) = ïíïï 0
(f x < 0),î

(7)ïì R(f x; n) (f x ³ 0),
R¢(f x; n) = ïíïï 0
(f x < 0),î

(8)IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control , vol. 56, no. 11, November 20092452

h ¢(x ) =¥¥

ò-¥H ¢(f x )e j2pf x df x Û H ¢(f x ) = ò-¥h¢(x)e -2pf xdx,
x

x

(9)
¥

ò0
¥
R¢(f x; n) = ò r ¢(x; n)e -2p f xdxx.
-¥
r ¢(x; n) =R¢(f x; n)e j2p f xx df x Û(10)

B. Principle of Lateral Displacement Estimation Using
the Lateral Phase
Let us define the complex correlation function
γ(Δm; m; n) between y(m; n) of (12) and y(m; n + 1) of(13) at lateral lag Δm · Δx asxMc

In the integration of (6), in a strict sense, the direct current component (at fx = 0) should be multiplied by 0.5.
However, in this study, the direct current component in
the measured signal s(x; n) was removed before applying
the Fourier transform to s(x; n). In this case, the Hilberttransform can be expressed by (6).
Similarly, the analytic signal y(x; n + 1) of s(x; n + 1) inthe (n + 1)th frame is expressed as follows; see Appendix
for derivation of (11):y(x; n + 1) =¥ò-¥H ¢(f x ) × R¢(f x; n) × D(f x )e -j2pf u (n)xx× e j 2p f xx df x
= h ¢(x ) * r ¢(x; n) *(11)¥

ò0 e j2pf (x -u (n))df x.
xx

In actual measurements, s(x; n) is sampled at the interval of scan lines Δx, and the sampled version of s(x; n)
is denoted by s(mΔx; n) ≡ s(m; n) (m = −M/2, −M/2 +1, …, −2, −1, 0, 1, 2, …, M/2), where M + 1 is the number of scan lines (lateral length of the scanned region
L = M · Δx). In such a discrete system, (6) and (11) aredenoted in the digital system as follows:y(m; n) = h ¢(m) * r ¢(m; n),(12)y(m; n + 1) = h ¢(m) * r ¢(m; n)*

1 - e j(2p/M )(M /2 +1)(m -(u x (n)/ Dx ))1 - e j(2p/M )(m -(u x (n)/ Dx ))

,

(13)

where (12) and (13) are obtained from (6) and (11), respectively, by replacing the integration, spatial frequency
fx, and lateral spatial position x by the summation, k/(MΔx) (discrete spatial frequency), and mΔx (discrete
spatial position), respectively. As shown in (12) and (13),the phase shift from y(m; n) to y(m; n + 1) actually depends on the lateral displacement ux(n) of the object between the nth and (n + 1)th frames. The phase shift ofcomplex signal g(x; n) of (3) obtained by the lateral modulation method is simply related to lateral displacement2ux(n) such as g *(x; n) × g(x; n + 1) = g(x; n) × e j2p f x 0u x (n),where * denotes complex conjugate. However, it is difficult
to relate the phase shift from y(m; n) to y(m; n + 1) to thelateral displacement ux(n) because h(x) cannot be assumed
to be a sinusoidal wave fluctuating at a known spatial
frequency of fx0 when a conventional beamformer is used.
This is a major difficulty for utilization of the lateral phasewith conventional beamformers. In the present study, a
new method was introduced to overcome this problem.

g(Dm; m; n) =

å y *(m + k; n) × y(m + Dm + k; n + 1),

k =-M c(14)where Mc determines the number of sampled points used
for calculating the correlation function.
The phase shift Δθ(m; n) from y(m; n) to y(m; n + 1)
induced by lateral displacement ux(n) between the nthand (n + 1)th frames can be obtained by setting laterallag Δm at 0 (conventional correlation technique [19]) as
follows:
Dq(m; n) = Ðg(Dm; m; n) Dm ® 0= Ðg(0; m; n).

(15)

In this study, as illustrated in Fig. 2, it was assumed
that there is a linear relationship between the lateraldisplacement ux(n) and the change in the lateral phaseΔθ(m; n), as expressed byu x(n) = a(n) × Dq(m; n)= a(n) × Ðg(0; m; n),(16)

where a(n) is a constant (corresponding to the slope of thelinear relationship), which linearly relates the phase shift
to the lateral displacement.
When Δm is set at 1, it can be considered thaty(m; n + 1) of (13) is artificially displaced by Δx (= an interval of scan lines) relative to y(m; n) of (12). Therefore,
the following relationship holds:
u x(n) + Dx = a(n) × Ðg(1; m; n).

(17)

By solving simultaneous equations consisting of (16)and (17), slope a(n) and lateral displacement ux(n) are
estimated as follows:â (Dm 1,Dm 2)(n) Dm ® 0, Dm ®1 =1

2Dx
Ðg(1; m; n) - Ðg(0; m; n)Dx,Ð ( g (0; m; n) × g(1; m; n) )(18)
Dxû x,(Dm 1,Dm 2)(n) Dm ® 0, Dm ®1 =1
2
Ð ( g *(0; m; n) × g(1; m; n) )=*× Ðg(0; m; n),(19)

where â (Dm 1,Dm 2)(n) and û x,(Dm 1,Dm 2)(n) are the slope andthe lateral displacement, respectively, which are estimatedhasegawa and kanai: phase-sensitive lateral motion estimator–phantom study

2453

the mean squared difference α(n) between the actual andmodel displacements is expressed as follows:
Dm 2a(n) =

å{ DmDx - a(n)Ðg¢(Dm; m; n) }2

Dm = Dm 1
Dm 2

=

å

{ Dm 2Dx 2 - 2DmDx × a(n)Ðg¢(Dm; m; n)

Dm = Dm 1+ a 2(n)Ðg ¢ 2(Dm; m; n) } .

(22)

To determine aˆ (Dm 1,Dm 2)(n), which minimizes the meansquared difference α(n), the partial derivative of (22) withrespect to a(n) is set to be zero:
Dm 2¶a= -2DxDmÐg ¢(Dm; m; n)¶aDm = Dm

å

1

Dm 2

+ 2a(n)å Ðg¢ 2(Dm; m; n)(23)

Dm = Dm 1= 0.

By solving (23), â (Dm 1,Dm 2)(n) is obtained as follows:
Dmâ (Dm 1,Dm 2)(n) =2
Dx åDm =DmÐg ¢(Dm; m; n)
Dm 1

Dm 2

åDm = Dm Ðg¢ 2(Dm; m; n)1

Dm

Fig. 2. Assumed relationship between phase shift and lateral displacement. Complex correlation function γ(Δm; m; n) is shown by γ(Δm).

by correlation functions at lateral lags from Δm1 toΔm2.
Although more computation is required, other correlation functions at different lateral lags can be used for
displacement estimation. As in (17), the following relationship holds:
u x(n) + Dm Dx = a(n)Ðg(Dm; m; n),

(20)

=

Dm 2
Dm = Dm 1å{ Ðg(Dm; m; n) - Ðg(0; m; n) }(21)

where Ðg ¢(Dm; m; n) = Ðg(Dm; m; n) - Ðg(0; m; n).
By considering the left and right sides of (21) to be
the actual and model lateral displacements, respectively,

2

.

(24)

By substituting (24) into (16), lateral displacement ux(n)is estimated as follows:û x,(Dm 1,Dm 2)(n)Dm=

2Dx åDm =Dm { Ðg(Dm; m; n) - Ðg(0; m; n) }
Dm 1

Dm 2

åDm = Dm { Ðg(Dm; m; n) - Ðg(0; m; n) }2

1× Ðg(0; m; n).

where it should be noted that ΔmΔx is an artificial displacement.
Using more than 2 correlation functions, lateral displacement ux(n) of a target can be estimated by the leastsquares method. To do that, let us consider the relationship that is obtained by subtracting (16) from (20) as
follows:
Dm Dx = a(n) { Ðg(Dm; m; n) - Ðg(0; m; n) }º a(n)Ðg ¢(Dm; m; n),2
Dx åDm =Dm { Ðg(Dm; m; n) - Ðg(0; m; n) }
Dm 1

(25)

In a subsequent section describing experiments using aphantom, accuracies in estimation of lateral displacements
achieved by the 2 estimators, i. e., the computationallyefficient version given by (19) and the estimator given by(25) consisting of 3 correlation functions at lags {Δm} of−1, 0, and 1 (Δm1 = −1, Δm2 = 1), are compared.
C. Experimental System
In this study, a cylindrical phantom (inner diameter:
8 mm; external diameter: 10 mm) made from silicone rubber (elastic modulus: 750 kPa) containing 5% carbon powder (by weight) was measured in the experimental system

2454

IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control , vol. 56, no. 11, November 2009

Fig. 3. Schematic of measurement system.

illustrated in Fig. 3. The radial motion ( = axial motion) of
the phantom was induced by changing the internal pressureusing a flow pump (pulse pressure: about 60 mmHg, theoretical resulting radial strain: about 4%). The longitudinalmotion ( = lateral motion) was simulated by moving a ultrasonic probe using an automatic stage. The maximum lateraldisplacements were controlled to be 0.1, 0.25, and 0.5 mmby the automatic stage. The stage was triggered by a signal
from the flow pump, which shows the beginning of ejection.
In ultrasonic measurements for the method proposedin this study, RF echoes from the phantom were acquired
at a frame rate of 286 Hz with a 10-MHz linear arrayprobe (UST-5545, Aloka, Tokyo, Japan) equipped withconventional ultrasonic diagnostic equipment (SSD-6500,
Aloka, Tokyo, Japan). The phantom was scanned in thelongitudinal direction at intervals Δx of 0.15 mm (46 scanlines), and RF echoes were sampled at 40 MHz at a 16-bit
resolution.
For the lateral modulation method, a scanner (its front
end is same as that of α-10 [Aloka, Tokyo, Japan]), which
was modified so that RF echoes received by each arrayelement could be acquired (frame rate: 289 Hz), was employed together with a 10-MHz linear array probe (UST5545, Aloka, Tokyo, Japan) [20]. With this system, planewaves were transmitted and the receive beamforming was
performed with the apodization and delay factors shownin Fig. 4 [16]. Figs. 5(a) and 5(b) show a beamformed
RF echo from a fine wire (diameter: 16 μm) placed about
20 mm away from the ultrasonic probe, and Fig. 5(c)shows the spectrum of the echo signal shown in Fig. 5(a)obtained by the 2-D Fourier transform. There were 72
scan lines at lateral intervals Δx of 0.2 mm, and the lateral modulation frequency fx0 was set at 0.89 mm−1.
III. Results of Experiments Using a Phantom
A. Estimation of Lateral Displacements by the
Proposed Method
Fig. 6 shows the procedure to obtain analytic signals{y(m; n)}. Fig. 6(a) shows RF echoes from the phantom
that was scanned in the longitudinal direction using the

Fig. 4. (a) Amplitude apodization and (b) time delay (for in-phase [solidline] and quadrature [dashed line] beamformers) values for creating thelateral modulation during receive beamforming.

linear array probe with a conventional beamformer. Envelopes {s′(m; n)} of the RF echoes were detected as shownin Fig. 6(c). The direct current component contained in
the envelope signal s′(m; n) has no phase information and,
therefore, it was removed. In this study, the envelope
signal without bias was used as echo amplitude s(m; n),which is given bys(m; n) = s ¢(m; n) - E m [ s ¢(m; n) ],(26)

where E[·] denotes the averaging operation with respect
to lateral position m · Δx. The Hilbert transform witha Tukey window shown in Fig. 6(b) was then applied tos(m; n) to obtain the analytic signal y(m; n). The real andimaginary parts of analytic signal y(m; n) were obtained asshown in Figs. 6(e) and 6(f).
The method proposed in Section II-B was applied toanalytic signal y(m; n) to estimate lateral displacementux(n). Twenty points of interest were assigned in the posterior wall at axial intervals of 50 μm along each scan line,
and lateral displacements {ux(n)} of these points were estimated. In Figs. 7(1) and 7(2), plots and vertical bars
show means and standard deviations of the maximum lateral displacements {û x (01),max(n)} and {û x (-11),max(n)} alongeach scan line estimated by the estimators given by (19)and (25), respectively, with 4 different sizes of kernels used

hasegawa and kanai: phase-sensitive lateral motion estimator–phantom study

2455

respectively. As shown in Fig. 7, standard deviations werereduced by the estimator consisting of 3 correlation functions given by (25).
Fig. 8 shows rms errors of the estimated lateral displacements obtained by the estimator given by (19) fromthe actual displacement, where the errors were evaluatedfrom all the estimates except for those obtained at 10
scan lines at each edge of the scanned region because theestimates at these scan lines were influenced by the shapeof the Tukey window (not flat). As shown in Fig. 8, the errors were reduced by increasing the size of the correlation
kernel, and similar errors resulted from the kernel sizes
greater or equal to 2.5 mm.
B. Comparison with Results Obtained by the Lateral
Modulation Method

Fig. 5. Radio frequency echo from a fine wire obtained by the lateralmodulation method. (a) Radio frequency echo Re{g(x; n)} obtained by inphase receive beamforming, where Re{·} means the real part. (b) Lateralprofiles at a depth indicated by the arrow in (a) obtained by in-phaseand quadrature receive beamforming. (c) Magnitude of spectrum of theecho signal Re{g(x; n)} shown in (a) obtained by the 2-D Fourier transform.

for calculating correlation function γ(Δm; m; n) defined by(14). The dashed lines in Fig. 7 show the actual assigned
maximum lateral displacements. In calculation of the correlation function, 2-D kernels were used, the axial size ofthe kernel being fixed to be the optimum value (0.5 μs)
determined in [18], and the lateral size of a kernel beingchanged. In Figs. 7(b), 7(c), 7(d), and 7(e), the lateralsizes of kernels were set at 0.6 mm (Mc = 2), 1.2 mm(Mc = 4) (see the multimedia file for an example under
this condition), 2.4 mm (Mc = 8), and 3.6 mm (Mc = 12),Fig. 9(a) shows RF echoes from the phantom obtained by the lateral modulation method (in-phasebeamformer), and Fig. 9(b) shows the lateral profilesof RF echoes at a depth indicated by the arrow in
Fig. 9(a) obtained by in-phase and quadrature beamformers. As in the displacement estimation by the
proposed method, changes in the phases of complex
signals {g(x; n)} obtained by in-phase and quadrature
beamformers were estimated by the correlation technique described in [16] to obtain lateral displacements{ux(n)} at 20 points of interest assigned along eachscan line at axial intervals of 50 μm.
As in Section III-A, the use of different lateral sizes ofkernels used for calculation of correlation functions was
examined. Lateral displacements estimated by applying
the motion estimator in [16] to complex signals {g(x; n)}obtained by the lateral modulation method are shown in
Fig. 10(1). In Figs. 10(1-a), (1-b), (1-c), and (1-d), lateraldisplacements {ux(n)} were estimated using kernel sizes of
0.6 mm (Mc = 2), 1.2 mm (Mc = 4), 2.4 mm (Mc = 8),
and 3.6 mm (Mc = 12) (see the multimedia file for an example under this condition), respectively. As in Fig. 10(1),lateral displacements obtained by applying the proposed
motion estimator to {g(x; n)} are shown in Fig. 10(2) for
different correlation kernel sizes. There were no improvements by using the proposed motion estimator together
with the lateral modulation method. As shown in Fig.
11, rms errors of the lateral displacements obtained by
the motion estimator in [16] were calculated from all theestimates shown in Fig. 10(1).
IV. Discussion
In this study, lateral displacements were estimated using the phases of complex signals generated by the Hilberttransform applied to ultrasonic echoes obtained by a conventional beamformer. The phase shift due to the lateraldisplacement was estimated by the conventional correlation technique, and a larger correlation kernel size in the2456

IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control , vol. 56, no. 11, November 2009

Fig. 6. (a) Radio frequency echoes from the phantom. (b) Tukey window used for frequency analysis of envelope signals of RF echoes. (c) Envelope
signals of RF echoes. (d) Envelope signals after bias removal. (e) Imaginary and (f) real parts of analytic signals {y(x; n)}.

lateral direction was found to achieve better accuracy. Errors in the estimated lateral displacements were reduced
by increasing the kernel size, and similar errors were obtained with kernel sizes greater or equal to 2.5 mm, whichroughly corresponds to 3 times the width at −20 dB ofthe maximum lateral ultrasonic field (about 0.8 mm). In aprevious study, the estimation of axial displacements, theoptimal kernel size in the axial direction corresponded toa pulse duration defined by the width at −20 dB of theenvelope of an ultrasonic pulse (about 0.4 mm) [18]. In theaxial displacement estimation, a smaller kernel relative tothe point spread function (0.4 mm) yielded good estimates
because there were 2 to 3 oscillations during a pulse duration. On the other hand, there was no oscillation in thelateral direction (or more exactly, one oscillation because
there was a peak in the lateral profile of the ultrasonicfield). Therefore, a kernel size of about 3 times the point
spread function was required to minimize the error.
To estimate the artery-wall motion using the phasesof ultrasonic echoes, it is necessary to avoid the aliasing
effect. Therefore, the frame rate fFR should be kept as
high as possible. In general, motion of the arterial wall in
the radial direction is larger than that in the longitudinal
direction, and the oscillation frequency in the axial direc-

tion (= radial direction) is higher than that in the lateral
direction (= longitudinal direction). Therefore, the aliasing limit for the radial motion (= axial motion) should beconsidered.
Basically, using the change in the axial phase Δϕ, theaxial velocity va is estimated as follows [19], [21]:
va =

c 0Df
f FR,
4p f 0(27)

where c0 and f0 are the speed of sound and the center
frequency of ultrasound, respectively. In addition, a typical maximal velocity of the carotid arterial wall in the radial direction is about 10 mm/s. The frame rate required
for the measurement of the axial motion at 10 mm/s isobtained by substituting va = 10 mm/s, Δϕ = π rad,c0 = 1500 m/s, and f0 = 10 MHz into (27) as follows: f−3
6
FR = 4 × (10 × 10 ) × (10 × 10 )/1500 ≈ 267 Hz. To
achieve a frame rate higher than 267 Hz at a fixed pulse
repetition frequency of 13156 Hz, the number of scan lines
was reduced to 46 (this is the lowest available number ofscan lines of the employed ultrasonic equipment).
As can be seen in Fig. 8, there is a trade-off between
the accuracy and the kernel size (spatial resolution).
When we want to estimate the global motion of an ob-hasegawa and kanai: phase-sensitive lateral motion estimator–phantom study

2457

Fig. 7. Means and standard deviation of lateral displacements obtained by the estimators given by (1) Eq. (19) and (2) Eq. (25). (a) Tukey windowused in the Hilbert transform. Lateral displacements estimated by the proposed method with correlation kernel sizes of (b) 0.6 mm, (c) 1.2 mm, (d)2.4 mm, and (e) 3.6 mm. The 3 horizontal dashed lines in each figure (b)–(e) show the actual assigned displacements.

ject, a large kernel can be used for higher accuracy. On
the other hand, it is necessary to use a smaller kernel
to obtain the spatial distribution of lateral displacement,
although the accuracy will be degraded. Therefore, themethod should be optimized depending on the purposes,
for example, a greater number of correlation functions inan estimator would reduce standard deviations, as shownin Fig. 7, at the expense of computational efficiency. In addition, further improvements would be required to achieve
better accuracy with a smaller kernel. In this study, alinear relationship between the lateral displacement andthe change in the phase of the complex signal obtainedby the Hilbert transform was assumed. Fig. 12(a) shows a
B-mode image of a fine wire (same as that in Fig. 5), and2458

IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control , vol. 56, no. 11, November 2009

Fig. 8. Root mean square errors of lateral displacements {ux,(0,1)(n)} obtained by the proposed estimator given by (19) are plotted as a functionof correlation kernel size in the lateral direction. The 3 curves show theerrors in cases of actual lateral displacements of 0.1, 0.25, and 0.5 mm.

Fig. 12(b) shows a lateral profile of envelopes of RF echoes
at a depth indicated by the arrow in Fig. 12(a). By applying the Fourier transform to the lateral profile shown in
Fig. 12(b), a power spectrum, which is shown by the solidline in Fig. 12(c), was obtained. As shown by the solidline in Fig. 12(c), in general, the direct current component
is largest (central spatial frequency is zero). Therefore,
it is difficult to use the lateral phase with conventional
beamformers. The dashed line in Fig. 12(c) shows a power
spectrum obtained by applying the Fourier transform to
the lateral profile after removing the direct current component based on (26). In this case, the power spectrum
is largest at a certain spatial frequency (≠ 0), and it canbe considered that the lateral profile is fluctuating at thecentral spatial frequency (≠ 0), as in the lateral modulation method. Therefore, a linear relationship was assumed
in the present paper. In addition, correlation functions ina range from −1 to 1, which corresponded to a range from−150 μm to 150 μm in the lateral direction were used for
estimation of lateral displacements. This region is muchsmaller than the size of the width at −20 dB of the maximum lateral ultrasonic field. Therefore, the assumption of
the linear relationship between the change in the lateralphase and the lateral displacement in such small range
was considered to be appropriate, even when the overall relationship was not perfectly linear. To use a greaternumber of correlation functions, identification of a better
function describing this relationship would be required.
In the present study, although it was assumed that there
was no deformation, the proposed method could estimate
the lateral displacements of a phantom under a specificdegree of deformation. However, a method, in which thedistortion of r(x; n) due to deformation of an object istaken into account, should be developed to improve theaccuracy, and such method should be validated under the
existence of various degrees of deformation.
In the present study, the results obtained by the proposed method were compared with those obtained by

Fig. 9. Radio frequency echo from the phantom obtained by the lateralmodulation method. (a) Radio frequency echo Re{g(x; n)} obtained byin-phase receive beamforming. (b) Lateral profiles of {g(x; n)} at a depth
indicated by the arrow in (a) obtained by in-phase and quadrature receive beamforming.

the lateral modulation method. The results obtained by
the lateral modulation method with the motion estimator described in [16] were worse than those obtained by
the proposed method, and the combination of the lateralmodulation method with the motion estimator proposedin the present paper did not achieve any improvements. It
could be considered from these results that some optimization might be necessary to obtain correctly modulated
fields. Therefore, it cannot be concluded at present thatthe proposed method is better than the lateral modulation method. Nevertheless, the proposed method would be
useful because it can be applied to the ultrasonic echoes
obtained by conventional beamformers.
V. Conclusions
In this study, a method was developed to estimate lateral displacements using the lateral phase, which can beapplied to ultrasonic echoes obtained by a conventional
beamformer. In the proposed method, complex signals
were generated by the Hilbert transform, and the phaseshift due to the lateral displacement was estimated by correlation-based estimators. The proposed method was vali-hasegawa and kanai: phase-sensitive lateral motion estimator–phantom study

2459

Fig. 10. Means and standard deviations of lateral displacements obtained by the lateral modulation method with (1) motion estimator given in [16]and (2) that given by (19). Means and standard deviations obtained with correlation kernel sizes of (a) 0.6 mm, (b) 1.2 mm, (c) 2.4 mm, and (d)3.6 mm in the lateral direction. The 3 horizontal dashed lines in each figure show the actual assigned displacements.

dated using a cylindrical phantom mimicking an artery.
As a result, the lateral displacements could be measured
with an error of 13.5% of the true displacement of 0.5 mm,
and the proposed method would be useful because it canbe applied to ultrasonic echoes obtained by conventional
beamformers.
Appendix
Hilbert Transform of Lateral
Amplitude Profile
As expressed by (6), analytic signal y(x; n) of lateralamplitude profile s(x; n) is obtained by the inverse Fouriertransform of S(fx; n) in the range of positive spatial frequencies as follows:¥ò0 H (f x ) × R(f x; n) × e j2pf x df x¥
= ò H ¢(f x ) × R¢(f x; n) × e j2p f x df x,
-¥

(28)

ìï H (f x ) (f x ³ 0),
H ¢(f x ) = ïíïï 0
(f x < 0),î

(29)ïì R(f x; n) (f x ³ 0),
R¢(f x; n) = ïíïï 0
(f x < 0).
î(30)y(x; n) =

x

x

where

IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control , vol. 56, no. 11, November 2009

2460

Fig. 11. Root mean square errors of lateral displacements estimated by
the lateral modulation method with the motion estimator in [16] areplotted as a function of correlation kernel size in the lateral direction.
The 3 curves show the errors in cases of actual lateral displacements of
0.1, 0.25, and 0.5 mm.

Eq. (28) can be expressed as follows:¥ ìïïí ¥ h ¢(x )e -j2p f x x 1dx ¥ r ¢(x ; n)e -j2p f x x 2dx ïüý
1
1ò
2
2
ò
ïïþ
-¥ ï
-¥
ïî -¥
j 2p f x x´e
df x ,(31)y(x; n) = òwhereh ¢(x ) = ò¥H ¢(f x )e j2p f xx df x Û H ¢(f x ) = ò-¥¥h ¢(x)e -j2p f x xdx,-¥¥

ò0 R¢(f x; n)e j2pf x df x Û¥R¢(f x; n) = ò r ¢(x; n)e -j2p f xdx.
-¥
r ¢(x; n) =Fig. 12. (a) B-mode image of a fine wire obtained by a conventional
beamformer. (b) Lateral profile of envelopes of received ultrasonic RF
echoes at a depth indicated by the arrow in (a). (c) Power spectra oflateral profiles of envelopes. The solid line shows power spectrum of theenvelope shown in (b). The dashed line shows power spectrum of theenvelope that was obtained by removing the direct current component
from (b) based on (26).

(32)

x

(33)

x

Analytic signal y(x; n + 1) in the (n + 1)th frame isexpressed as follows:

By defining ξ = ξ1 + ξ2, (31) is modified as shown in (34),
see next page.
As can be seen in (34), the integration with respect to ξcorresponds to the Fourier transform of h′(x)*r′(x; n) and,
therefore, the integration with respect to fx corresponds
to the inverse Fourier transform of F[h′(x)*r′(x; n)], whereF[∙]means the Fourier transform. Thus, y(x; n) is expressedas y(x; n) = h′(x)*r′(x; n), which corresponds to (6).
Similarly, analytic signal y(x; n + 1) in the (n + 1)th
frame is obtained as described below. Spatial distributionof reflection coefficient r(x; n + 1) in the (n + 1)th frame
is expressed as r(x; n + 1) = r(x − ux(n); n). Therefore,
complex spectrum R(fx; n + 1) is expressed as
R(f x; n + 1) = R(f x; n) × e -j2p f xu x (n). Then, let us expressR′(fx; n + 1) as R′(fx; n + 1) = R′(fx; n) · P(fx; n), wherey(x; n + 1) =¥¥ò-¥H ¢(f x ) × R¢(f x; n + 1) × e j2pf x df xxò-¥H ¢(f x ) × R¢(f x; n) × P(f x; n) × e j2pf x df x
¥ ì
¥
ï ¥
= ò í ò h ¢(x 1)e -j2p f x dx 1ò r ¢(x 2; n)e -j2p f x dx 2
ïï -¥
-¥ î
-¥
=

x

x 2

x 1

´

¥

ïü

ò-¥p(x 3; n)e -j2pf x dx 3 ïýïïþe j2pf x df x,
x 3

(35)

(36)

where
¥

ò-¥P(f x; n)e j2pf x df x¥
= ò e -j2p f u (n)e j2p f x df x.
0p(x; n) =xxïìe -j2p f xu x (n) (f x ³ 0),P(f x; n) = ïíïï 0
(f x < 0),
î

x

x

x

By defining ξ′ = ξ1 + ξ2, (36) is modified as follows:(37)

hasegawa and kanai: phase-sensitive lateral motion estimator–phantom studyy(x; n) =

=

¥ ï
ì

¥

¥ ì
ï

¥ ì
ï

¥ ì
ï

¥

¥

x 1

üï

¥

1

üï

x

=

ò-¥ ïíîïï ò-¥ {h¢(x) * r ¢(x; n)}e -j2pf xdx ïýþïïe j2pf x df x.
¥ ï
ì

¥

¥ ï
ì

¥ ï
ì

¥ ì
ï

¥

x

x

üï

(34)

x

x

¥

ïü

ò-¥ ïíïïî ò-¥h¢(x¢) * r ¢(x¢; n)e -j2pf x dx¢ò-¥p(x - x¢; n)e -j2pf (x -x )dx ïýïïþe j2pf x df x
x

¢

x

ïü

¥

¢

ïü

x

=

ò-¥ ïíïïî ò-¥ ïíïïî ò-¥ {h¢(x¢) * r ¢(x¢; n)} p(x - x¢; n)dx¢ ïýïïþe -j2pf xdx ïïïþýe j2pf x df x=ò-¥ íïîïï ò-¥ {h¢(x) * r ¢(x; n) * p(x; n)}e -j2pf xdx ïýþïïe j2pf x df x .

¥

x

x 1

¥

¢

´

ò-¥r ¢(x¢ - x 1; n)e -j2pf (x -x )dx¢´

ò-¥p(x 3; n)e -j2pf x dx 3 ïýþïïe j2pf x df x (38)

x

¥

x 3

1

üï

x

ò-¥ { ò-¥h¢(x¢) * r ¢(x¢; n)e -j2pf x¢dx¢´xò-¥ ïíîïï ò-¥ ïíîïï ò-¥h¢(x 1) × r ¢(x - x 1; n)dx 1 ïýþïïe -j2pf xdx ïýþïïe j2pf x df xò-¥ ïíïïî ò-¥h¢(x 1)e -j2pf x dx 1¥üïò-¥ ïíïïî ò-¥h¢(x 1)e -j2pf x dx 1ò-¥r ¢(x - x 1; n)e -j2pf (x -x )dx ïýïïþe j2pf x df x=y(x; n + 1) =y(x; n + 1) =

¥ ì
ï

2461

¥

x

ò-¥p(x 3; n)e -j2pf x dx 3 }e¥x 3

j 2p f x x

df x .

Again, by defining ξ = ξ′ + ξ3, (38) is modified asshown in (39), see above.
As in (34), the integration with respect to ξ corresponds
to the Fourier transform of h′(x)*r′(x; n)*p(x; n) and, therefore, the integration with respect to fx corresponds to theinverse Fourier transform of F[h′(x)*r′(x; n)*p(x; n)]. Thus,y(x; n + 1) is expressed as follows:y(x; n + 1) = h ¢(x ) * r ¢(x; n) * p(x; n)= h ¢(x ) * r ¢(x; n) *¥ò0e j2p f x (x -u x (n))df x.

(40)

Eq. (40) corresponds to (11).
Acknowledgments
We are grateful to Panasonic Co., Ltd., and Panasonic
Shikoku Co., Ltd., for providing the cylindrical phantom
used in this study.
References[1]	R. T. Lee, A. J. Grodzinsky, E. H. Frank, R. D. Kamm, and F.
J. Schoen, “Structure-dependent dynamic mechanical behavior of

üï

x

x

x

(39)

fibrous caps from human atherosclerotic plaques,” Circulation, vol.
83, no. 5, pp. 1764–1770, 1991.
[2] H. M. Loree, A. J. Grodzinsky, S. Y. Park, L. J. Gibson, and R. T.
Lee, “Static circumferential tangential modulus of human atherosclerotic tissue,” J. Biomech., vol. 27, no. 2, pp. 195–204, 1994.
[3] J. M. Meinders, P. J. Brands, J. M. Willigers, L. Kornet, and A. P.
G. Hoeks, “Assessment of the spatial homogeneity of artery dimension parameters with high frame rate 2-D B-mode,” Ultrasound Med.
Biol., vol. 27, no. 6, pp. 785–794, 2001.
[4]	A. P. G. Hoeks, C. J. Ruissen, P. Hick, and R. S. Reneman, “Transcutaneous detection of relative changes in artery diameter,” Ultrasound Med. Biol., vol. 11, no. 1, pp. 51–59, 1985.
[5] T. Länne, H. Stale, H. Bengtsson, D. Gustafsson, D. Bergqvist, B.
Sonesson, H. Lecerof, and P. Dahl, “Noninvasive measurement of
diameter changes in the distal abdominal aorta in man,” UltrasoundMed. Biol., vol. 18, no. 5, pp. 451–457, 1992.
[6] P. J. Brands, A. P. G. Hoeks, M. C. M. Rutten, and R. S. Reneman,“A noninvasive method to estimate arterial impedance by means ofassessment of local diameter change and the local center-line bloodflow velocity using ultrasound,” Ultrasound Med. Biol., vol. 22, no.
7, pp. 895–905, 1996.
[7]	C. L. de Korte, E. I. Céspedes, A. F. W. van der Steen, and C. T.
Lanée, “Intravascular elasticity imaging using ultrasound: Feasibility studies in phantoms,” Ultrasound Med. Biol., vol. 23, no. 5, pp.
735–746, 1997.
[8]	O. Bonnefous, “Blood flow and tissue motion with ultrasound for
vascular applications,” Comptes Rendus Acad. Sci. IV, vol. 2, pp.
1161–1178, 2001.
[9] H. Hasegawa, H. Kanai, Y. Koiwa, and N. Chubachi, “Noninvasive
evaluation of Poisson’s ratio of arterial wall using ultrasound,” Electron. Lett., vol. 33, no. 4, pp. 340–342, 1997.
[10]	M. Cinthio, A. R. Ahlgren, T. Jansson, A. Eriksson, H. W. Persson, and K. Lindström, “Evaluation of an ultrasonic echo-tracking
method for measurements of arterial wall movements in two dimensions,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 52, no.
8, pp. 1300–1311, Aug. 2005.
[11] H. Ribbers, S. Holewijn, J. D. Blankensteijn, and C. L. de Korte,“Non-invasive two dimensional elastography of the carotid artery,”in Proc. 2005 IEEE Ultrasonics Symp., 2005, pp. 1113–1116.
[12] J. Ophir, I. Céspedes, H. Ponnekanti, Y. Yazdi, and X. Li, “Elastography: A quantitative method for imaging the elasticity of biologicaltissues,” Ultrason. Imaging, vol. 13, pp. 111–134, 1991.
[13] I. Céspedes and J. Ophir, “Reduction of image noise in elastography,” Ultrason. Imaging, vol. 15, pp. 89–102, 1993.
[14] E. E. Konofagou, J. Ophir, T. Varghese, and F. Kallel, “A newelastographic method for estimation and imaging of lateral displacements, lateral strains, corrected axial strain and Poisson’s ratios intissues,” Ultrasound Med. Biol., vol. 24, no. 8, pp. 1183–1199, 1998.

2462

IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control , vol. 56, no. 11, November 2009[15] J. Ophir, B. Garra, F. Kallel, E. Konofagou, T. Krouskop, R.
Righetti, and T. Varghese, “Elastographic imaging,” UltrasoundMed. Biol., vol. 26, pp. S23–S29, 2000.
[16] J. A. Jensen and P. Munk, “A new method for estimation of velocity
vectors,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 45,
no. 3, pp. 837–851, May 1998.
[17] X. Chen, M. J. Zohdy, S. Y. Emelianov, and M. O’Donnell, “Lateral speckle tracking using synthetic lateral phase,” IEEE Trans.
Ultrason. Ferroelectr. Freq. Control, vol. 51, no. 5, pp. 540–550, May
2004.
[18] H. Hasegawa and H. Kanai, “Reduction of influence of variation incenter frequencies of RF echoes on estimation of artery-wall strain,”IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 55, no. 9, pp.
1921–1934, Sep. 2008.
[19]	C. Kasai, K. Namekawa, A. Koyano, and R. Omoto, “Real-time twodimensional blood flow imaging using an autocorrelation technique,”IEEE Trans. Sonics Ultrason., vol. SU-32, no. 3, pp. 458–464, 1985.
[20] H. Hasegawa and H. Kanai, “Simultaneous imaging of artery-wallstrain and blood flow by high frame rate acquisition of RF signals,”IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 55, no. 12, pp.
2626–2639, Dec. 2008.
[21] H. Kanai, M. Sato, Y. Koiwa, and N. Chubachi, “Transcutaneous
measurement and spectrum analysis of heart wall vibrations,” IEEE
Trans. Ultrason. Ferroelectr. Freq. Control, vol. 43, no. 5, pp. 791–
810, Sep. 1996.

Hideyuki Hasegawa was born in Oyama, Japan, in 1973. He received the B.E. degree fromTohoku University, Sendai, Japan, in 1996. He received the Ph.D. degree from Tohoku University
in 2001. He is presently an associate professor in
the Graduate School of Biomedical Engineering,Tohoku University. His main research interest is
medical ultrasound, especially diagnosis of atherosclerosis based on measurements of motion andmechanical properties of the arterial wall. Dr. Hasegawa is a member of
the IEEE, the Acoustical Society of Japan, the Japan Society of Ultrasonics in Medicine, and the Institute of Electronics, Information andCommunication Engineers.

Hiroshi Kanai was born in Matsumoto, Japan,
on November 29, 1958. He received a B.E. degreefrom Tohoku University, Sendai, Japan, in 1981,
and M.E. and Ph.D. degrees, also from Tohoku
University, in 1983 and in 1986, respectively, bothin electrical engineering. From 1986 to 1988, he
was with the Education Center for Information
Processing, Tohoku University, as a research associate. From 1990 to 1992, he was a lecturer inthe Department of Electrical Engineering, Facultyof Engineering, Tohoku University. From 1992 to
2001, he was an associate professor in the Department of Electrical Engineering, Faculty of Engineering, Tohoku University. Since 2001, he has
been a professor in the Department of Electronic Engineering, Graduate
School of Engineering, Tohoku University. Since 2008, he has been also a
professor in the Department of Biomedical Engineering, Graduate Schoolof Biomedical Engineering, Tohoku University. His present interests arein transcutaneous measurement of heart wall vibrations and their spectrum analysis for diagnosis of myocardium and cross-sectional imaging ofelasticity around atherosclerotic plaque with transcutaneous ultrasoundfor tissue characterization of the arterial wall.
Dr. Kanai is a member of the IEEE, the Acoustical Society of Japan,
the Institute of Electronics Information and Communication Engineeringof Japan, the Japan Society of Ultrasonics in Medicine, Japan Society of
Medical Electronics and Biological Engineering, and the Japanese Circulation Society. Since 1998, he has been a member of Technical Program
Committee of the IEEE Ultrasonic Symposium.