Refraction-corrected backscatter tensor imaging of excised porcine ventricular myocardium Open Access

Backscatter tensor imaging (BTI) is an ultrasound imaging technique that estimates the orientation of microstructure in fibrous tissues by detecting variation with direction of the spatial coherence of received echoes.¹ The spatial coherence of echoes from the tissue is proportional to the alignment of the linear array with the fiber direction at the ultrasound focal spot.² Coherent plane wave compounding³ (PWC) is employed in BTI to synthesize a focus throughout the imaging region, thereby reducing the number of transmits required compared to using focused waves to obtain the images.

Benchtop and in vivo experimental configurations often comprise two or more horizontal layers. In these configurations, the ultrasound pulses in transmit and receive are subjected to distinct propagation speeds as well as refraction through the interfaces between layers. Refraction-correction for ultrasound imaging of layered regions has been applied to soft tissue^4,5 and bone.⁶ 

In this letter, we report BTI of excised right- and left-ventricular porcine myocardium using measured sound speeds to account for refraction of the transmitted and reflected pulses through the fluid-tissue interface. Refraction-correction significantly increases both the average spatial coherence and the anisotropy of the spatial coherence, i.e., contrast between the along- and across-fiber directions, which indicate superior beamforming and a larger effect of the microstructure on the backscattered echoes, respectively. Thus, refraction-correction significantly improves the sensitivity of BTI compared to the assumption of a homogeneous imaging region.

An image is first obtained using a single un-steered plane wave of two reference reflectors with a measured separation distance of 6.34 mm, submerged in phosphate buffered saline (PBS), using a commercial linear array transducer (L22-14v, Verasonics, Kirkland, WA) controlled by a programmable ultrasound platform (Vantage, Verasonics). The transmitted waveform consists of a four-cycle pulse with center frequency of 18 MHz, sampled at 62.5 MHz. The reference reflector near to the imaging array is the upper surface of a piece of flat aluminum that rests on the bottom of the tank that contains the fluid. The far reference reflector is the bottom of the tank, which is composed of POM (Delrin). A representative B-mode image of the reference imaging configuration is shown in Fig. 1(b). The time-of-flight of echoes received from each reference reflector are used to calculate the speed of sound in the PBS as described in Sec. 2.3.

BTI scans are performed on one sample each of excised right- and left-ventricular porcine myocardium. Samples were cut from the corresponding ventricular free wall with approximate dimensions of 20 × 20 mm and thinned to thickness h < 5 mm. With care taken to ensure that the position of the aluminum reference reflector is maintained, the myocardium is positioned between the imaging array and the aluminum reference reflector [Figs. 1(a) and 1(c)]. The sample is suspended using a small amount of tension with a biaxial testing machine (BioTester, Cellscale, Waterloo, Canada) to ensure that the surface of the sample is approximately parallel with the linear array as is necessary for BTI and application of Eqs. (1) and (3) below. The BTI acquisitions are modeled after the protocol of Papadacci et al.:¹ a rotational scan of the linear array about the imaging (z) axis is performed by a motorized rotation stage (Velmex, Rochester, NY) that is controlled by the programmable ultrasound platform. Images are obtained in steps of 5° over 360°. Each image in the rotational scan employs coherent plane wave compounding (PWC)³ with 41 plane wave transmissions and steering angles ranging from −20° to 20° in steps of 1°.¹ PWC is employed instead of traditional focused transmissions to reduce the number of transmits needed for transmural BTI.¹ In addition to radio frequency (RF) signals corresponding to echoes from the myocardium, in each image, echoes from the near reference reflector are obtained. Arrival times of echoes from the tissue and the near reference reflector in this configuration are used to compute the sound speed in the myocardium as described in Sec. 2.3. Beamformed RF data from the myocardium is used to estimate the transmural fiber orientation as described in Sec. 2.4.

The linear imaging array of ultrasound transducers is located along the x axis and excites waves traveling in the positive z direction (Fig. 2). The interface between the PBS and myocardium is represented by a horizontal layer at z = z₁. The location of a scatterer within the myocardium is given by (x, z). It is assumed that the transmitted pulse propagates through the imaging region as a planar wave and that the echoes return to the array as spherically spreading wavefronts.

A plane wave with steering angle α relative to the x axis is transmitted by the linear array. The time-of-flight of the transmitted wave to the scatterer at $(x,z>z1)$ is

The first term in Eq. (1) corresponds to travel time from the array to the interface and is identical to Eq. (4) in Ref. 3. The second term in Eq. (1) accounts for change in propagation direction (refraction) and sound speed in the second layer. Note that the expression under the square root in Eq. (1) can also be written as $cos2(α′)/c22$⁠, where the plane wave propagates at an angle $α′$ relative to the x axis after refraction through the interface [Fig. 2(b)].

The receive delays are found from the time-of-flight of an acoustic ray emanating from the scatterer location at (x, z) and intersecting with the ultrasound linear array at $(x0,z=0)$⁠. The lateral intersection point x₁ of such a ray with the interface at z₁ is found from numerical solution of a quartic equation,⁷ 

which follows from the law of refraction. The time-of-flight from the scatterer to the linear array is then

where the two terms are the acoustic path length divided by the sound speed in each layer.

Thus, the appropriate delay to apply in a delay-and-sum coherent PWC scheme in a two-layered region of interest is

Received RF data are delayed according to Eq. (4) and then compounded by summation over the set of plane wave steering angles.

The sound speeds in the PBS and the myocardium sample are calculated using a time-of-flight based approach.^8,9 RF data from the reference configuration images [Fig. 1(b)] are used to calculate the sound speed c₁ in the PBS as

where L is the distance between reference reflectors, measured beforehand with calipers to be L = 6.34 mm, and t₀ and t_A are the arrival times of echoes received at the linear array from the far and near reference reflectors, respectively. Echo arrival times t₀ and t_A from the reference reflectors are determined by finding the points in local search windows where the echo envelope amplitude reaches 50% of the pulse maximum. Arrival times t_A are calculated using each of the signals received on the 31 elements with –5.5 mm $<x<−2.5$ mm, and times t₀ are calculated using each of the signals received on the 31 elements with 2.5 mm $<x<5.5$ mm. Equation (5) is evaluated for each pair of t_A and t₀ (497 pairs), and the mean value is used in evaluation of Eq. (6) and in the BTI calculations.

After the reference measurements, the myocardium is introduced as shown in Fig. 1(c), and the sound speed c₂ in the myocardium is calculated following the procedure outlined in Ref. 9. Echo times are obtained from A-lines (envelope of beamformed RF signals, not logarithmically compressed), which are calculated as functions of time for x = 0 in each of the 72 images of the rotational scan by assuming a constant sound speed equal to 1540 m/s. Echo times calculated in this way were found to be robust to this choice of sound speed in the range of 1500–1560 m/s. The arrival time t_D from the reference reflector with the myocardium present was determined by finding the point where the echo envelope amplitude reaches 50% of the maximum. This timing criterion was found to yield the best accuracy and consistency in calculation of $tA−tD$ as evaluated by visual inspection of the delayed reference pulses. Echo arrival times t_B and t_C corresponding to the upper and lower myocardium surfaces, respectively, were defined as the point at which the A-line amplitude reached 25% of the local maximum in a search window approximately corresponding to the surface. This definition of the surface echo timing was chosen based on the resulting consistency between computed values of $tC−tB$⁠, i.e., tissue thickness, and the visual identification of the tissue surface depths in B-mode images. The mean value of $tA−tD$ and $tC−tB$ over all 72 acquisitions of the rotational scan is used to compute the speed of sound in the myocardium using^8,9

The formulation of Ref. 9 also returns the thickness of the tissue sample, $h=(c2/2)(tC−tB)$⁠, which was used as validation for the algorithm to determine t_B and t_C by comparison with the observed thickness in B-mode images.

The interface $z1=c1tB/2$ between PBS and myocardium is found from the time-of-flight calculation. Ultrasound RF data from the myocardium are delayed and compounded according to Eq. (4) (two-layer case), as well as the baseline case in which it is assumed that both PBS and myocardium have equal sound speeds of 1540 m/s. Spatial coherence of the backscattered echoes from a location in the myocardium is quantified by the coherence factor C of the beamformed data, which is computed as^1,10

where $s(xk,t)$ is the time-delayed and compounded RF signal corresponding to the kth element of the receive aperture (before summation over the receive aperture takes place), N is the number of elements in the receive aperture, and the range $T1<t<T2$ defines a temporal window centered at the focal time (⁠$T2−T1≈0.19$μs, or approximately 3.4 cycles at the RF center frequency, in this work). The coherence factor C is computed on the imaging axis x = 0 through the depth of the tissue in a moving averaging window centered at x = 0 and with dimensions of 1.5 × 0.25 mm in the x and z directions to compensate for heterogeneities that appear as coherent specular scatterers but are unrelated to the fiber direction in the tissue.¹ 

Variations of interest in C with array orientation θ at a given depth of tissue thickness are due to the fiber direction at that depth, with maxima and minima in C expected to occur when the linear array is oriented parallel and perpendicular to the fiber direction, respectively.² Thus, the variations of interest in C should have sinusoidal dependence with period of 180°, and the fiber angle at each depth is determined by a cosine fit to C at each z,

where $θfib(z)$ is the fiber angle. Angles θ and $θfib$ have units of radians in Eq. (8).

Accuracy of BTI is assessed by comparison with the ground-truth fiber angle determined from histological sectioning (Sec. 2.5). Sensitivity of BTI is evaluated for both the baseline and two-layer beamforming strategies by assessing the average coherence factor $A0(z)$⁠, the coefficient of determination R² of the cosine fit, and the fractional anisotropy

A larger value of average coherence A₀ indicates superior beamforming of echoes from the myocardium.¹¹ Higher values of R² and fractional anisotropy (FA) indicate a stronger effect of the microstructure on the backscattered coherence and, therefore, increased sensitivity of BTI for identification of the fiber angle through the tissue thickness.

Transmural fiber angle variation derived from BTI is validated against the ground truth of optical microscopy after histological sectioning. Samples were fixed in 10% neutral buffered formalin and then sectioned and stained with hematoxylin and eosin (H&E). Digital images of each section are taken with a trinocular microscope (1.08 μm/pixel, field of view 2.76 mm × 2.08 mm), and the fiber angle is determined with the OrientationJ/DistributionJ toolbox¹² in ImageJ¹³ with local window set to 40 pixels (43.2 μm), minimum coherency to 2%, and minimum energy to 2%. Similar to previous studies,^14,15 the fiber angle and dispersion are reported as the circular mean and circular standard deviation of the computed pixel-level local orientations.

Sound speeds in the PBS and myocardium during the BTI scan of the RV sample were $c1=1481±4$ m/s and $c2=1540±4$ m/s, respectively. Sound speeds in the PBS and myocardium during the BTI scan of the left-ventricular (LV) sample were $c1=1479±3$ m/s and $c2=1516±5$ m/s, respectively. Sample thicknesses were 4.03±0.04 and 2.67±0.04 mm for the RV and LV samples, respectively.

Results of the BTI calculation in the baseline and two-layer (refraction-corrected) cases are presented in Table 1 for both RV and LV samples. Detailed results are presented in Fig. 3 for the RV tissue scan. While not presented or discussed in detail here, the corresponding results for the sample of LV myocardium are qualitatively identical. Normalized coherence factor $C/A0$ as a function of array rotation angle θ and depth are presented in Figs. 3(a) and 3(b) for the baseline and two-layer cases, respectively. Depth is expressed in terms of percent thickness through the tissue, where the upper and lower surfaces of the tissue were identified as described in Sec. 2.3, and the lower surface of the tissue is the epicardium. Regions of high and low coherence correspond to the orientations θ of the linear array that are parallel and perpendicular to the fiber direction at that depth, respectively. Increased sharpness of the bands of high and low coherence in the two-layer case is apparent in the refraction-corrected case [Fig. 3(b)] compared to baseline [Fig. 3(a)].

Accuracy of the fiber angle estimation is presented in Figs. 3(c) and 3(d). Ultrasound-derived fiber orientation is in good agreement with histology for 0%–80% thickness for both baseline and two-layer cases. The region of tissue from 80% to 100%, indicated by the gray region in Figs. 3(c) and 3(d), is a region of rapidly varying and non-uniform microstructure unique to porcine RV epicardium.¹⁶ The rapid change in microstructure is seen in Fig. 3(c) by the abrupt jump in fiber angle of approximately 70° at 80% thickness. BTI is not effective in that region of the tissue thickness due to the rapidly changing and disorganized tissue microstructure. The LV sample exhibited smooth transmural variation in fiber orientation, with BTI able to estimate the fiber orientation throughout the thickness.

The main result of the present work is in the significant increase in sensitivity metrics for BTI throughout the tissue thickness when accounting for the slower sound speed in the PBS and for refraction through the fluid-tissue interface. Computed coherence C of the backscattered echoes from 50% thickness is shown versus array angle θ in Figs. 3(e) and 3(f) for the baseline and two-layer cases, respectively. Comparison of the cosine fit (red curves) in Figs. 3(e) and 3(f) demonstrates the increased average coherence A₀ and amplitude of coherence modulation A₁ in the two-layer case compared to baseline. Larger values for the average coherence A₀ indicate increased accuracy of beamforming, and larger values of A₁ indicate increased sensitivity of BTI to the fiber structure when the slower sound speed of the PBS and refraction are accounted for. Larger values of A₀ and A₁ in the two-layer case are found throughout the tissue thickness apart from the epicardial layer.

The coefficient of determination R² for the cosine fit is shown versus tissue thickness in Fig. 3(g). Outside of the epicardial layer, the value of R² is higher in the two-layer case (red curve) compared to baseline (blue curve), indicating that the fibrous microstructure accounts for more of the variation in coherence with probe angle when the slower sound speed in the PBS and refraction are accounted for. FA of the coherence is presented in Fig. 3(h). Higher values of FA are obtained in the two-layer case compared to baseline, indicating increased contrast of spatial coherence and increased sensitivity of BTI when the slower sound speed of the PBS and refraction are accounted for. Higher anisotropy of coherence allows for identification of the fiber orientation in tissues with less organized fiber structure, e.g., tissue with disarrayed fibers or a distribution of fiber orientations at a given location in the tissue thickness.² 

The same metrics were computed for the sample of LV myocardium and are presented along with those for the RV sample in Table 1. Root mean square error (RMSE) values compared to histology are similar to those reported for BTI in porcine LV tissues with a two-dimensional (2D) array.¹⁷ Significant increases in the average transmural A₀, R², and FA are seen in both samples throughout the tissue thickness (values for the RV sample are calculated ignoring the epicardial layer from 80% to 100% thickness). Average A₀ and FA are notably lower in this study compared to those of Papadacci et al.,¹ who report $A0≈0.5$ and $FA≈0.4$ in porcine LV myocardium. The difference is possibly due to the higher frequency used in this study (18 MHz here versus 6 MHz in Ref. 1), which is more susceptible to reduced coherence resulting from off-axis specular scattering from small heterogeneities.

This study demonstrates the benefits of accounting for the distinct sound speeds in a layered imaging region during estimation of transmural fiber orientation using BTI. For example, by using this approach, the sensitivity may be increased for BTI performed in applications for which the sample cannot contact the imaging array or be embedded in a gel, such as biaxial mechanical testing.

This work was supported by National Institutes of Health (NIH) Grant Nos. F32HL160051 and T32HL129964.