1 November 2011 Toward soft-tissue elastography using digital holography to monitor surface acoustic waves
Author Affiliations +
J. of Biomedical Optics, 16(11), 116005 (2011). doi:10.1117/1.3646211
Measuring the elasticity distribution inside the human body is of great interest because elastic abnormalities can serve as indicators of several diseases. We present a method for mapping elasticity inside soft tissues by imaging surface acoustic waves (SAWs) with digital holographic interferometry. With this method, we show that SAWs are consistent with Rayleigh waves, with velocities proportional to the square root of the elastic modulus greater than 2-40 kPa in homogeneous tissue phantoms. In two-layer phantoms, the SAW velocity transitions approximately from that of the lower layer to that of the upper layer as frequency is increased in agreement with the theoretical relationship between SAW dispersion and the depth-dependent stiffness profile. We also observed deformation in the propagation direction of SAWs above a stiff inclusion placed 8 mm below the surface. These findings demonstrate the potential for quantitative digital holography-based elastography of soft tissues as a noninvasive method for disease detection.
Li, Mohan, Sanders, and Oldenburg: Toward soft-tissue elastography using digital holography to monitor surface acoustic waves



Because many illnesses are accompanied by a change in tissue elastic properties, medical practitioners often use palpation to “feel” the presence of disease. This underlies the motivation behind the development of biomedical imaging modalities that contrast soft-tissue elasticity,1 which has particular relevance in diseases of the liver,2 abdomen,3 blood vessels,4 and breast.5 Current methods for elastography include acoustic radiation force impulse imaging,6, 7 sonoelastography,8 and magnetic resonance (MR) elastography,9 the latter two of which work by monitoring low-frequency (typically, <1 kHz) elastic waves inside the tissue. In comparison to ultrasonic and magnetic resonance techniques, optical detection techniques have the potential to be lower cost and provide real-time imaging over a larger tissue surface area. In particular, coherent optical techniques, such as optical coherence elastography (OCE),10, 11, 12 detect nanometer tissue displacements and thus have the potential to provide improved sensitivity and resolution. Recent methods in OCE have employed sinusoidal actuation for higher specificity13 and for fitting the tissue response to elastic wave models14, 15 in analogy to the ultrasonic and MR methods. However, although OCE maps internal tissue deformations, its depth penetration is limited to only a few millimeters.16 In order to circumvent this limitation to direct optical measurement of internal tissue deformation, we propose to probe the depth-dependent tissue stiffness by measuring the elastic waves only on the tissue surface. This method is analogous to nondestructive testing techniques in materials analysis17, 18, 19, 20 and seismology,21 where internally mapped stiffness properties are inferred from surface acoustic waves (SAWs).

In fact, the use of SAWs for elastography of soft tissues is a growing area of interest,22, 23 but many of these experiments have been accomplished using a laser Doppler vibrometer, which cannot rapidly assess the entire tissue surface. Early studies using stroboscopic photography, while lacking nanoscale motion sensitivity clearly reveal the radial propagation of SAWs in human tissue, under point-source excitation.24 In comparison, holographic interferometry is capable of imaging nanoscale deformation over a large sample surface and has been widely used for vibrometry in materials analysis.25, 26 This technique has also been employed on biological tissues, where it was suggested that the method may be suitable for elasticity measurements.27 Related work using speckle pattern imaging demonstrates the detection of skin lesions based on their distinct mechanical properties.28 These techniques may also be suitable for breast cancer detection, because most malignant breast phenotypes are associated with an increased elastic modulus.29 The potential for holographic detection of breast abnormalities was demonstrated in human patients under quasi-static deformation.30 Under sinusoidal excitation, a more recent study of breast tissue phantoms revealed variations in the surface displacement as a function of the position of a tumorlike inhomogeneity.31 Although these results are promising, currently only empirical results have been demonstrated. Furthermore, because tumor location and breast size are highly variable, a versatile method with the ability to image at different depths is desirable.

In this paper, we describe a digital holographic interferometry method for elastography based on the quantitative description of SAWs as Rayleigh waves on an elastic half-space. Although this is only a first step toward 3-D reconstruction of elastograms, it provides needed intuition as to the nature of SAWs. In contrast to previous work,25, 26, 27, 30, 31 we employ a continuous, rather than a pulsed (i.e., stroboscopic), laser source, which is enabled by the increasing availability of high-speed digital cameras to capture fast tissue motion. We also developed a phase-shifting algorithm in the time domain for use with our on-axis holography setup, instead of the traditional Fourier transform or spatial phase-shifting techniques, which require careful control of the fringe spacing via off-axis holography.25, 32 A distinguishing feature of this algorithm is that it directly provides the phase of the mechanical wave (modulo π) using a priori knowledge of the SAW frequency, without the need to shift the reference beam phase. However, the key difference between this and previous work is that we measure the frequency-dependent velocity (i.e., the velocity dispersion of SAWs). This is important because the depth-dependent attenuation of the out-of-plane motion of SAWs is proportional to the SAW wavelength.33, 34 By sweeping the wavelength of SAWs, we therefore probe the tissue at varying depths, which enables reconstruction of the depth-dependent tissue elasticity. In this paper, we study the feasibility of this system for quantitative elastography using homogeneous and two-layer soft-tissue phantoms. We also demonstrate empirical results in phantoms with a local inhomogeneity to investigate the potential of this holographic interferometry technique for breast cancer detection.




Tissue-Phantom Preparation

Tissue phantoms with optical and mechanical properties comparable to that of soft tissues were used in these experiments, similar to those used previously.35 For the homogeneous phantoms, they were fabricated from a mixture of pure non–cross-linking polydimethylsiloxane (PDMS) fluid (50-cSt viscosity, ClearCo, Inc., Bensalem, Pennsylvania), and a cross-linking room-temperature PDMS base and its associated curing agent (Sylgard 184 A and B, respectively, Dow Corning, Inc. Midland, Wisconsin). The Young's moduli of the phantoms were adjusted over the range from 2 to 40 kPa by changing the ratio of the pure PDMS and the PDMS base over the range from 16:1 to 5:1, whereas the ratio between the PDMS base and its curing agent was held constant at 10:1. The relation between the mixture ratio and Young's modulus was estimated based on previous work.36 Titanium dioxide (TiO2) powder [Sigma-Aldrich Co., LLC (St. Louis, Missouri), No. 224227, mean size 1 μm, <5 μm] was embedded at a concentration of 4 mg/g in the tissue phantoms to function as tissuelike optical scatterers for holographic imaging. The phantom solutions were mixed thoroughly on an ultrasonicator for several tens of minutes at room temperature and then poured into a 126-mm-diam cylindrical container, with a thickness of ∼30 mm. The phantom dimensions were chosen to be large compared to the attenuation length of Rayleigh waves to eliminate significant boundary reflections. To further eliminate the boundary reflections, a few millimeters of the phantom adhering to the sidewall were removed. All phantoms were cured in an oven at 80°C for one day and subsequently at room temperature for another day. To prepare the two-layer phantoms, the same process was repeated for the top layer. The solution mixture was poured onto the cross-linked lower layer, and the whole phantom was cured as for the homogeneous phantoms. To construct a phantom with a local inhomogeneity, we prepared two homogeneous phantoms with the same mixing ratio and stacked them together, between which a thin metal ring with a 6.6 mm o.d., 3 mm i.d., and thickness of 0.5 mm was sandwiched.


Digital Holographic Interferometry

Figure 1 shows the in-line, image-plane holographic interferometry system for monitoring SAWs. A 633-nm helium-neon laser of 5 mW power (Thorlabs) was divided into sample and reference arms by a 92:8 beam splitter, where the majority of the optical power was directed to the sample arm. Two lenses in the sample arm acted as a 2.5× beam expander with the lens spacing adjusted so that the output of the last lens was slightly divergent, illuminating the phantom with a beam nominally 30 mm in diameter. A 2–in.-diam mirror (downward reflecting mirror in Fig. 1) directed the beam downward onto the horizontal surface of the phantom and subsequently directed the diffusely reflected light from the phantom through a 200-mm focal length, 2–in.-diam imaging lens, which imaged the light onto the sensor of a fast CMOS camera (Photron SA3). The phantom was mounted on a vertical stage with which the height was adjusted to focus the image. The camera frame rate was adjusted between 1 and 7.5 kHz to sample images at a rate ≥10× the frequency of the SAW. At these frame rates, the corresponding pixel array size varied between 1024 × 1024 and 256 × 256, respectively. The pixel size is 17 μm, and the imaging area for the full 1024 × 1024 pixel array was 16 × 16 mm2. Because there are highly spatially modulated speckles in the image, an adjustable aperture was placed close to the imaging lens. The aperture size was adjusted until the speckles could be distinguished by the camera sensor. The reference arm beam was passed through a spatial filter (i.e., the combination of a short focal length lens and a pinhole) and a 150-mm focal-length positive lens (f4 in Fig. 1), before being directed into the camera. The distance from the pinhole to the lens f4 was adjusted so that interference fringes between the reference and sample could be identified over the whole measurement area. Because the optical path of an edge ray from the sample to the camera is greater than that of a central ray, the output light from the lens f4 was chosen to be divergent. Ideally, the optical path difference (OPD) between edge rays and central rays should be balanced for the two arms. A neutral density filter was inserted in the reference arm to optimize the interference signal within the dynamic range of the camera. A retroreflector was used in the reference arm to adjust the OPD between the two arms so that it was within the coherence length of the laser, ∼200 mm.

Fig. 1

Schematic diagram of holographic elastography experimental setup. NDF: neutral density filter, f1 = 40 mm, f2 = 100 mm, f3 = 5 mm; f4 = 150 mm, and f5 = 200 mm.


The image of the phantom interacted with the reference beam, forming a local holographic image. As SAWs propagated across the phantom surface, surface displacement along the axis that bisects the propagation directions of the incident and detected scattered fields (approximately the axis normal to the phantom surface) induced a phase shift in the interferogram due to the changing OPD. As such, this system was predominantly sensitive to out-of-plane sample motions. For generation of SAWs, a piezoelectric transducer (Kinetic Ceramics, Inc., Hayward, California) with a travel range of 20 μm was placed into contact with the phantom immediately to one side of the illumination area. The transducer contact area could be regarded as a point because its diameter is much less than the SAW wavelengths in our experiments. The transducer was driven by a high-voltage power supply controlled by a function generator. Typically, the transducer was initiated at a fixed frequency and the system was allowed to reach a steady state (∼1−2 s), after which the camera recorded the holographic images in real time. On-board image data were transferred to the computer for further processing, as described later.


SAW Phase Extraction

The SAWs propagate outward from the point of actuation, which cause temporally modulated phase shifts in the holographic image. When the phantom's thickness is much greater than the SAW wavelength, we expect the SAW to be characterized as a Rayleigh wave.33 In that case, the displacement in the z direction (out of plane) can be expressed as [TeX:] $\Delta z\left({\vec r,t} \right) = A\cos ({\omega t + \vec k \cdot \vec r})$ Δzr,t=Acos(ωt+k·r) , where A is the SAW amplitude (assumed to be slowly varying in [TeX:] $\vec r$ r , roughly as the inverse square root of [TeX:] $\left| {\vec r} \right|$ r ), [TeX:] $\vec r$ r is the position in the (x,y) surface plane, ω is the angular frequency of the SAW, and [TeX:] $\vec k$ k is the spatial frequency of the SAW. The intensity recorded on the camera can then be written as follows:


[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} I\left({\vec r,t} \right) = I_{{\rm DC}} + I_{{\rm AC}} \cos [ {\varphi _{{\rm SAW}} \cos ({\omega t + \vec k \cdot \vec r}) + \phi ({\vec r})}], \end{equation}\end{document} Ir,t=IDC+IACcos[ϕSAWcos(ωt+k·r)+φ(r)],
where φSAW = 4πA/λ is the amplitude of the optical phase shift induced by the SAW, [TeX:] $\phi \left({\vec r} \right)$ φr is the optical phase difference between the reference and sample fields in the absence of the SAW, I DC is the sum of the reference and sample beam intensities, I AC is the cross-spectral density of the reference and sample beams, and λ is the laser wavelength.

To extract the phase of the SAW, that is, [TeX:] $\beta _{{\rm SAW}} \left({\vec r,t} \right)\break = \omega t + \vec k \cdot \vec r$ βSAWr,t=ωt+k·r , in the limit of small SAW amplitudes (Δz ≪ λ) and sufficiently fast temporal sampling (ωΔt ≪ 2π, where Δt is the time between successive frames), it is possible to use a temporal phase-shifting algorithm similar to Refs. 37 and 38, such that


[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{eqnarray} \beta _{{\rm SAW}} ({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t}) &=& \tan ^{ - 1} \left[\! {\frac{{I({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t {-} \Delta t}) {-} I({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t {+} \Delta t})}}{{I({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t {-} \Delta t}) {+} I({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t {+} \Delta t}) {-} 2I({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t})}}}\right.\hspace*{-5pt}\nonumber\\ && \times\left. {\tan \left({\frac{{\omega \Delta t}}{2}} \right)} \vphantom{\left[\! {\frac{{I\left({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t {-} \Delta t} \right) {-} I\left({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t {+} \Delta t} \right)}}{{I\left({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t {-} \Delta t} \right) {+} I\left({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t {+} \Delta t} \right) {-} 2I\left({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t} \right)}}}\right.}\right], \end{eqnarray}\end{document} βSAW(r,t)=tan1I(r,tΔt)I(r,t+Δt)I(r,tΔt)+I(r,t+Δt)2I(r,t)×tanωΔt2,

where a group of three successive frames at times t – Δt, t, and t + Δt are processed to produce the SAW phase image at time t, and [TeX:] $\mathord{\buildrel{\lower3pt \hbox{\scriptscriptstyle\rightharpoonup}}\over r}$ r is the (x,y) pixel position in the frame. Briefly, we can understand why this works by rewriting the phase-shifting algorithm in terms of time derivatives of I and using the small angle approximation,


[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{eqnarray} \beta _{{\rm SAW}} ({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t}) &=& \tan ^{ - 1} \left[ {\frac{{I'({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t})}}{{I''({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t})}} \cdot \frac{2}{{\Delta t}}\tan \left({\frac{{\omega \Delta t}}{2}} \right)} \right]\nonumber\\ & \approx& \tan ^{ - 1} \left({\frac{{I'({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t})}}{{I''({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t})}}\omega } \right). \end{eqnarray}\end{document} βSAW(r,t)=tan1I(r,t)I(r,t)·2ΔttanωΔt2tan1I(r,t)I(r,t)ω.
If we then perform a Taylor expansion of I in Eq. 1 about φSAW (in the limit of small SAW amplitudes), then we obtain


[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{eqnarray} I\left({\vec r,t} \right) &\approx& I_{{\rm DC}} + I_{{\rm AC}} \cos [ {\phi \left({\vec r} \right)} ]\nonumber\\ && -\, I_{{\rm AC}} \varphi _{{\rm SAW}} \sin [ {\phi \left({\vec r} \right)}]\cos ({\omega t + \vec k \cdot \vec r}), \end{eqnarray}\end{document} Ir,tIDC+IACcos[φr]IACϕSAWsin[φr]cos(ωt+k·r),
the time derivatives of which become


[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} \begin{array}{l} I'({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t}) \approx \omega I_{{\rm AC}} \varphi _{{\rm SAW}} \sin [ {\phi ({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r} })}]\sin ({\omega t + \mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over k} \cdot \mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r} }), \\[8pt] I''({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r},t}) \approx \omega ^2 I_{{\rm AC}} \varphi _{{\rm SAW}} \sin [ {\phi ({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r} })} ]\cos ({\omega t + \mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over k} \cdot \mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}} \over r} }). \\ \end{array} \end{equation}\end{document} I(r,t)ωIACϕSAWsin[φ(r)]sin(ωt+k·r),I(r,t)ω2IACϕSAWsin[φ(r)]cos(ωt+k·r).

Plugging these back into Eq. 3 reveals how the extracted phase βSAW using the phase-shifting algorithm in Eq. 2 is equivalent to the phase of the mechanical wave, [TeX:] $\omega t + \vec k \cdot \vec r$ ωt+k·r . It should be noted that our method is nontraditional in that the reference beam is not phase shifted and the optical phase is not obtained directly; instead, we employ a priori knowledge of the mechanical wave frequency to obtain the mechanical phase. This has the advantage of skipping the unnecessary evaluation of the optical phase, in favor of directly assessing the mechanical phase. One consequence is that we do not know at any instant in time whether the SAW is in its upswing or downswing and, therefore, we can only measure βSAW modulo π. In our phase images, therefore, we represent βSAW = –π/2 to π/2 as a gray-scale value from black to white, respectively.

Figure 2a (Video 1) shows an example interference pattern measured by the camera, [TeX:] $I({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}}\over r},t})$ I(r,t) . In the video, the evenly distributed fringes aligned from lower right to upper left are background fringes caused by the protection glass in front of the camera sensor due to multiple reflections (Fabry–Perot). These are stationary and do not affect the measurement result using the temporal phase-shifting algorithm. The SAW phase image, [TeX:] $\beta _{{\rm SAW}} ({\mathord{\buildrel{\lower3pt\hbox{\scriptscriptstyle\rightharpoonup}}\over r},t})$ βSAW(r,t) , which was extracted from Fig. 2a according to Eq. 3, is shown in Fig. 2b (Video 2). It reveals a SAW propagating away from the actuator position located outside of the image to the right. Note that the phase-wrap position from βSAW = π/2 to –π/2 provides a highly contrasted way to visualize the SAW wavefront, which we also verified with simulations. We then obtain the SAW phase velocity by tracking a specific wavefront in the image series. Although we found human inspection to be sufficient for SAW velocity measurement within the scope of this work, in future work we plan to develop an automated method for SAW velocity extraction.

Fig. 2

Images of a phantom during mechanical excitation by a piezoelectric transducer with a driving frequency of 180 Hz. The camera frame rate is 2 kHz, and the measurement area is 10 × 10 mm2. (a) Holographic images recorded by the camera (Video 1) and (b) the associated extracted SAW phase according to the methods in Sec. 2.3 (Video 2). (Video 1, QuickTime, 356 KB); (Video 2, QuickTime, 472 KB)


It is important to note that modulation of the speckle is the dominant effect giving rise to the SAW phase signal. This is because the phantoms were designed so that the diffuse scattering from the TiO2 micropowder dominates over the specular reflection from the surface in order to more accurately mimic the optical signal from real tissues. This experimental configuration can actually be thought of as a specific case of digital holography known as electronic speckle pattern interferometry.37, 39, 40


Results and Discussion


Surface Acoustic Waves in Homogeneous Phantoms with Varying Young's Moduli

Our hypothesis is that SAWs on a soft-tissue–like, homogeneous, elastic half-space will behave as Rayleigh waves.33 consistent with other observations.23 Here, we study this phenomenon using digital holographic interferometry. We generated tissue phantoms with Young's moduli ranging from 2 to 40 kPa to match the range of elastic moduli encountered in soft tissue such as the human breast.29 According to elastic wave theory, the phase velocity of Rayleigh waves in a homogeneous half-space depends on the Young's modulus (E) and density of the material (ρ) according to23


[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} c_{\rm R} \approx \frac{1}{{1.05}}\sqrt {\frac{E}{{3\rho }}}, \end{equation}\end{document} cR11.05E3ρ,
assuming the Poisson's ratio of the material is 0.5 [a reasonable assumption for PDMS (Ref. 41)].

To produce SAWs, we excited the phantoms with a 150-Hz sinusoidal signal on the piezoelectric transducer. This choice of frequency was such that the corresponding SAW wavelength (and, thus, penetration depth) was smaller than the phantom depth (30 mm). This produces a surface wave with constant phase velocity due to the elastic homogeneity of the material. Using the SAW phase-extraction method previously described, we quantified the wave speed as a function of Young's modulus and plotted the results in Fig. 3. Each data point was obtained by averaging eight measurements for a single sample. The solid line is the theoretical curve calculated with Eq. 6 using ρ = 960 kg/m3. The dashed line is a least-squares best-fit curve with one fitting parameter, k 0, such that c R = k 0√E. The fitting value of k 0 = 0.525 is within 6% of the theoretical value of 0.560. The R 2 value of the fit was 0.99. For additional validation, we measured the SAW speeds for one homogeneous phantom while sinusoidally driving the transducer at different frequencies (75−300 Hz) and found that the measured velocity remained constant with a standard deviation of ∼7.5% of the mean. This demonstrates that the induced out-of-plane SAW motions measured by holographic interferometry are well described as Rayleigh waves.

Fig. 3

Variation in SAW speed as a function of the Young's modulus in homogeneous tissue phantoms follows the expected relationship of Eq. 6.



Surface Acoustic Wave Velocity Dispersion in Two-Layer Phantoms

Because the penetration depth of Rayleigh waves is proportional to the wavelength, it is known that a variation in the material's elasticity versus depth gives rise to dispersion in the SAW velocity.17, 18, 19, 20, 21 We hypothesize that SAWs in soft tissues will exhibit a velocity according to Eq. 6 given by the average material elasticity E within the surface depth probed by the wave. In other words, we expect the long-wavelength (low-frequency) SAWs to probe the tissue response over a large depth range, while the short-wavelength (high-frequency) SAWs should probe only the upper tissue layers.

In order to test this hypothesis, we generated two two-layer phantoms, one with a soft upper layer and stiff lower layer, and a second with a stiff upper layer and soft lower layer. For the first phantom, the Young's modulus E of the upper layer was 7 kPa [corresponding to a theoretical velocity of 1.4 m/s according to Eq. 6], and the lower layer was 22 kPa (corresponding to a velocity of 2.6 m/s). The top layer was 4 mm thick, and the bottom layer was 10 mm thick. The measured SAW speed at different driving frequencies is shown in Fig. 4a, where the dashed lines indicate the theoretical velocities in a homogenous medium composed of the upper or lower layer. For the second phantom, the Young's modulus E of the upper layer was 38 kPa (corresponding to the velocity of 3.4 m/s) and the lower layer was 11 kPa (corresponding to the velocity of 1.8 m/s). The top layer was 6 mm thick, and the bottom layer was 30 mm thick. The measured SAW speed at different driving frequencies is shown in Fig. 4b.

Fig. 4

SAW speed variation with driving frequency for (a) a phantom with upper softer layer and (b) a phantom with upper stiffer layer.


From the two plots, it can be seen that the SAW velocity gradually transitions from close to the velocity of the lower layer to that of the upper layer with the increase of the driving frequency. This expected behavior is observed in both cases, with the stiffer layer either being on top or underneath, such that a positive dispersion is indicative of decreasing stiffness versus depth, and a negative dispersion that of increasing stiffness versus depth. In future work, this measured dispersion profile can be used to quantitatively predict the depth-dependent elasticity either in the perturbative regime18, 19 or by treating the exact solution as an inverse problem.17, 20

Qualitatively, we can compare the transition wavelength between lower and upper wave velocities to the layer thickness. This is obtained at the point at which the slope dv/df is the maximum and is observed to occur at approximately 100 and 450 Hz in each of the phantoms, respectively. The associated SAW wavelengths, using the velocities measured at each of these frequencies, are 18 and 7 mm, respectively. In the first phantom with a 4-mm-thick soft upper layer, the transition wavelength of 18 mm is much larger than the layer thickness, suggesting significant decay of the wave has reduced its effective probing depth. In the second phantom with a 6-mm-thick stiff upper layer, the transition wavelength was ∼7 mm. Depending on Poisson's ratio, out-of-plane Rayleigh wave motion is reduced to ∼20–33% of its surface value at one wavelength in depth.34 We might therefore expect the transition to occur when the wavelength is on the order of the upper layer thickness. However, comparing the transition wavelength with upper layer thickness in the two phantoms suggests that the softer phantom material exhibits greater decay than the stiffer material. This suggestion is also supported by the observation that the SAW velocity in the phantom in Fig. 4b does not quite reach that of the soft lower layer at low frequency. If the attenuation in the soft layer is too rapid, then the wave will tend to bias toward a velocity given by the stiffer upper layer. Further investigation is needed to determine whether the softer phantom material exhibits more dissipative loss, which might be expected because the softer phantom contains a lower concentration of cross-linking PDMS (which is dominantly elastic) relative to the non–cross-linking PDMS fluid (which is dominantly viscous).

In experiments, we also note that the observable range of SAW frequencies (frequencies where SAWs had sufficient amplitude to provide a detectable phase variation) were limited in practice. No SAWs were observed at <20 Hz, while at >80 Hz, SAWs were observed in every experiment. We found the high-frequency limit to lie between 300 Hz for phantoms with a soft upper layer and 800 Hz for phantoms with a stiff upper layer. Improvements in phase sensitivity or actuator sweep may increase the detectable frequency range, in order to increase the depth penetration. Thus, we expect the spatial resolution for elastographic reconstruction to be fundamentally limited by the SAW bandwidth and signal-to-noise ratio. To speed the data acquisition, a chirped (frequency-swept) SAW waveform can be applied to rapidly measure the dispersion curve; the rate of chirp can be slowed to provide the desired frequency resolution, depending on the needed spatial resolution.


Surface Acoustic Wave Propagation in Phantoms with a Local Inhomogeneity

To investigate the potential of this technique for detecting local inhomogeneities, we prepared an otherwise homogeneous soft phantom (7 kPa) with an embedded thin metal ring as described in Sec. 2.1. The embedding depth was 8 mm below the surface. Figure 5 shows the phase reconstruction of the SAW propagation at frequencies of 65 Hz (Video 3) and 89 Hz (Video 4) using a frame rate of 1 kHz, and 200 Hz (Video 5) using a frame rate of 2 kHz. These frequencies nominally correspond to wavelengths of 22, 16, and 7 mm, respectively. The measurement area was 16 × 16 mm2 at 1 kHz and 10 × 10 mm2 at 2 kHz. The effect of the ring is clearly observed as a deformation of the SAW propagation, and the shape of this deformation appears to change at different frequencies due to the varying penetration depth. In particular, at the frequency of 200 Hz, the impact of the ring is minimized because its wavelength of 7 mm is shorter than the ring depth of 8 mm, while the impact is larger at the shorter wavelengths. The impact appears to be largest at 89 Hz, although significant wavefront deformation is observed in the upper left corner of the images at 65 Hz as well. This demonstrates that this technique is sensitive to local inhomogeneities at least as deep as 8 mm. The capability to identify local inhomogeneities is useful because tumors and other disease states are often characterized by a palpable, local elastic abnormality.

Fig. 5

SAW propagation in a phantom with a local inhomogeneity at frequencies of 65 Hz (a) (Video 3), 89 Hz (b) (Video 4), and 200 Hz (c) (Video 5). The image at 200 Hz has a smaller field of view due to the increased frame rate. The dashed circle indicates the lateral position of the metal ring. (Video 3, QuickTime, 728 KB); (Video 4, QuickTime, 878 KB); (Video 5, QuickTime, 986 KB)




This paper demonstrates the potential for a new elastography technique for soft tissue by monitoring the out-of-plane motion of SAWs with digital holography. In homogeneous and two-layer soft-tissue phantoms, SAWs display the characteristics of Rayleigh waves. SAW propagation is also clearly affected in response to a stiff inclusion placed 8 mm beneath the tissue phantom surface. As such, the frequency-dependent velocity of SAWs contains information about the depth-dependent tissue elasticity (Fig. 4), whereas the lateral propagation direction of SAWs contains information about the lateral profile of the tissue elasticity (Fig. 5). By deconstructing this system into simplified phantoms, we have provided needed intuition as to how SAWs behave in response to elastic inhomogeneities.

Before this system can be used in living subjects, motion artifacts (e.g., from cardiopulmonary function) will need to be addressed. Previous experiments using modulated contrast methods in vivo with phase-resolved optical coherence tomography42 suggest that bandpass filtering around the excitation frequency may be sufficient to reject these types of motion artifacts.

We note that current depth-profiling methods using SAWs17, 18, 19, 20 have been performed in elastic media, and future work is needed to extend these for use in viscoelastic media (where the mechanical properties of tissue may be a function of frequency).43 It may be possible to reduce the viscoelastic parameter space by characterizing the scaling of the tissue frequency dependence with appropriate models such as that in Ref. 44 to enable depth profiling in viscoelastic tissues. In future work, we hope that this will lead to reconstruction methods based on elastic scattering theory to convert the observed SAW behavior into 3-D tissue elastograms.45 This purely optical technique is noninvasive, takes advantage of the high displacement sensitivity afforded by interferometry, and provides depth-dependent elasticity information potentially several centimeters into soft tissue. In particular, holographic interferometry allows us to assess a large surface area simultaneously. Our hope is that it may eventually lead to a new clinical tool for quantitative palpation.


We thank Vasilica Crecea, at the University of Illinois at Urbana-Champaign, for providing calibration data for the Young's modulus of the silicone tissue phantoms. This work was supported by startup funds from the University of North Carolina at Chapel Hill.


1.  J. Ophir, I. Céspedes, H. Ponnekanti, Y. Yazdj, and X. Li, “Elastography: a quantitative method for imaging the elasticity of biological tissues,” Ultrason. Imaging 13, 111–134 (1991). 10.1016/0161-7346(91)90079-W Google Scholar

2.  L. Sandrin, B. Fourquet, J. Hasquenoph, S. Yon, C. Fournier, F. Mal, C. Christidis, M. Ziol, B. Poulet, F. Kazemi, M. Beaugrand, and R. Palau, “Transient elastography: a new noninvasive method for assessment of hepatic fibrosis,” Ultrasound Med. Biol. 29(12), 1705–1713 (2003). 10.1016/j.ultrasmedbio.2003.07.001 Google Scholar

3.  B. J. Fahey, R. C. Nelson, D. P. Bradway, S. J. Hsu, D. M. Dumont, and G. E. Trahey, “In vivo visualization of abdominal malignancies with acoustic radiation force elastography,” Phys. Med. Biol. 53, 279–293 (2008). 10.1088/0031-9155/53/1/020 Google Scholar

4.  C. Schmitt, G. Soulez, R. L. Maurice, M. Giroux, and G. Cloutier, “Noninvasive vascular elastography: toward a complementary characterization tool of atherosclerosis in carotid arteries,” Ultrasound Med. Biol. 33(12), 1841–1858 (2007). 10.1016/j.ultrasmedbio.2007.05.020 Google Scholar

5.  D. B. Plewes, J. Bishop, A. Samani, and J. Sciarretta, “Visualization and quantification of breast cancer biomechanical properties with magnetic resonance elastography,” Phys. Med. Biol. 45, 1591–1610 (2000). 10.1088/0031-9155/45/6/314 Google Scholar

6.  A. Y. Iyo, “Acoustic radiation force impulse imaging: a literature review,” J. Diagn. Med. Sonog. 25, 204–211 (2009). 10.1177/8756479309335741 Google Scholar

7.  K. Nightingale, M. S. Soo, R. Nightingale, and G. Trahey, “Acoustic radiation force impulse imaging: in vivo demonstration of clinical feasibility,” Ultrasound Med. Biol. 28(2), 227–235 (2002). 10.1016/S0301-5629(01)00499-9 Google Scholar

8.  L. S. Taylor, B. C. Porter, D. J. Rubens, and K. J. Parker, “Three-dimensional sonoelastography: principles and practices,” Phys. Med. Biol. 45, 1477–1494 (2000). 10.1088/0031-9155/45/6/306 Google Scholar

9.  R. Muthupillai, D. J. Lomas, P. J. Rossman, J. F. Greenleaf, A. Manduca, and R. L. Ehman, “Magnetic resonance elastography by direct visualization of propagating acoustic strain waves,” Science 269, 1854–1857 (1995). 10.1126/science.7569924 Google Scholar

10.  J. M. Schmitt, “OCT elastography: imaging microscopic deformation and strain of tissue,” Opt. Express 3(6), 199–211 (1998). 10.1364/OE.3.000199 Google Scholar

11.  R. C. Chan, A. H. Chau, W. C. Karl, S. Nadkarni, A. S. Khalil, N. Iftimia, M. Shishkov, G. J. Tearney, M. R. Kaazempur-Mofrad, and B. E. Bouma, “OCT-based arterial elastography: robust estimation exploiting tissue biomechanics,” Opt. Express 12(19), 4558–4572 (2004). 10.1364/OPEX.12.004558 Google Scholar

12.  S. J. Kirkpatrick, R. K. Wang, and D. D. Duncan, “OCT-based elastography for large and small deformations,” Opt. Express 14(24), 11585–11597 (2006). 10.1364/OE.14.011585 Google Scholar

13.  B. F. Kennedy, X. Liang, S. G. Adie, D. K. Gerstmann, B. C. Quirk, S. A. Boppart, and D. D. Sampson, “In vivo three-dimensional optical coherence elastography,” Opt. Express 19(7), 6623–6634 (2011). 10.1364/OE.19.006623 Google Scholar

14.  A. L. Oldenburg and S. A. Boppart, “Resonant acoustic spectroscopy of soft tissues using embedded magnetomotive nanotransducers and optical coherence tomography,” Phys. Med. Biol. 55, 1189–1201 (2010). 10.1088/0031-9155/55/4/019 Google Scholar

15.  X. Liang, A. L. Oldenburg, V. Crecea, E. J. Chaney, and S. A. Boppart, “Optical micro-scale mapping of dynamic biomechanical tissue properties,” Opt. Express 16(15), 11052–11065 (2008). 10.1364/OE.16.011052 Google Scholar

16.  J. M. Schmitt, “Optical coherence tomography (OCT): a review,” IEEE J. Sel. Top. Quant. 5(4), 1205–1215 (1999). 10.1109/2944.796348 Google Scholar

17.  C. Glorieux, W. Gao, S. E. Kruger, K. V. de Rostyne, W. Lauriks, and J. Thoen, “Surface acoustic wave depth profiling of elastically inhomogeneous materials,” J. Appl. Phys. 88(7), 4394–4440 (2000). 10.1063/1.1290457 Google Scholar

18.  T. L. Szabo, “Obtaining subsurface profiles from surface-acoustic-wave velocity dispersion,” J. Appl. Phys. 46(4), 1448–1454 (1975). 10.1063/1.321793 Google Scholar

19.  B. R. Tittman, L. A. Ahlberg, J. M. Richardson, and R. B. Thompson, “Determination of physical property gradients from measured surface wave dispersion,” IEEE Trans. Ultrason. Ferr. UFFC-34(5), 500–507 (1987). 10.1109/T-UFFC.1987.26975 Google Scholar

20.  G. A. Gordon and B. R. Tittmann, “Surface acoustic wave determination of subsurface structure,” Nondestr. Test. Eval. 11, 21–32 (1994). 10.1080/10589759408952814 Google Scholar

21.  D. G. Harkrider, “Surface waves in multilayered elastic media I. Rayleigh and Love waves from buried sources in a multilayered elastic half-space,” Bull. Seismol. Soc. Am. 54(2), 627–679 (1964). Google Scholar

22.  T. J. Royston, H. A. Mansy, and R. H. Sandler, “Excitation and propagation of surface waves on a viscoelastic half-space with application to medical diagnosis,” J. Acoust. Soc. Am. 106(6), 3678–3686 (1999). 10.1121/1.428219 Google Scholar

23.  X. Zhang and J. F. Greenleaf, “Estimation of tissue's elasticity with surface wave speed (L),” J. Acoust. Soc. Am. 122(5), 2522–2525 (2007). 10.1121/1.2785045 Google Scholar

24.  H. E. von Gierke, H. L. Oestreicher, E. K. Franke, H. O. Parrack, and W. W. von Wittern, “Physics of vibrations in living tissues,” J. Appl. Physiol. 4(12), 886–900 (1952). Google Scholar

25.  G. Pedrini, S. Schedin, and H. J. Tiziani, “Pulsed digital holography combined with laser vibrometry for 3D measurements of vibrating objects,” Opt. Laser Eng. 38, 117–129 (2002). 10.1016/S0143-8166(02)00005-2 Google Scholar

26.  C. Trillo, D. Cernadas, Á. F. Doval, C López, B. V. Dorrío, and J. L. Fernández, “Detection of transient surface acoustic waves of nanometric amplitude with double-pulsed TV holography,” Appl. Opt. 42(7), 1228–1235 (2003). 10.1364/AO.42.001228 Google Scholar

27.  S. Schedin, G. Pedrini, and H. J. Tiziani, “Pulsed digital holography for deformation measurements on biological tissues,” Appl. Opt. 39(16), 2853–2857 (2000). 10.1364/AO.39.002853 Google Scholar

28.  S. J. Kirkpatrick, R. K. Wang, D. D. Duncan, M. Kulesz-Martin, and K. Lee, “Imaging the mechanical stiffness of skin lesions by in vivo acousto-optical elastography,” Opt. Express 14(21), 9770–9779 (2006). 10.1364/OE.14.009770 Google Scholar

29.  A. Samani, J. Zubovits, and D. Plewes, “Elastic moduli of normal and pathlogical human breast tissues: an inversion-technique-based investigation of 169 samples,” Phys. Med. Biol. 52, 1565–1576 (2007). 10.1088/0031-9155/52/6/002 Google Scholar

30.  J. Woisetschlager, D. B. Sheffer, C. W. Loughry, K. Somasundaram, S. K. Chawla, and P. J. Wesolowski, “Phase-shifting holographic interferometry for breast cancer detection,” Appl. Opt. 33(22), 5011–5015 (1994). 10.1364/AO.33.005011 Google Scholar

31.  M. del S. Hernández-Montes, C. Pérez-López, and F. M. Santoyo, “Finding the position of tumor inhomogeneities in a gel-like model of a human breast using 3-D pulsed digital holography,” J. Biomed. Opt. 12(2), 024027 (2007). 10.1117/1.2717515 Google Scholar

32.  L. Mertz, “Real-time fringe-pattern analysis,” Appl. Opt. 22(10), 1535–1539 (1983). 10.1364/AO.22.001535 Google Scholar

33.  L. Rayleigh, “On waves propagating along the plane surface of an elastic solid,” Proc. London Math. Soc. S1-17(1), 4–11 (1885). 10.1112/plms/s1-17.1.4 Google Scholar

34.  K. F. Graff, Wave Motion in Elastic Solids, Dover, New York (1991). Google Scholar

35.  A. L. Oldenburg, F. Jean-Jacques Toublan, K. S. Suslick, A. Wei, and S. A. Boppart, “Magnetomotive contrast for in vivo optical coherence tomography,” Opt. Express 13(17), 6597–6614 (2005). 10.1364/OPEX.13.006597 Google Scholar

36.  V. Crecea, A. L. Oldenburg, X. Liang, T. S. Ralston, and S. A. Boppart, “Magnetomotive nanoparticle transducers for optical rheology of viscoelastic materials,” Opt. Express 17(25), 23114–23122 (2009). 10.1364/OE.17.023114 Google Scholar

37.  G. Pedrini, H. J. Tiziani, and Y. Zou, “Digital double pulse-TV-holography,” Opt. Laws Eng. 26, 199–219 (1997). 10.1016/0143-8166(95)00113-1 Google Scholar

38.  M. Hipp, J. Woisetschläger, P. Reiterer, and T. Neger, “Digital evaluation of interferograms,” Measurement 36, 53–66 (2004). 10.1016/j.measurement.2004.04.003 Google Scholar

39.  J. M. Huntley, G. H. Kaufmann, and D. Kerr, “Phase-shifted dynamic speckle pattern interferometry at 1 kHz,” Appl. Opt. 38(31), 6556–6563 (1999). 10.1364/AO.38.006556 Google Scholar

40.  O. Løkberg, “Electronic speckle pattern interferometry,” Phys. Technol. 11, 16–22 (1980). 10.1088/0305-4624/11/1/303 Google Scholar

41.  D. Armani, C. Liu, and N. Alum, “Re-configurable fluid circuits by PDMS elastomer micromachining,” in Proc. of IEEE Conf. on Micro Electro Mechanical Systems (MEMS), pp. 222–227 (1999). Google Scholar

42.  R. John, R. Rezaeipoor, S. G. Adie, E. J. Chaney, A. L. Oldenburg, M. Marjanovic, J. P. Haldar, B. Sutton, and S. A. Boppart, “In vivo magnetomotive optical molecular imaging using targeted magnetic nanoprobes,” Proc. Natl. Acad. Sci. 107(18), 8085–8090 (2010). 10.1073/pnas.0913679107 Google Scholar

43.  R. Sinkus, K. Siegmann, T. Xydeas, M. Tanter, C. Claussen, and M. Fink, “MR elastography of breast lesions: understanding the solid/liquid duality can improve the specificity of contrast-enhanced MR mammography,” Magn. Reson. Med. 58, 1135–1144 (2007). 10.1002/mrm.21404 Google Scholar

44.  F. Can Meral, T. J. Royston, and R. L. Magin, “Surface response of a fractional order viscoelastic halfspace to surface and subsurface sources,” J. Acoust. Soc. Am. 126(6), 3278–3285 (2009). 10.1121/1.3242351 Google Scholar

45.  D. Colton and R. Kress, Inverse Acoustic and Electromagnetic Scattering Theory, Springer, New York (1993). Google Scholar

© 2011 Society of Photo-Optical Instrumentation Engineers (SPIE)
Shiguang Li, Karan D. Mohan, William W. Sanders, Amy L. Oldenburg, "Toward soft-tissue elastography using digital holography to monitor surface acoustic waves," Journal of Biomedical Optics 16(11), 116005 (1 November 2011). https://doi.org/10.1117/1.3646211

Back to Top