|Home | About | Journals | Submit | Contact Us | Français|
The speed at which shear waves propagate in tissue can be used to quantify the shear modulus of the tissue. As many groups have shown, shear waves can be generated within tissues using focused, impulsive, acoustic radiation force excitations, and the resulting displacement response can be ultrasonically tracked through time. The goals of the work herein are two-fold: first, to develop and validate an algorithm to quantify shear wave speed from radiation force-induced, ultrasonically-detected displacement data that is robust in the presence of poor displacement signal-to-noise ratio (SNR), and second, to apply this algorithm to in vivo datasets acquired in human volunteers in order to demonstrate the clinical feasibility of using this method to quantify the shear modulus of liver tissue in longitudinal studies. The ultimate clinical application of this work is non-invasive quantification of liver stiffness in the setting of fibrosis and steatosis.
In the proposed algorithm, time to peak (TTP) displacement data in response to impulsive acoustic radiation force outside the region of excitation (ROE) are used to characterize the shear wave speed of a material, which is used to reconstruct the material’s shear modulus. The algorithm is developed and validated using finite element method (FEM) simulations. Using this algorithm on simulated displacement fields, reconstructions for materials with shear moduli (μ) ranging from 1.3–5 kPa are accurate to within 0.3 kPa, while stiffer shear moduli ranging from 10–16 kPa are accurate to within 1.0 kPa. Ultrasonically tracking the displacement data, which introduces jitter in the displacement estimates, does not impede the use of this algorithm to reconstruct accurate shear moduli.
Using in vivo data acquired intercostally in 20 volunteers with body mass indices (BMI) ranging from normal to obese, liver shear moduli have been reconstructed between 0.9 and 3.0 kPa, with an average precision of ±0.4 kPa. These reconstructed liver moduli are consistent with those reported in the literature (μ = 0.75–2.5 kPa) with a similar precision (±0.3 kPa). Repeated intercostal liver shear modulus reconstructions were performed on 9 different days in 2 volunteers over a 105 day period, yielding an average shear modulus of 1.9 ± 0.50 kPa (1.3–2.5 kPa) in the first volunteer, and 1.8 ± 0.4 kPa (1.1–3.0 kPa) in the second volunteer. The simulation and in vivo data to date demonstrate that this method is capable of generating accurate and repeatable liver stiffness measurements and appears promising as a clinical tool for quantifying liver stiffness.
Impulsive acoustic radiation force excitations can be generated in tissue at remote, focused locations, and the resulting dynamic tissue response can be monitored using ultrasonic displacement tracking methods. The rate at which tissue responds to an impulsive excitation, including the speed at which shear waves propagate away from the region of excitation (ROE), can be measured to quantify the tissue’s shear modulus, as originally proposed by Sarvazyan et al. (1998). In contrast to elastographic strain images (Ophir et al., 1991; Hall et al., 2003; Greenleaf et al., 2003) and Acoustic Radiation Force Impulse (ARFI) images (Nightingale et al., 2006) that show relative structural stiffness compared with adjacent tissues, the ability to quantify an absolute tissue modulus will be useful for different clinical applications. This approach will allow disease processes that involve the stiffening or softening of tissue without large scale structural changes, such as liver fibrosis and steatosis, to be monitored longitudinally to determine when to initiate/cease treatment protocols and to stage disease progression and resolution.
One of the fundamental difficulties with reconstructing shear moduli from shear wave speeds is generating shear waves in vivo. Systems that use external mechanical excitation (static or dynamic) are challenged in coupling the excitation into the organ/structure of interest, especially if the tissue is deep within the body (Sandrin et al., 2003). The use of focused acoustic energy circumvents this challenge by providing mechanical excitation directly to the focal region of the acoustic beam and generating shear waves directly in the tissue of interest. Impulsive acoustic radiation force excitations generate shear waves within tissues (Sarvazyan et al., 1998). The speed at which these shear waves propagate away from the ROE is related to the shear modulus and density of the tissue; therefore, measuring this shear wave speed facilitates estimation of the tissue’s shear modulus.
There are several methods available to measure the speed of shear waves using dynamic displacement data. Inversion of the Helmholtz equation to reconstruct shear wave speed from displacement data has been used in ARFI and supersonic imaging (Nightingale et al., 2003; Bercoff et al., 2004). One drawback to this reconstruction method is having to perform second-order differentiation of displacement data in space and time. Jitter associated with ultrasonically tracking these displacement fields (Walker and Trahey, 1995; Pinton and Trahey, 2006; Palmeri et al., 2006) necessitates that significant filtering operations be performed on the displacement data. These filtering operations are computationally intensive and better suited for offline, rather than real-time, processing. An alternative approach to estimating shear wave speed is to use time-of-flight measurements, where the shear wave position is characterized as a function of time. Shear waves can be tracked using correlation-based algorithms (McLaughlin and Renzi, 2006a,b), or as described herein, times to peak (TTP) displacement outside the ROE can be used to estimate shear wave speed.
This manuscript presents an imaging system capable of generating and monitoring radiation force-induced shear waves in human liver in vivo, along with a robust algorithm for reconstructing shear wave speeds from ultrasonically detected displacements to quantify shear moduli. The Background section presents the work of Sarvazyan et al. (1998) that motivated the development of this algorithm, in addition to discussing the general mechanics surrounding shear wave propagation and the use of Helmholtz and Eikonal methods to reconstruct shear wave speed. The Methods section outlines the implementation of the new algorithm, along with describing the simulation and experimental setups used to quantify the accuracy and precision of the algorithm in the clinical context of quantifying liver stiffness. The Results sections shows the algorithm applied to simulation data, experimental phantom data, and in vivo human liver data. The human studies were performed in 20 volunteers, with the repeatability of this shear modulus reconstruction approach studied in two volunteers over a 105 day period. Additional approaches to optimize this algorithm, along with this algorithm’s limitations, are explored in the Discussion section.
Liver disease is among the ten major causes of death in the United States (Foundation, 2006). Chronic hepatitis and cirrhosis are diseases that progress over several decades, and hepatic fibrosis staging is the key factor in determining liver health for the majority of liver diseases. Fibrosis staging is currently accomplished by a single needle core biopsy, typically performed without image guidance. Liver biopsies are typically not well tolerated, can be associated with complications, and are thus generally performed only for initial diagnosis and clinical treatment endpoints. While core biopsy is considered to be the current gold standard for evaluating liver health, its accuracy is limited by a small sample size, with misdiagnosis in fibrosis staging reported in 20–40% of cases (Ratziu et al., 2005; Regev et al., 2002; Bedossa et al., 2003). There are several treatments being studied to arrest or reverse liver fibrosis; however, these developments are hindered by the lack of an inexpensive, noninvasive method for monitoring hepatic fibrosis. There are strong and immediate clinical needs to develop an accurate and efficient method to noninvasively and longitudinally characterize liver stiffness in an outpatient setting, and there is an immediate need for a noninvasive imaging modality capable of monitoring disease progression during the course of treatment trials.
In linear, isotropic, elastic solids, the speed of shear wave propagation (cT ) is related to shear modulus (μ) and density (ρ) by:
Equation 1 provides a relationship between shear wave speed and shear modulus; however, there are two significant challenges to using this relationship to characterize the modulus of soft tissue: (1) generating shear waves within tissues in vivo, and (2) reconstructing cT from measured displacement fields.
Generating shear waves within tissues can be accomplished by coupling external mechanical sources through the skin into the organ of interest or generating the shear wave within tissues using acoustic radiation force. The FibroScan® system (EchoSens, Paris, France) uses an external vibrator to generate shear waves in tissue and has successfully quantified differences in liver stiffness as correlated with fibrosis stage (Sandrin et al., 2003). Similar approaches of external shear wave excitation have also been used in MR-based elastography techniques (Huwart et al., 2006; Roiviere et al., 2006). While these findings are promising, such setups can be challenged in their ability to couple enough energy through the skin and subcutaneous fat to generate adequate shear wave displacements within organs such as the liver, especially in obese patients. External mechanical excitation sources can also be limited by the ribs when trying to reach more superior and lateral regions of the liver.
Some of these challenges can be overcome with the use of focused acoustic radiation force excitations where mechanical excitation occurs along the acoustic wave propagation path and within the focal region of the acoustic beam. These radiation force excitations generate shear waves directly in the tissue of interest. Acoustic radiation force is applied to absorbing and/or reflecting materials in the propagation path of an acoustic wave. This phenomenon is caused by a transfer of momentum from the acoustic wave to the propagation medium. The spatial distribution of the radiation force field (i.e., the region of excitation, or ROE) is determined by both the acoustic excitation parameters and the tissue properties. In soft tissues, where the majority of attenuation results from absorption (Parker, 1983; Christensen, 1988), the following equation can be used to determine radiation force magnitude (Torr, 1984; Nyborg, 1965):
where F [dyn (1000 cm)−3], or [kg s−2 cm−2], is acoustic radiation force (in the form of a body force), Wabsorbed [W (100 cm)−3] is the power absorbed by the medium at a given spatial location, α[m s−1] is the speed of sound in the medium, a [cm−1] is the absorption coefficient of the medium, and I [W cm−2] is the temporal average intensity at a given spatial location. The spatial extent of the ROE varies with focal characteristics and tissue attenuation; however, it is always distributed within the geometric shadow of the active transmit aperture and is typically most energetic within the focal region of the acoustic beam.
After shear waves are generated or coupled into the organ of interest, and the tissue displacement response has been measured, there are multiple methods that can be used to estimate the shear wave speed. One method to estimate shear wave speed from dynamic tissue displacement data () involves algebraic inversion of the Helmholtz equation (Oliphant et al., 2001; Bercoff et al., 2004):
This method has been successfully applied when ultrasound or magnetic resonance imaging has been used to track tissue displacements (Oliphant et al., 2001; Bercoff et al., 2004; Sandrin et al., 2002; Nightingale et al., 2003), but as Equation 3 indicates, second-order spatial and temporal derivatives of displacement are required for the shear wave speed reconstruction. Given the jitter that exists when ultrasonically estimating displacement, appreciable filtering and smoothing of displacement data must be performed prior to processing (Oliphant et al., 2001; Bercoff et al., 2004). The benefit of this method is that no a priori assumption about shear wave propagation direction needs to be made when analyzing a given region of dynamic displacement data.
Ideally, shear wave speeds would be reconstructed from high signal-to-noise ratio (SNR), three dimensional displacement data using inversion of the Helmholtz equation, allowing for good spatial resolution. However, in our experience, the typically low SNR displacement data obtained with ultrasonic displacement tracking in a single imaging plane yield inaccurate results when relying on second-order differentiation of displacement data. This has motivated the development of alternate methods to reconstruct shear wave speeds in response to impulsive radiation force excitations.
Time-of-flight methods track the position of shear waves through time and correlate their space/time coordinates to estimate shear wave speeds. McLaughlin et al. (2006a) have implemented such an approach using correlation methods on displacement datasets to determine the position of the shear wave, and shear wave speeds are estimated by inverting Eikonal equations that rely on first-order differentiation of shear wave positions through time (McLaughlin and Renzi, 2006a,b).
The Lateral TTP algorithm developed herein is a time-of-flight method. In order to make the algorithm robust in the presence of noise, the following assumptions are made: (1) homogeneity of the region adjacent to the ROE, (2) shear wave propagation exclusively in the lateral direction (perpendicular to the ROE axis of symmetry), and (3) negligible dispersion over the analyzed region.
In this algorithm (described in detail in the Methods Section, referencing Figures 1–4), shear wave position is estimated by quantifying the TTP displacement at laterally-offset positions outside of the ROE. Calculations are performed over the depth of field (DOF) of the focused radiation force excitation to satisfy the assumption of shear wave propagation in the lateral direction. Using TTP displacements to estimate shear wave speed assumes that the shear wave’s peak displacement propagates at the shear wave group velocity through the analyzed region, which is the case for purely elastic or mildly dispersive media. Linear regressions are then performed on the TTP displacement data versus lateral position. The R2 value of the linear regression and the associated 95% confidence interval are used as goodness of fit metrics to exclude datasets corrupted by motion, noise, or other artifacts. The inverse slopes of the remaining lines represent the shear wave speeds and can be used to estimate the shear moduli of the material as a function of depth. More details about the implementation of this algorithm will follow in the Methods section.
The Lateral TTP algorithm (Figures 1–4) was applied to radiation force-generated, ultrasonically-tracked axial displacement data monitored through time at laterally-offset locations in the imaging plane relative to a fixed excitation location. In order to satisfy the assumption of uniform shear wave propagation parallel to the lateral dimension, the axial extent of the data utilized to estimate the shear wave speed was confined to be within the depth of field (DOF) of the excitation beam, as demonstrated in Figure 1. The DOF was defined by 8(F/#)2λ, where F/# and λ represent the focal configuration and wavelength of the excitation beam, respectively. The dimensionless excitation beam f-number (F/#) was defined as , where z was the excitation focal depth, and d was the electronically active aperture width of the excitation beam. At each lateral location, the DOF was subdivided into 0.5 mm increments for analysis, with the displacement averaged over ±0.25 mm at each depth to reduce jitter/noise. No axial averaging was performed on FEM model datasets without simulated tracking. This DOF is indicated by the horizontal dashed lines in the simulated TTP displacement data in Figure 1. The TTP displacement was estimated from these displacement through time data after upsampling to 50 kHz using low-pass interpolation.
The rate at which TTP displacement changes with lateral position was evaluated using linear regression. Regressions were performed starting one excitation beamwidth from the center of the ROE and extended over a lateral range where the peak displacements remained above 1 μm. The 1 μm threshold was chosen to ensure a peak displacement estimate SNR of at least 2 dB, as demonstrated in Figure 2. The inverse slopes of these regression lines, with goodness of fit metrics exceeding a threshold (R2 > 0.8, 95% CI < 0.2), represent the material’s local shear wave speeds. These specific goodness of fit metrics were applied to all of the datasets presented throughout this manuscript. The interrogated material’s density was assumed to be 1.0 g/cm3, and the material’s shear modulus was then estimated using Equation 1. Young’s moduli (E) could be estimated by E = 2(1 + ν) μ = 3μ, where ν = 0.5 is an incompressible material’s Poisson’s ratio, but only values for shear modulus (μ) are quoted herein.
This procedure is graphically demonstrated in Figure 3 using the simulation data sampled at 10 kHz. (Note that these data would be upsampled to 50 kHz before being processed with the Lateral TTP algorithm.) When two adjacent time steps happened to yield the same peak displacement values, the smaller time step was chosen.
Three-dimensional Finite Element Method (FEM) models of the dynamic response of elastic media to impulsive acoustic radiation force excitations were used to study the accuracy of the proposed method in reconstructing shear moduli ranging from 1.3–16 kPa. These shear moduli represent those reported for healthy through cirrhotic livers (Foucher et al., 2006; Sandrin et al., 2003). These models have been previously validated to accurately simulate shear waves that are generated in response to impulsive acoustic radiation force excitations in elastic media (Palmeri et al., 2005). Table 1 outlines the simulated excitation beam configuration.
The impact of ultrasonically tracking the axial components of simulated displacement fields was also characterized using previously validated methods (Palmeri et al., 2006). These simulations were performed under noise- and physiologic-motion-free conditions in purely elastic media and were used to characterize the accuracy and precision of the proposed Lateral TTP algorithm. The tracking beam configuration that was simulated is outlined in Table 1.
A modified Siemens SONOLINE™ Antares scanner (Siemens Medical Solutions USA, Inc., Ultrasound Division, Issaquah, WA, USA) was utilized for all experiments. Different transducers and system parameters were used for each experiment, depending upon the depth of the region of interest (ROI), and are specified in the experiment-specific sections that follow.
For all experiments, custom beam sequences were programmed into the scanner and either radio-frequency (RF) or the real and imaginary (IQ) components of the RF tracking data were stored for offline processing using custom written algorithms (described below) on a Linux Beowulf cluster. IQ data were acquired using 4:1 parallel receive mode, where 4 receive beams are acquired for each tracking transmit beam (Dahl et al., 2007), while RF data were acquired in conventional receive mode (1:1) for the gelatin phantom study. Each interrogation consisted of a reference tracking pulse (conventional B-mode pulse) followed by high intensity pushing pulse and a series of tracking pulses to track the displacement and full recovery of the material after excitation. Further details of this procedure are covered by Dahl et al. (Dahl et al., 2007).
The size of the ROI for estimating shear wave speed is large relative to the ROE and cannot be adequately characterized with only 4 tracking beams after a single excitation. To fully characterize the shear wave propagating across the ROI, nine interrogations of the reference:push:tracking sequence were fired, where the excitation remained in the same location (at a lateral position of 0 in all of the images herein), but the tracking beams were offset at greater lateral positions from the excitation with each repetition. Displacements were estimated using either Loupas’ method on IQ data or normalized cross-correlation with a 1.5 λ kernel with 99% overlap on the RF data, as detailed by Pinton et al. (Pinton et al., 2006).
For all in vivo datasets, a linear motion filter with a temporal range that varied as a function of lateral offset from the ROE was applied to remove physiologic motion from the displacement data. The temporal range was chosen to coincide with when a shear wave was expected to travel through a given lateral position (xlat) in the ROI. For each xlat, the earliest time for the first non-zero displacement due to the propagating shear wave (tstartMin) was computed as:
where Texc was the duration of the excitation, 3 * BW was three times the −6 dB lateral beamwidth of the excitation beam in the DOF that approximates the spatial extent of the propagating shear wave, cTmax was the maximum shear wave speed that is expected in the material being imaged, and cTmin was the minimum expected shear wave speed. Three times the −6 dB lateral beamwidth of the excitation beam was empirically chosen as a conservative estimate to estimate > 95% of the shear wave’s spatial extent. The latest time that a shear wave was expected to pass through a given lateral position (tstopMax) was computed as:
where TswMax was the maximum expected shear wave period, and tstartMax was the latest time when a shear wave would have started passing through xlat The residual displacement was measured at tstopMax, when the tissue should have fully recovered from the excitation. This residual displacement at tstopMax was used to subtract a linear displacement artifact through time based on the displacement at tstartMin, before which no displacement due to the propagating shear wave is expected. This filter removes transducer and physiologic motion, in addition to low frequency artifacts that can arise from transducer heating and power supply variations.
Phantom measurements were performed using a gelatin-based, tissue-mimicking phantom (Hall T.J. et al., 1997) using the VF10-5 array and a calibrated tissue-mimicking phantom (Computerized Imaging Reference Systems, Inc., Norfolk, VA) using the PH4-1 and VF10-5 array setups, using the parameters detailed in Table 2. These studies were used to characterize the precision of the algorithm empirically and to demonstrate the algorithm’s independence on the array/focal configuration used to generate and track the shear waves.
The feasibility of reconstructing hepatic shear moduli with the Lateral TTP algorithm was also demonstrated in human volunteers with written consent under approval from the Duke University Medical Center IRB (#9328-06-12). An intercostal imaging approach was used in all volunteers, with imaging performed between the ninth and tenth ribs, slightly anterior to the mid-axillary line. Six independent measurements were made during an imaging session, where a measurement consists of the the transducer being placed in the intercostal space and the patient performing a full inspiration breath hold. Electrocardiogram (ECG) triggering was not utilized in this study, but could be added in future studies to potentially reduce cardiac motion artifacts. One of two individuals performed all of the volunteer scanning in these studies.
Data were acquired using the procedure outlined above with the PH4-1 array with the Antares scanner, using the imaging parameters outlined in Table 2. For each independent measurement, a total of nine excitation pulses were utilized, allowing monitoring of shear wave propagation to the right over a 6 mm lateral region, with the excitation beam centered at 0 mm laterally. The radiation force excitation power was chosen to achieve adequate displacement magnitudes (10–20 μm) within the ROE, while minimizing tissue heating. Based upon hydrophone measurements in water, using linear extrapolation of small signal derated fields (α = 0.7 dB/cm/MHz) (NCRP, 2002), the in situ spatial peak pulse average intensity of the excitation pulses is estimated to be 1470 W/cm2. For these shear wave sequences (nine repeated excitations in a single location), validated FEM simulations indicate that the total cumulative temperature rise in liver associated with this sequence occurs at the focal point, and is 0.25°C for less than 0.2 seconds of tissue insonification without taking into account perfusion effects that would lower heating estimates (Palmeri et al., 2004). This heating is less than the 6°C of heating that is accepted for diagnostic ultrasonic imaging per the US FDA (NCRP, 2002).
Displacement tracking was performed using 4:1 parallel receive, where 4 receive beams were acquired for each tracking beam (Dahl et al., 2007). Displacements were estimated using the Loupas algorithm on IQ data (Pinton et al., 2006) and motion filtering was performed as outlined in the earlier Motion Filtering section.
Simulation data were used to characterize the accuracy and precision of the Lateral TTP algorithm for elastic materials with known shear moduli that were not corrupted by noise and/or motion artifacts. Figure 4(a) shows the TTP displacement estimates as a function of lateral position extending away from the ROE at the excitation focal depth of 20 mm in elastic media with shear moduli of 1.33 and 2.83 kPa. The tracked simulation data were compiled over 20 independent speckle realizations. (The complete configurations of the simulated excitation and tracking beams are provided in Table 1.) Figure 4(b) shows the reconstructed shear moduli at depths within the DOF of the excitation beam. As indicated on the figure, both the non-tracked and tracked FEM data can be reconstructed within ±0.03 kPa of the 1.33 kPa material and ±0.08 kPa of the 2.83 kPa material for this demonstrative speckle realization.
The modulus reconstructions using the simulation data were then extended to a greater range of moduli that may be encountered when characterizing diseased (fibrotic) livers (Foucher et al., 2006; Sandrin et al., 2003). The shear modulus reconstructions for both raw FEM displacement data and ultrasonically-tracked FEM displacement data are shown in Figure 5. The error bars associated with the non-tracked FEM displacement data represent variations in the reconstructed moduli over the DOF, while the error bars associated with the tracked FEM data represent variations in the reconstructed moduli over 20 independent speckle realizations at the focal depth (20 mm).
The Lateral TTP algorithm was applied to ARFI shear wave data obtained in a gelatin tissue-mimicking phantom (Hall T.J. et al., 1997; Palmeri et al., 2005). Two-dimensional (2D) images of reconstructed shear moduli, as shown in Figure 6, were made using the Lateral TTP algorithm. The increased spatial resolution of the Lateral TTP algorithm was achieved by reducing the lateral domain over which the linear regressions were performed on the TTP versus lateral position data. This implementation of the Lateral TTP algorithm is similar to what would be achieved with differentiation of such data in Eikonal methods (McLaughlin and Renzi, 2006a,b), if the lateral extent of the linear regression was reduced to fewer points, approaching a local spatial derivative.
To demonstrate the independence of the Lateral TTP algorithm on the excitation focal configuration and associated tracking configuration, a single location in a calibrated CIRS tissue-mimicking phantom was characterized 12 times using both the PH4-1 and the VF10-5 arrays. Both arrays were positioned so their lateral foci (F/2 focal configuration) occurred at the same phantom location, at a depth of 20 mm from the phantom surface. The VF10-5, which has an elevation focus near 19 mm, was operating at 5.7 MHz with a lateral focus at 20 mm. The PH4-1, which has an elevation focus near 70 mm, was operated at 2.2 MHz with a lateral focus at 37.5 mm (offset by a 17.5 mm water path standoff to image the same location as the more shallowly focused VF10-5 array). The PH4-1 characterization yielded a reconstructed shear modulus of 1.7 ± 0.2 kPa, while the VF10-5 yielded a reconstructed shear modulus of 1.6 ± 0.1 kPa.
While the simulations and phantoms are good controls to test the accuracy and precision of the Lateral TTP algorithm in reconstructing moduli, performing such reconstructions in vivo presents additional challenges, such as physiologic motion and the penetration of acoustic energy through skin and fat.
To demonstrate the feasibility of using the Lateral TTP algorithm to longitudinally monitor liver stiffness in humans, two studies were conducted: (1) 20 human volunteers were imaged intercostally to reconstruct the shear moduli of their livers, and (2) two volunteers (7 and 8 from the 20 volunteer study) had their liver stiffness reconstructed 9 times over a 105 day period to evaluate the repeatability of such measurements. An intercostal imaging approach was chosen for these studies because it characterizes the same right posterior lobe of the liver that is typically biopsied by clinicians. Volunteers without known liver disease were chosen for these studies, but no definitive clinical studies (e.g., liver function tests or coagulation studies) were performed to confirm normal liver health. There were 10 male and 10 female volunteers with a mean age of 37 years (25–71 years old) and BMIs ranging from 19.7–30.5 (obese).
Figure 7 shows a demonstrative data set from Volunteer 16 (25 year old female, BMI 20). Figure 7(a) shows the B-mode image with the ROI used for the shear wave speed characterization outlined by the yellow box. The motion-filtered displacement through time data at the focal depth are shown in Figure 7(b); these data represent one of the depth increments processed over the entire DOF of the excitation beam. Figure 7(c,d) shows the times to peak displacement as a function of lateral position for all of the 0.5 mm increments over the DOF of the excitation beam before (c) and after (d) the goodness of fit metrics were applied to the data. Figure 7(e) shows the reconstructed shear moduli over the 6 trials that were performed during the study, where each trial consisted of a new breath hold and repositioning of the transducer in the same intercostal space.
Figure 8(a) shows the reconstructed shear moduli from all data (red) that had regressions exceeding the goodness of fit metrics over the 6 trials (mean ± one standard deviation) in the 20 human volunteers. The mean precision over each of the individual trials is shown in blue for each volunteer. Figure 8(b) shows the repeated shear modulus reconstruction at 9 different time points over a 105 day period in two of the volunteers. As with the 20 volunteer study, each time point in the repeatability study consisted of 6 independent trials. (The absence of an estimate for Volunteer 7 on the first day was because the volunteer was not available for imaging that day, not because of an absence of valid shear modulus reconstructions.) Figure 8(c) shows the reconstructed shear moduli from the 20 volunteers as a function of their BMI with the vertical dashed lines indicated the transitions between normal, overweight (25–30), and obese (> 30) BMIs.
The data presented herein demonstrate the feasibility of using acoustic radiation force imaging methods to non-invasively quantify liver stiffness in vivo. The reconstructed shear moduli shown in Figures 7 and and88 are consistent with those reported for healthy human liver as determined by external excitation methods (e.g., the FibroScan® system (Sandrin et al., 2003; Foucher et al., 2006), assuming a relation of ). In a study by Rouviere et al. using MR elastography of the liver, the shear stiffness in healthy subjects was found to be 2.0 ± 0.3 kPa (Roiviere et al., 2006); a study by Huwart et al., also using MR elastography, found the shear stiffness of liver to be 2.24 ± 0.23 kPa in healthy subjects. The reconstructed hepatic shear moduli from these MR elastography studies agree well with the data in our human volunteers (Figures 7 and and8).8). In Figure 8(a), the error bars (red) over the six trials are greater than the mean standard deviations over the DOF for a single interrogation (blue) because they include the variability between the 6 different interrogations, including differences in the regions of liver that were interrogated. While we are assuming the liver is homogeneous over the ROIs associated with a given interrogation, there may be heterogeneities between different interrogations depending on transducer location and depth of breath hold.
Figure 8(b) shows that for the same two volunteers, reconstructed shear moduli varied between 1–3 kPa at 9 different time points over a 105 day period. There are several factors that may be affecting this variability, including when imaging is performed relative to eating, blood pressure, operator variability, and actual heterogeneity of stiffness within the liver. These factors are being investigated in greater detail in ongoing studies. While the number of obese volunteers imaged in this study was limited, obesity did not pose a problem to performing shear modulus reconstructions in vivo (Figure 8(c)).
The SNR of the displacement estimates used to characterize shear wave propagation in vivo greatly impacts the success of the algorithm (Figure 2). Jitter in the displacement estimates increases with decreasing tracking beam frequency and is generally between 0.2–1 μm for the experimental setup presented herein (Walker and Trahey, 1995). Clearly, the larger the signal (i.e., displacement due to radiation force excitation), the better the SNR (Figure 2). In general, for all datasets, the number of good estimates (R2 > 0.8, 95% CI < 0.2) increases with increasing displacement magnitude (SNR). In our experience, it is advantageous to generate peak focal displacements of at least 10 μm within the ROE and to restrict the lateral range over which linear regressions are performed to displacement magnitudes of at least 1 μm. It has also been observed that over repeated measurements at a given location, the peak focal displacement for a given transmit power can vary considerably, presumably due to differences in transducer coupling and differences in acoustic attenuation and aberration caused by tissue in the propagation path of the acoustic excitation.
Both the displacement magnitude and the thermal response of the tissue to acoustic excitation are linearly related to the in situ intensity of the acoustic pulse and the attenuation of the tissue (Equation 2), which are impacted by many factors, including: transmit beam parameters (frequency, F/#, and power), tissue absorption at the focus, and attenuation and aberration of the acoustic beam through intervening tissue. Limits on the maximum energy included in an interrogating beam arise from both the thermal response (for diagnostic ultrasonic imaging, total tissue temperature increase must be less than 6°C (Herman and Harris, 2002)) and from system power output limitations. In the experimental implementation presented herein, the tissue is repeatedly excited in one spatial location in order to monitor shear wave propagation throughout a large lateral field of view (FOV). The use of 4:1 parallel receive tracking allows four lateral locations to be tracked for each excitation, reducing the total number of excitations by a factor of 4. This parallel receive implementation allows for either a reduction of the heat generated in the tissue or quadrupling the power in each excitation without any additional tissue heating. The tradeoff between lateral FOV size (the region over which shear wave propagation is monitored) and acoustic energy within each excitation (as determined both by tissue heating and system power output considerations) is currently under investigation.
A potential challenge for methods using acoustic radiation force to quantify elastic moduli in human liver in vivo are people who have “poor ultrasonic image quality”. Acoustic waves pass through varying amounts of skin, connective tissue and fat before entering the liver, and such heterogeneous propagation paths can be associated with phase aberration, varying attenuation, and other effects that distort the distribution and decrease the magnitude of the acoustic energy within the focal region. For such cases, a potential solution would be to use lower frequency acoustic excitations, which are less susceptible to aberrations and near field energy loss.
The simulation and phantom experiment data (Figures 5 and and6)6) demonstrate the ability for the Lateral TTP algorithm to accurately reconstruct the shear modulus of a material when the assumptions of material homogeneity and negligible dispersion are satisfied. These datasets are free of noise and physiologic motion, and they allow the theoretical performance of the algorithm to be evaluated over the range of shear moduli that are expected in healthy and diseased livers. While the agreement between the reconstructed and actual moduli is good, the variance increases with stiffness. This effect is due to the fixed tracking beam pulse repetition frequency (PRF) that dictates the temporal sampling of the shear wave. As demonstrated in Figure 3, the more compliant medium (μ = 1.33 kPa) has a slower shear wave speed and greater sampling per shear wavelength than the stiffer medium (μ = 8 kPa). This results in increased variance in the determination of the time of peak displacement in the stiffer medium. Clearly, increasing the PRF of the tracking beams would benefit estimates in stiffer media; however, as with conventional B-mode imaging, the maximum PRF is dictated by the focal depth and the propagation speed of the imaging pulse.
One of the benefits of the Helmholtz algorithm, as implemented by several groups (Oliphant et al., 2001; Bercoff et al., 2004; Nightingale et al., 2003), is that a priori knowledge of the shear wave’s propagation direction is not necessary. The Lateral TTP algorithm, as implemented in this manuscript, has assumed that shear waves are propagating parallel to the lateral dimension (perpendicular to the central axis of the excitation), thus limiting the axial extent of its application to the DOF of the excitation beam, as demonstrated in Figure 1. Incorporation of geometric compensations for shear wave propagation direction in the near field are under investigation, but decreased displacement magnitudes away from the focal zone will make reconstructions at these depths challenging. As demonstrated by the comparison of the VF10-5 and PH4-1 arrays in the phantom experiments, the configuration of the array used to generate the acoustic radiation force and track the resulting displacements should not impact the Lateral TTP algorithm, as long as analysis is confined to the DOF.
In the current study, highly focused excitation beams (F/1.5 and F/2) were used to generate larger displacements because of the greater number of transmit elements. The depths over which this lateral propagation assumption is valid could be expanded by using larger F/# excitations to extend the DOF at the expense of using fewer excitation elements and generating less displacement. Another approach would be to implement sequences that interrogate multiple axial locations, either sequentially or in a rapid fire mode, as is performed in supersonic imaging, to simulate a ‘line source’ of shear waves, although the latter approach can be limited by system power constraints (Bercoff et al., 2004).
In addition to assuming that the direction of shear wave propagation is known, the Lateral TTP algorithm also assumes that significant shear wave dispersion does not occur over the ROI that would cause distortion of the shear wave’s shape and would reduce the correlation between the TTP displacement and the mean energy of the shear wave. Such dispersion would cause the TTP displacements to trend nonlinearly with respect to lateral position, making the linear regression assumptions invalid. Nonlinear TTP displacement relationships do not appear to pose a significant challenge in liver modulus quantification, but if present, this challenge could be reduced by restricting the regression domains, as was done in Figure 6, or using a nonlinear fit to quantify the dispersion. Such trends would appear in the current implementation of the algorithm as linear fits with poor goodness of fit metrics. As demonstrated in Figure 7(c–d), the goodness of fit metrics allow for TTP data that are corrupted by jitter, noise, and/or physiologic motion to be removed from the analysis without any user intervention.
The use of the linear motion filter applied over a dynamic temporal range that varies for each lateral position allows the underlying linear motion assumption to be restricted to a subdomain of the total temporal data acquired, making the assumption more valid. The use of higher-order filters (quadratic) has not yielded different results in liver data to date.
Overall, the human data presented herein demonstrate the feasibility of using acoustic radiation force methods to quantify liver moduli clinically, using an intercostal imaging approach that coincides with the location that liver biopsies are currently performed. While these results are promising, larger studies are necessary to optimize the imaging parameters and to confirm the clinical utility of these methods for diagnosing liver pathology.
Measurement of times to peak displacements at laterally offset positions within the depth of field of an acoustic radiation force excitation allows for accurate estimation of shear wave speeds and reconstruction of shear moduli in homogeneous elastic media. This approach, termed the Lateral TTP algorithm, has been successfully validated in simulation and phantoms, and has been demonstrated in vivo in 20 human livers. The Lateral TTP algorithm does not require second-order temporal and spatial differentiation of displacement data, as is done in Helmholtz reconstructions, but does rely on a presumed direction of shear wave propagation and minimal shear wave dispersion over the region of interest. Simulation studies indicate that shear moduli of healthy livers (μ = 0.8–3.0 kPa) can be reconstructed with a precision of 0.3 kPa, while more fibrotic tissues (μ = 10–15 kPa) can be reconstructed with a precision of 1.0 kPa. Liver modulus reconstructions have been successfully performed in vivo in human volunteers, with shear moduli (0.8–3.0 kPa) consistent with those reported in the literature. The ability to generate at least 10 μm displacements in the liver leads to increased valid shear wave reconstructions in vivo, and the application of linear regression goodness of fit metrics allows for data corrupted by noise and/or physiologic motion to be excluded from the analysis. Obesity (BMI > 30) was not an obstacle to performing shear modulus reconstructions in vivo in the volunteer study. These results demonstrate the feasibility of using radiation force methods to non-invasively quantify liver stiffness, which may be used clinically to correlate with liver fibrosis and to longitudinally monitor disease progression and aid in treatment decisions.
This work was supported by NIH grants R01 EB002132, and R01 CA114075. The authors would also like to thank Dr. Gregg Trahey for his valuable insights.
Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.