Estimated Viscoelasticity of Voigt tissue having
Studies have found that shear moduli, having the dynamic range of several orders of magnitude for various biological tissues , are highly correlated with the pathological statues of human tissue such as livers [2, 3]. The shear moduli can be investigated by measuring the attenuation and velocity of the shear wave propagation in a tissue region. Many efforts have been made to measure shear wave propagations induced by different types of force, which include the motion force of human organs, external applied force , and ultrasound radiation force .
In past 15 years, ultrasound radiation force has been successfully used to induce tissue motion for imaging tissue elasticity. Vibroacoustography (VA) uses bifocal beams to remotely induce vibration in a tissue region and detect the vibration using a hydrophone . The vibration center is sequentially moved in the tissue region to form a two-dimensional image. Acoustic Radiation Force Imaging (ARFI) uses focused ultrasound to apply localized radiation force to small volumes of tissue for short durations and the resulting tissue displacements are mapped using ultrasonic correlation based methods . Supersonic shear image remotely vibrates tissue and sequentially moves vibration center along the beam axis to create intense shear plan wave that is imaged at a high frame rate (5000 frames per second) . These image methods provide measurements of tissue elasticity, but not the viscosity.
Because of the dispersive property of biological tissue, the induced tissue displacement and the shear wave propagation are frequency dependent. Tissue shear property can be modeled by several models including Kelvin-Voigt (Voigt) model, Maxwell model, and Zener model . Voigt model effectively describes the creep behavior of tissue, Maxwell model effectively describes the relaxation process, and the Zener model effectively describes both creep and relaxation but it requires one extra parameter. Voigt model is often used by many researchers because of its simplicity and the effectiveness of modeling soft tissue. Voigt model consists of a purely viscous damper and a purely elastic spring connected in parallel. For Voigt tissue, the tissue motion at a very low frequency largely depends on the elasticity, while the motion at a very high frequency largely depends on the viscosity . In general, the tissue motion depends on both elasticity and viscosity, and estimates of elasticity by ignoring viscosity are biased or erroneous.
Back to the year of 1951, Dr. Oestreicher published his work to solve the wave equation for the Voigt soft tissue with harmonic motions . With assumptions of isotropic tissue and plane wave, he derived equations that relate the shear wave attenuation and speed to the elasticity and viscosity of soft tissue. However, Oestreicher’s method was not realized for applications until the half century later.
In the past ten years, Oestreicher’s method was utilized to quantitatively measure both tissue elasticity and viscosity. Ultrasound vibrometry has been developed to noninvasively and quantitatively measure tissue shear moduli [10-16]. It induces shear waves using ultrasound radiation force [5, 6] and estimates the shear moduli using shear wave phase velocities at several frequencies by measuring the phase shifts of the propagating shear wave over a short distance using pulse echo ultrasound [10-16]. Applications of the ultrasound vibrometry were conducted for viscoelasticities of liver , bovine and porcine striated muscles [17, 18], blood vessels [12, 19-21], and hearts . A recent
One of potential applications of the ultrasound vibrometry is to characterize shear moduli of livers. The shear moduli of liver are highly correlated with liver pathology status [24, 25]. Recently, the shear viscoelasticity of liver tissue has been investigated by several research groups [23, 26-28]. The most of these studies applied ultrasound radiation force in liver tissue regions, measured the phase velocities of shear wave in a limited frequency range, and inversely solved the Voigt model with an assumption that liver local tissue is isotropic without considering boundary conditions. Because of the boundary conditions, shear wave propagations are impacted by the limited physical dimensions of tissue. Studies shows that considerations of boundary conditions should be taken for characterizing tissue that have limited physical dimensions such as heart , blood vessels [19-21], and liver , when ultrasound vibrometry is used.
2. Shear wave propagation in soft tissue and shear viscoelasticity
The shear wave propagation in soft tissue is a complicated process. When the tissue is isotropic and modeled by the Voigt model, the phase velocity and attenuation of the shear wave propagation in the tissue are associated with tissue viscoelasticity. Oesteicher documented the detailed derivations of the solution of the sound wave equation for Voigt tissue . We extended the solution to other models  for the applications of ultrasound vibrometry . In this section, we provide the simplified descriptions of the shear wave propagation in tissue modeled by Voigt model, Maxwell model, and Zener model.
Assuming that a harmonic motion produces the shear wave that propagates in a tissue region, the phase velocity
The phase velocity is associated with the tissue property, which can be found by solving the wave equation with a tissue viscoelasticity model. For a small local region, the wave is approximated as a uniform plane wave, which has a simple form in isotropic medium:
where S is the phasor notation of the displacement of the time-harmonic field of the shear wave, z is the wave propagation distance which is perpendicular to the direction of the displacement of the shear wave, and the complex wave number is
The solution of (2) is a standard solution of a homogeneous wave equation:
Although attenuation coefficient α = –ki carries information of the complex modulus of tissue, the phase measurement is often more reliable because it is relatively independent to transducers and measurement systems. The phase velocity is the speed of the wave propagating at a constant phase, which is a solution of :
The complex wave number
which describes the relationship between stress and strain in the Voigt tissue. The Voigt model consists of an elastic spring μ1 and a viscous damper μ2 connected in parallel, which represents the same strain in each component as shown in Figure 1.
The relation between stress σ and strain ε of the Maxwell tissue is:
For a harmonic motion, (9) becomes:
The elasticity μ1 and viscosity μ2 are two constants and independent to the frequency.
A numerical example of phase velocity of Voigt tissue is shown in Figure 2. Equation (11) shows that cs(ω) increases at the rate of square root of the frequency and there is no the upper limit for cs(ω). As shown in the Figure 2, the phase velocity is determined by both elasticity and viscosity. Ignoring the viscosity introduces errors and biases for elasticity estimates. However, examining the velocities at the extreme frequencies is useful for understanding the model and obtaining initial values for numerical solutions of μ1 and μ2. In tissue characterization applications, μ1 is often in the order of a few thousands and μ2 is often less than 10. Thus, when the wave frequency is very low (less than a few Hz),
When the frequency is very high (higher than a few tens of kHz),
A broad frequency range is needed to accurately estimate both μ1 and μ2. (12) and (13) are only useful for estimating initial values for the numerical solutions of (11) with measured velocities, and they should not be used for final estimates.
Equation (7) can be used for other models for the plane shear wave having a single frequency. The Maxwell model consists of a viscous damper η and an elastic spring E connected in series, which represents the same stress in each component, as shown in Figure 3.
The relation between stress σ and strain ε of the Maxwell tissue is:
For a harmonic motion, (14) becomes:
which is the complex shear modulus of the Maxwell model. Unlike the Voigt model, real and imaginary components of (15) are functions of the frequency. When the frequency is fixed, the complex modulus is a function of and E. Substituting (15) into (7), the shear wave speed in Maxwell medium can be found from (6):
A numerical example of phase velocity of Maxwell tissue is shown in Figure 4. Note that cs(ω) gradually increases to a limit that is proportional to the square root of the elasticity. As shown in the Figure 4, the phase velocity is determined by both elasticity and viscosity. However, examining the velocities at the extreme frequencies is useful for understanding the model and obtaining initial values for numerical solutions of E and η. for a very large ω, for a very small ω,
The Zener model adds an additional elastic spring, having the elasticity of
Equation (18) shows that for a very large ω, for a very small ω,
3. Ultrasound vibrometry
Ultrasound vibrometry has been developed to induce shear wave in a tissue region, measure phase velocity of the shear wave, and calculate the tissue viscoelasticity based on (11), or (16), or (18). The basics of the ultrasound vibrometry are described in details in references [11-17, 32]. Ultrasound vibrometry induces tissue vibrations and shear waves using ultrasound radiation force and detects the phase velocity of the shear wave propagation using pulse-echo ultrasound.
From the solution of the wave equation, equation (5) can be represented by a harmonic motion at a location,
where T is the period of the push pulses shown in Figure 9 and the modulation index is:
Operating on the in-phase and quadrature components
A phase constant can be added to the local oscillator of the demodulator  to avoid zeros in I. The signal extracted by (23) is proportional to the displacement of a harmonic motion induced by the push pulses.
Thus, the motion velocity in slow time can be obtained,
which is proportional to the velocity of the tissue harmonic motion for
The slow-time signal s(t,k) represents the tissue motion at a particular location, its amplitudes and phases change over distances are described by (5). The measurements of amplitudes and phases at two locations are used to calculate attenuation and phase velocity. As shown in (1), the phase velocity is related to the frequency and inversely related to the phase difference Δ
There are several methods to estimate the phases of slow-time signals: Fourier transform, correlation method, and Kalman filter . The estimated phase of the slow-time signal at a location include some phase constants due to the tissue location t and different pulse k, and phase
Ultrasound vibrometry is developed to induce the shear wave described by (19) and detect the phase shift described by (26) for characterizing the tissue shear property using (1) and (11), (14), and (16). Ultrasound virbometry uses interleaved periodic pulses to induce shear wave and detects the phase velocity of shear wave propagation using pulse-echo ultrasound. Figure 8 shows an application setup of the ultrasound vibrometry. An ultrasound transducer transmits push beams to a tissue region to induce vibrations and shear waves. The push beams are periodic pulses that have a fundamental frequency
There are different variations of the excitation pulses beside the on-off binary pulses: continuous waves , non-uniform binary pulses , and composed pulses or Orthogonal Frequency Ultrasound Vibrometry (OFUV) pulses [30, 31]. The OFUV pulses can be designed to enhance higher harmonics to compensate the high attenuations of high harmonics. The OFUV pulses have multiple binary pulses in one period of the fundamental period [30, 31]. Other variations of the ultrasound vibrometry include consideration of background motion and boundary conditions that require more complicated models of tissue motions  and wave propagations .
4. Finite element simulation of shear wave propagation
Simulations using Finite Element Method (FEM) were conducted to understand the shear wave propagation in tissue. The simulation tool is COMSOL 4.2. The simulated tissue region is a two-dimensional axisymmetric finite element model of a viscoelastic solid with a dimension of 100 mm × 100 mm, as shown in Figure 10. The size of domain Ω1 is 100 mm × 80 mm. The domain is divided to 25,371 mesh elements and the average distance between adjacent nodes is 0.95 mm. The schematic diagram shown in Figure 10 includes simulation domains (Ω1, Ω2, Ω3) and boundaries (B1,B2). A line source (with a length of 60 mm) in the left of the solid represents as an excitation source of the shear wave.
All domains had the same material property of the Voigt tissue and all boundaries were set free to avoid reflections. The material parameters were: density of 1055 kg/m3, Poisson’s ratio of 0.499, and Voigt rheological model of the viscoelasticity model. The Voigt model was converted and represented in the form of Prony series. The store modulus and loss modulus were calculated using frequency response analysis for demonstrating the conversion of the Prony series. The complex shear modulus of the Voigt model is the same as (8):
where elasticity modulus
Transient analysis was used and the time step for solver was one eightieth of the time period of the shear wave. Uniform plane shear wave was produced by oscillating the line source with ten cycles of harmonic vibrations in the frequency range from 100 Hz to 400 Hz with a maximal displacement in the order of tens of micrometers. The displacements of the shear wave were recorded for post-processing at 8 locations, 1 mm apart, along a straight line that is normal to the line source. The phases of the wave were estimated by the Kalman filter and the average phase shifts were estimated using a linear fitting method . The estimates of shear wave velocity and viscoelasticity are shown in Table 1.
|Shear Wave Velocity (m/s)||Viscoelasitcity Estimation|
The shear wave velocities in red represent the theoretical values of wave speeds in Voigt tissue. The estimates of the speeds and viscoelasticity moduli of three simulations are shown by three sets of the measurement. Their average values are close to the theoretical values as shown in Figure 11, except the elasticity
5. Experiment system and results
Experiments were conducted for evaluating ultrasound vibrometry. The diagram of an experiment system is shown in Figure 12. This system mainly consists of a transmitter to produce the ultrasound radiation force and a receiver unit using a SonixRP system. Two arbitrary signal generators were utilized to generate the system timing and excitation waveform. The waveform was amplified by a power amplifier having a gain of 50 dB to drive an excitation transducer for inducing vibrations in a tissue region. The SonixRP system was applied to detect the vibration using pulse-echo mode with a linear array probe. The SonixRP is a diagnostic ultrasound system packaged with an Ultrasound Research Interface (URI). It has some special research tools which allow users to perform flexible tasks such as low-level ultrasound beam sequencing and control. The center frequency of the excitation transducer was 1 MHz. The center frequency of the linear array probe was 5 MHz and the sampling frequency of SonixRP was 40 MHz. The excitation transducer and detection transducer were fixed on multi-degree adjustable brackets and were controlled by three-axis motion stages.
The picture of experiment system setup is shown in Figure 13. The left lobe of a SD rat liver was embedded in gel phantom and placed in water tank. Before experiment, the SonixRP URI was run first to preview the internal structure of the liver. In the interface shown in Figure 14, the B-mode image and RF signal of a selected scan line were displayed together to help users selecting test points inside the liver tissue. The positions of the excitation transducer and the detection probe were adjusted to focus on two locations in the liver at the same vertical depth.
Computer programs based on the software development kit (SDK) of SonixRP were developed for detecting the vibrations and shear wave propagation. The programs defined a specific detection sequencing and timing that repeatedly transmit pulses to a single scan line and repeatedly receive the echoes with a PRF of 2 KHz. The timing of the excitation and detection pulses is shown in Figure 15. The pulse repetition frequency of the excitation pulses was 100 Hz.
An example of the typical fast-time RF ultrasound signal acquired by the SonixRP is shown in Figure 16. Figure 16a shows the echo through the entire liver tissue region, while Figure 16b shows the echo around the focus point (75 mm in depth) in the liver tissue.
The vibration of shear wave at a location was extracted from I and Q channels using the I/Q estimation algorithm described by equation (23). Figure 17a shows the vibration displacement and Figure 17b shows the spectral amplitude of the vibration.
The extracted displacement signal
Figure 19 shows the phase velocities at different harmonics and the fitting curves of three models: Voigt, Maxwell, and Zener models. The fitting values are shown in Table 2. As shown by the figure and table, the Voigt model and Zener model fit the measurements of the phase velocity of the liver tissue better than the Maxwell model for this liver.
|Voigt Model, μ1, μ2, fitting error||4.10 kPa||1.51 Pa·s||0.019 m/s|
|Maxwell Model, E, η, fitting error||7.18 kPa||4.27 Pa·s||0.143 m/s|
|Zener Model, E1, E2, η, fitting error||4.07 kPa, 45.9 kPa||1.47 Pa·s||0.020 m/s|
The second experiment was conducted to demonstrate the impact of boundary conditions. Because boundary conditions play very important roles in wave propagation,
The shear wave speeds were measured from 100 Hz to 800 Hz over a distance up to 5 mm away from the center of the radiation force application. Figure 20 shows the estimates of the shear wave speeds. Each error bar was the standard deviation of 30 estimates from five data sets of repeated measurements and six distances (1 to 4 mm, 1 to 5 mm, etc). The estimates from 100 Hz to 400 Hz were almost identical for the binary excitation pulses and the OFUV excitation pulses. Because the estimate errors using binary excitation pulses were too high for the frequency beyond 400 Hz, the estimates at 4.9 mm were based on the OFUV method. Figure 20 represents the trend of our experiment results that the shear wave speed in the superficial liver tissue is generally higher than that in the deep tissue. The results should be carefully examined. One of the possibilities is that we think it is caused by the liver capsule as we have verified it with Finite Element (FE) simulations, and another possibility is that the shear wave speeds of the gelatin are between 3 to 4 m/s from 100 to 800 Hz, higher than that in the liver tissue.
The estimates of shear wave speeds at deep tissue of 4.9 mm and superfical tissue of 0.4 mm were used to numerically solve for the shear moduli of the three models. The curves generated by the models were compared with the measurements. As shown in Figures 21a and 21b, we find that the Voigt model may not always suitable for modeling liver shear viscoelasticity, at least for
Table 3 shows the estimated shear moduli of different models with two different frequency ranges at two different depths in liver tissues based on our experiment data. Each modulus is an average of 30 estimates from 5 data sets and 6 distances. All elasticity has the unit of kPa and all viscosity has the unit of Pa·s. The fitting errors (m/s) are the deifferences between the measurements and calculated shear wave speeds using the models. The changes represent the variations of the estiamtes from one frequency range to another. The statistics are not conclusive because of the small number of samples. But this study indicates the variations of estimates and importance of the selection of tissue viscoelasticity models.
|2.48, 2.00, 0.152||3.71, 1.46, 0.204||50%, 27%|
|10.7, 2.50, 0.043||11.7, 2.36, 0.048||10%, 6%|
|0.578, 9.033, 2.85, 0.028||1.34, 9.843, 2.56, 0.0569||132%, 9%, 10%|
|2.74, 1.35, 0.108||3.59, 0.791, 0.151||31%, 41%|
|5.68, 2.82, 0.016||5.90, 2.70, 0.021||4%, 4%|
|1.49, 4.20, 2.44, 0.015||1.70, 4.25, 2.19, 0.018||14%, 1%, 10%|
The third experiment was conducted to demonstrate the effectiveness of the ultrasound vibrometry to characterize the injury of liver tissue. Table 4 shows that the measured shear moduli of the livers thermally damaged by a microwave oven using different amount of cooking time (3, 6, 9, and 12 seconds). All estimates were from the superficial tissue region. It shows that the shear wave speeds estimated in the superficial tissue region are effective for indicating the damage levels of the livers. The errors are the standard deviations of the differences between the measurements and calculated speeds of the models. The Zener model provides the best curve fitting with the minimum fitting error.
|3 sec.||6 sec.||9 sec.||12 sec.|
Shear moduli have very high dynamic ranges and are highly correlated with the pathological statues of human tissue. The solutions of the wave equation with constitutional models of tissue viscoelasticity show that the shear moduli of tissue can be estimated by measuring the phase velocity and attenuation of shear wave propagation in the tissue. However, it is a challenge to effectively and remotely generate vibrations and shear waves in a tissue region. It is also a challenge to measure shear wave because shear wave attenuates very fast as the propagation distance increases.
In the past fifteen years, the use of pulsed and focused ultrasound beams has been demonstrated as an effective method to remotely induce localized vibrations and shear waves in a tissue region. Several useful technologies have been developed for characterizing tissue viscoelasticity: Vibroacoutography, ARFI, Supersonic imaging, and ultrasound vibrometry, etc.
The ultrasound vibrometry is only technique that quantitatively estimates both tissue elasticity and viscosity. We found that the estimates of tissue elasticity by ignoring the viscosity are erroneous. Shear phase velocity are frequency dependent because the dispersive property of the biological tissue. Therefore, regardless of the usefulness of the viscosity, accurate estimates of tissue elasticity require the inclusion of the viscosity in the tissue models, as indicated by the solutions of the wave equation with three viscoelasticity models.
The ultrasound vibrometry transmits periodic push pulses to induce vibrations and shear waves in a tissue region, and detects the shear wave propagation using the pulse-echo ultrasound. The push pulses and detection pulses are interleaved so that one array transducer can be used for the applications of both pulses. The application of the array transducer allows the detection over a distance so that the phase velocities of several harmonics can be measured for calculating shear moduli.
Accurate estimates of shear moduli require an extended frequency range over an extended distance. The current technology is only effective for a few hundred Hertz in the frequency and a few mm in the distance away from the center of the radiation force applied. Shear wave having a high frequency attenuates very quickly as distance increases. Other vibration methods such as OFUV may be worth to explore.
We found that the shear wave speeds of livers are location dependent or dispersive in locations. Our experiment results indicate that the shear moduli estimated from a superficial tissue region and from a deep tissue region can be significantly different. Boundary conditions play a very important role in shear wave propagation and its phase velocity. The solution of the wave equation with boundary conditions should be considered for a tissue region that has a limited physical size. Some studies in this area have been done for myocardium and blood vessel walls.
The measurements of the ultrasound vibrometry are based on the assumption that tissue under the test is isotropic, which is not true for most tissues. Nevertheless, the measurements may be useful in clinical practices, which need to be evaluated
Limited by the extensive contents in this chapter, we do not discuss the application of the Kalman filter in this work. The Kalman filter has great potential to include more complicated tissue models and motion models that are not fully explored yet, at least are not publically reported yet. On the other hand, Fourier transform and correlation method are also effective tools to calculate phases of the slow-time signals, if the motion model is simply sinusoidal.
Our experiments demonstrate that the ultrasound vibrometry can be readily implemented by using commercial medical ultrasound scanners with minimum alterations. Our experiment results also demonstrate that the ultrasound vibrometry is effective to characterize the stiffness and injury levels of livers.
We find that the Zener model fit the shear wave speeds of the livers better than the Voigt model and Maxwell model in almost all cases that include different frequency ranges, different locations, and different tissue conditions. Our study also indicates that the Voigt model is sensitive to the change of the observation frequency. Measurements at higher frequencies should be included when the Voigt model is used. In this case, the OFUV is useful to enhance the higher frequency components of the shear waves. The Zener model and Maxwell model appear to be less impacted by the frequency changes with our experiment data.
Tissue pathological statues are related to tissue shear moduli, which can be estimated by measuring the phase velocity of shear wave propagation in a tissue region. Ultrasound vibrometry is an effective tool to quantitatively measure tissue elasticity and viscosity. Ultrasound vibrometry induces vibrations in a tissue region using pulsed and focused ultrasound radiation force and detects the shear wave propagation using pulse-echo ultrasound. Experiment results demonstrate the effectiveness of the ultrasound vibrometry for characterizing tissue stiffness and liver damages.
This research was supported in part through a grant from National Institute of Health (NIH) of USA with a grant number of EB002167 and a grant from Natural Science Foundation of China (NSFC) with a grant number of 61031003.