System and Method for Non-Invasive Determination of Pressure in a Biological Compartment

Abstract

Systems and methods are provided for a mechanical analysis for determining a kinetic value, such as pressure, from ultrasonically measured Lamb wave speed in a pressurized biological compartment. Acoustoclastic analysis of Lamb waves may be used to determine pressure in a biological compartment, such as with a spherical pressure-vessel analysis. Non-invasive pressure estimation is provided in biological compartments which is relevant to a range of biological structures.

Description

STATEMENT REGARDING FEDERALLY SPONSORED RESEARCH

N/A

BACKGROUND

Biological compartments that sustain internal fluid pressure are ubiquitous structures in mammalian organisms. Various scaling relationships between anatomical structures of this character in different species can be understood in relation to their mechanical role as pressurized compartments. In relation to human health and medicine, the pressures that such compartments are exposed to are often vital indicators of their functioning. Examples, where pressure measurements hold clinical significance in biological compartments include the diaphragm and lungs, large blood vessels, the heart, the bladder, and body cavities.

Pressure measurements are typically collected by invasive procedures that use catheters instrumented with pressure transducers. Because of the potential for discomfort and morbidity is high in such procedures, there is great interest toward non-invasive assessments.

Due to the sensitivity of shear wave speed to tissue loading [11], wave-based elastographic techniques (i.e. shear wave and Lamb wave elastography) are increasingly being investigated as a noninvasive means of inferring loads sustained by soft tissues. Recently, a variety of investigations have correlated internal fluid pressures with transverse wave speeds at the walls of a variety of biological compartments. However, these investigations are empirical and the underlying mechanics that produce such associations has received limited treatment. An analysis that relates elastic wave speed to pressure in such structures would need to incorporate analyses for acoustoelastic wave propogation in the presence of finite deformation in a nonlinear and hyperelastic solid; wave dispersion associated with the wall thickness of the compartment; and pressure-vessel analysis for highly deformed compartments. None of these criteria have been applied toward estimation of pressure in a biological compartment. Thus there remains a need for a non-invasive method for measuring pressure, such as in a bladder.

SUMMARY OF THE DISCLOSURE

The present disclosure addresses the aforementioned drawbacks by providing systems and methods for a mechanical analysis for determining pressure from ultrasonically measured Lamb wave speed in a pressurized biological compartment. In some configurations, the method includes acoustoelastic analysis of Lamb waves with spherical pressure-vessel analysis. The systems and methods may provide for non-invasive pressure estimation in biological compartments which is relevant to a range of biological structures.

In one configuration, a method is provided for determining a kinetic value of a tissue volume. The method includes, with a transducer, generating Lamb waves in a tissue wall of the tissue volume using radiation force. The tissue volume is formed by a tissue wall that spatially separates a fluid material from a rigid material. The method also includes forming ultrasonic echo data by detecting ultrasonic energy reflected by multiple locations along the tissue volume that is subject to the radiation force. The method also includes determining frequency-resolved wave speed data for the generated Lamb waves from the ultrasonic echo data and determining pressure of the tissue volume based upon the determined frequency-resolved wave speed data.

In one configuration, a system is provided for determining a kinetic value of a tissue volume. The system includes a transducer configured to generate Lamb waves in a tissue wall of the tissue volume using radiation force and detect ultrasonic energy reflected by multiple locations along the tissue volume that is subject to the radiation force to form ultrasonic echo data. The tissue volume is formed by the tissue wall that spatially separates a fluid material from a rigid material. The system also includes a computer system configured to: determine frequency-resolved wave speed data from the generated Lamb waves from the ultrasonic echo data; and determine pressure of the tissue volume based upon the determined frequency-resolved wave speed data.

The foregoing and other aspects and advantages of the present disclosure will appear from the following description. In the description, reference is made to the accompanying drawings that form a part hereof, and in which there is shown by way of illustration a preferred embodiment. This embodiment does not necessarily represent the full scope of the invention, however, and reference is therefore made to the claims and herein for interpreting the scope of the invention. Like reference numerals will be used to refer to like parts from Figure to Figure in the following description.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is a block diagram of an ultrasound imaging system used in accordance with the present disclosure.

FIG. 2 is a block diagram of a transmitter that forms part of the ultrasound imaging system of FIG. 1.

FIG. 3 is a block diagram of a receiver that forms part of the ultrasound imaging system of FIG. 1.

FIG. 4A is a flow chart of non-limiting example steps for of a method to determine kinetic values from ultrasonically measured Lamb wave speed in a pressurized biological compartment.

FIG. 4B is a flow chart of non-limiting example steps for of a method to determine pressure from ultrasonically measured Lamb wave speed in a pressurized biological compartment.

FIG. 5A is a diagram of the basic principal of ultrasound bladder vibrometry (UBV) of the present disclosure.

FIG. 5B illustrates in vivo human B-mode of the bladder with a sketch of the UBV excitation beam.

FIG. 6A depicts a non-limiting example compartment modeled as a hyperelastic spherical pressure vessel with Lamb waves propagating along the circumference.

FIG. 6B depicts the acoustic pulse producing Lamb waves propagating forward and backward away from the pulse of FIG. 6A.

FIG. 6C depicts a 2D-FFT based method used to estimate frequency dependent phase velocities from the waves of FIG. 6B.

FIG. 6D depicts a non-limiting example curve fitting applied to the median wave speeds from the acquisitions collected at each volume of FIG. 6A.

FIG. 7 shows graphs of non-limiting examples of vessel pressures plotted along with the vessel pressure measured by a pressure gauge as a function of bladder volume.

FIG. 8 is a correlation plot of measurements from a non-limiting example set of bladders.

FIG. 9 is a graph of a non-limiting example signal and noise region of interest (ROI).

DETAILED DESCRIPTION

Systems and methods are provided for a mechanical analysis for determining kinetic values, such as pressure, stress, and the like, from ultrasonically measured Lamb wave speed in a pressurized biological compartment. In some configurations, the method includes acoustoelastic analysis of Lamb waves with spherical pressure-vessel analysis. In a non-limiting example, bladder pressure may be determined with varied filling volumes and elasticities. The systems and methods may provide for non-invasive pressure estimation in biological compartments which is relevant to a range of biological structures.

“Kinetic value” is a quantity relating to the state of the forces acting on the material such as pressure, stress and/or force.

A mechanical analysis may be used that directly relates the pressure in a spherical pressure-vessel to measurements of the dispersion curve of a Lamb wave. Pressure may be determined from UBV measurements collected from the biological compartment, such as the bladder.

Conventional methods rely upon simple empirical association or curve fitting for the purpose of pressure estimation in biological compartments. The systems and methods in accordance with the present disclosure may be generalized to compartments of varied characteristics (e.g. elasticity, size, and the like) as long as the characteristics are included within the mechanistic analysis. The systems and methods in accordance with the present disclosure also include greater accuracy than conventional methods. Deviations in individual bladders suggest that additional effects may be affecting the measurement. For a mechanistic analysis, greater accuracy may be attained by examining the simplifying assumptions in the analysis and adapting the analysis and measurements to incorporate additional characteristics. Non-limiting examples of characteristics that may be modified include adjusting a spherical compartment to a non-spherical geometry that may alter the relation between wall-stress and internal pressure.

Referring to FIG. 1, a transducer array 23 includes a plurality of separately driven elements 11 that each produce a burst of ultrasonic energy when energized by a pulse produced by a transmitter 13. The ultrasonic energy reflected back to the transducer array 23 from the subject under study is converted to an electrical signal by each transducer element 11 and applied separately to a receiver 9 through a set of switches 15. The transmitter 13, receiver 9 and the switches 15 are operated under the control of a digital controller 19 responsive to the commands input by the human operator. A complete scan is performed by acquiring a series of echoes in which the switches 15 are set to their transmit position, the transmitter 13 is gated on momentarily to energize each transducer element 11, the switches 15 are then set to their receive position, and the subsequent echo signals produced by each transducer element 11 are applied to the receiver 9. The separate echo signals from each transducer element 11 are combined in the receiver 9 to produce a single echo signal that is employed to produce a line in an image on a display system 17.

The transmitter 13 drives the transducer array 23 such that the ultrasonic energy produced is directed, or steered, in a beam or pulse. A B-scan can therefore be performed by moving this beam through a set of angles from point-to-point rather than physically moving the transducer array 23. To accomplish this, the transmitter 13 imparts a time delay (Ti) to the respective pulses 20 that are applied to successive transducer elements 11. If the time delay is zero (Ti=0), all the transducer elements 11 are energized simultaneously and the resulting ultrasonic beam is directed along an axis 24 normal to the transducer face and originating from the center of the transducer array 23. As the time delay (Ti) is increased, the ultrasonic beam is directed downward from the central axis 24 by an angle θ.

A sector scan is performed by progressively changing the time delays Ti in successive excitations. The angle θ is thus changed in increments to steer the transmitted beam in a succession of directions. When the direction of the beam is above the central axis 24, the timing of the pulses 7 is reversed.

Referring still to FIG. 1, the echo signals produced by each burst of ultrasonic energy emanate from reflecting objects located at successive positions (R) along the ultrasonic beam. These are sensed separately by each segment 11 of the transducer array 23 and a sample of the magnitude of the echo signal at a particular point in time represents the amount of reflection occurring at a specific range (R). Due to the differences in the propagation paths between a focal point P and each transducer element 11, however, these echo signals will not occur simultaneously and their amplitudes will not be equal. The function of the receiver 9 is to amplify and demodulate these separate echo signals, impart the proper time delay to each and sum them together to provide a single echo signal that accurately indicates the total ultrasonic energy reflected from each focal point P located at range R along the ultrasonic beam oriented at the angle θ.

To simultaneously sum the electrical signals produced by the echoes from each transducer element 11, time delays are introduced into each separate transducer element channel of the receiver 9. In the case of the linear array 23, the delay introduced in each channel may be divided into two components, one component is referred to as a beam steering time delay, and the other component is referred to as a beam focusing time delay. The beam steering and beam focusing time delays for reception are precisely the same delays (Ti) as the transmission delays described above. However, the focusing time delay component introduced into each receiver channel is continuously changing during reception of the echo to provide dynamic focusing of the received beam at the range R from which the echo signal emanates.

Under the direction of the digital controller 19, the receiver 9 provides delays during the scan such that the steering of the receiver 9 tracks with the direction of the beam steered by the transmitter 13 and it samples the echo signals at a succession of ranges and provides the proper delays to dynamically focus at points P along the beam. Thus, each emission of an ultrasonic pulse results in the acquisition of a series of data points that represent the amount of reflected sound from a corresponding series of points P located along the ultrasonic beam.

By selecting proper time delays, echoes from multiple focused locations can be received to measure vibration information from several points of the tissue. The limitation of the lateral resolution of the transducer for two closely located points can be improved by assigning different transmitting codes for different locations.

The display system 17 receives the series of data points produced by the receiver 9 and converts the data to a form producing the desired image. For example, if an A-scan is desired, the magnitude of the series of data points is merely graphed as a function of time. If a B-scan is desired, each data point in the series is used to control the brightness of a pixel in the image, and a scan comprised of a series of measurements at successive steering angles (0) is performed to provide the data necessary for display of an image.

Referring to FIG. 2, the transmitter 13 includes a set of channel pulse code memories that are indicated collectively at 900. Each pulse code memory 900 stores a bit pattern 902 that determines the frequency of the ultrasonic pulse 904 that is to be produced. This bit pattern is read out of each pulse code memory 900 by a master clock and applied to a driver 906 that amplifies the signal to a power level suitable for driving the transducer 23. In the example shown in FIG. 2, the bit pattern is a sequence of four “1” bits alternated with four “0” bits to produce a 5 megahertz (“MHz”) ultrasonic pulse 904. The transducer elements 23 to which these ultrasonic pulses 904 are applied respond by producing ultrasonic energy.

As indicated above, to steer the transmitted beam of the ultrasonic energy in the desired manner, the pulses 904 for each of the N channels must be produced and delayed by the proper amount. These delays are provided by a transmit control 908 that receives control signals from the digital controller 19. When the control signal is received, the transmit control 908 gates a clock signal through to the first transmit channel 900. At each successive delay time interval thereafter, the clock signal is gated through to the next channel pulse code memory 900 until all the channels to be energized are producing their ultrasonic pulses 904. Each transmit channel 900 is reset after its entire bit pattern 902 has been transmitted and the transmitter 13 then waits for the next control signal from the digital controller 19.

Referring to FIG. 3, the receiver 9 is comprised of three sections: a time-gain control section 100, a beam forming section 101, and a mid processor 102. The time-gain control section 100 includes an amplifier 105 for each of the N=128 receiver channels and a time-gain control circuit 106. The input of each amplifier 105 is connected to a respective one of the transducer elements 11 to receive and amplify the echo signal that it receives. The amount of amplification provided by the amplifiers 105 is controlled through a control line 107 that is driven by the time-gain control circuit 106. As is well known in the art, as the range of the echo signal increases, its amplitude is diminished. As a result, unless the echo signal emanating from more distant reflectors is amplified more than the echo signal from nearby reflectors, the brightness of the image diminishes rapidly as a function of range (R). This amplification is controlled by the operator who manually sets TGC linear potentiometers 108 to values that provide a relatively uniform brightness over the entire range of the sector scan. The time interval over which the echo signal is acquired determines the range from which it emanates, and this time interval is divided into by the TGC control circuit 106. The settings of the potentiometers are employed to set the gain of the amplifiers 105 during each of the respective time intervals so that the echo signal is amplified in ever increasing amounts over the acquisition time interval.

The beam forming section 101 of the receiver 9 includes N=128 separate receiver channels 110. Each receiver channel 110 receives the analog echo signal from one of the TGC amplifiers 105 at an input 111, and it produces a stream of digitized output values on an I bus 112 and a Q bus 113. Each of these I and Q values represents a sample of the echo signal envelope at a specific range (R). These samples have been delayed in the manner described above such that when they are summed at summing points 114 and 115 with the I and Q samples from each of the other receiver channels 110, they indicate the magnitude and phase of the echo signal reflected from a point P located at range R on the steered beam (0).

Referring still to FIG. 3, the mid processor section 102 receives the beam samples from the summing points 114 and 115. The I and Q values of each beam sample is a 16-bit digital number that represents the in-phase and quadrature components of the magnitude of the reflected sound from a point (R,0). The mid processor 102 can perform a variety of calculations on these beam samples, where choice is determined by the type of image to be reconstructed.

For example, a conventional ultrasound image may be produced by a detection processor 120 that calculates the magnitude of the echo signal from its I and Q components:

$M=\sqrt{I^2+Q^2}$

The resulting magnitude values output at 121 to the display system 17 result in an image in which the magnitude of the reflected echo at each image pixel is indicated.

The present disclosure may be implemented by a mechanical property processor 122 that forms part of the mid-processor 102. As will be explained in detail below, this processor 102 receives the I and Q beam samples acquired during a sequence of measurements of the subject tissue and calculates a mechanical property (i.e., thickness) of the tissue.

Referring to FIG. 5A and FIG. 5B, a focused ultrasound beam 605 may be used to induce an acoustic radiation force that excites waves 630 in the wall of a compartment 608 and a pulse-echo detection beam 606 may be used to measure the motion at several points along a line of propagation in order to track tissue motion through successive B-modes. The line of propagation being, for example, a bladder wall. Though the following description is made with respect to the bladder wall, other elastic tissue volumes may likewise be considered. The B-mode scanning can be performed by, for example, a linear array transducer 602 with ultrasound system 600 and a data processing module for kinetic value estimation 601, as shown in FIG. 5A. The detection beam 606 can be transmitted with a pulse repetition frequency (PRF). FIG. 5B illustrates a B-mode cross-section of a non-limiting example compartment.

Referring to FIG. 4A, a flow chart is provided setting forth example steps 500 of a method to determine a kinetic value of a biological compartment, such as pressure, tensile stress, and the like. The process begins at process block 502 with the formation of acoustic radiation force excitation signals that are applied, at process block 504, to a subject to, thereby, apply vibration pulses to a tissue wall of a tissue of interest.

At process block 506, a pulse-echo ultrasound detection beam is used to measure the motion at several points along a line of propagation, the line of propagation being, for example, a bladder wall, in order to track tissue motion through successive B-modes. Cross-spectral analysis of the received echoes can then be used to calculate the wall motion as a function of time. Because an impulse contains mostly the frequency components up to the inverse of the time length of the impulse, many phase velocities from a single impulse push can be extracted.

As shown at process block 508, a two-dimensional fast Fourier transform (2D-FFT) of the wall motion may be used to calculate the change of Lamb wave velocity as a function of frequency, or the Lamb wave dispersion. In other words, Fourier-space analysis of the motion may be used as a function of time that yields the k-space whose coordinates are frequency, f, and wave number, k.

Since the wave velocity c=k/f, the phase velocity at each frequency can be calculated at process block 510 by searching for peaks at the given frequency and dividing by the wave number coordinate by the frequency coordinate for the given peak. The dispersion data may be fit to the pre-stressed anti-symmetric Lamb wave model at process block 512, to determine a kinetic value (such as pressure, stress, and the like), as shown at process block 518. As indicated by process block 514, the digital controller determines whether the last frequency has been measured when a harmonic excitation approach is used. If not, at process block 516, another frequency is selected and process blocks 506 through 512 are repeated at each desired prescribed frequency.

Referring still to FIG. 4A, an ultrasound deformation estimation may be performed where ultrasound imaging frames are acquired at process block 520. The resulting image may be segmented at process block 522. The geometry of the compartment and the wall thickness may be determined at process block 524. Large-strain deformation of the compartment may be determined at process block 526, which may be fed into process block 512.

Referring now to FIG. 4B, a flow chart is provided setting forth non-limiting example steps 550 of a method to determine pressure from ultrasonically measured Lamb wave speed in a pressurized biological compartment. An ultrasound system may be used to form an acoustic radiation force excitation signal at step 552. Using the radiation force excitation signal, vibration pulses may be applied to a compartment with tissue of interest of a subject at step 554. Lamb waves may be generated therefrom in the wall of the compartment at step 556. Any necessary corrections for the geometry of the compartment may be determined at step 558, such as adjusting for the curvature of a bladder wall, and the like.

Axial particle velocities in the compartment wall may be determined using phase-based autocorrection at step 560. Frequency-resolved wave speed of the Lamb waves may be determined using a 2D Fourier transformation of the particle velocities at step 562. Pressure of the compartment may be determined using an optimization curve fitting of an expression that includes wave speed and compartment characteristics at step 564. An expression may include a dispersion relation that describes wave speed as a function of pressure. In a non-limiting example, the expression may be of the form of eq. (26) below. Characteristics may include compartment geometry, solid and fluid wave densities, and the like. A signal to noise ratio (SNR) may optionally be determined for the acquired data at step 566 and any data below a threshold SNR value may be excluded from the analysis at step 568.

Expressions of stress and deformation in terms of principle stress and stretches and associated stretch invariants may be used in a mechanical analysis. The spectral decomposition of the deformation tensor F can be expressed as

$\begin{array}{cc}F=RU=VR=\sum _{i=1}^3{\lambda }_i{v}^{\left(i\right)}\otimes u^{\left(i\right)}=\sum _{i=1}^3{\lambda }_i{\mathrm{Ru}}^{\left(i\right)}\otimes u^{\left(i\right)}& \left(1\right)\end{array}$

Here, R represents the orthogonal tensor associated with rigid body rotation, U is the right stretch tensor and, λ_i denotes the ith principle stretch and u⁽ⁱ⁾ and v⁽ⁱ⁾ denotes associated unit vector of the principle stretch in the respective material and spatial (e.g. undeformed and deformed) configurations of the medium. In some configurations, v⁽ⁱ⁾=Ru⁽ⁱ⁾. ⊗ denotes a diadic product of two vectors. To remove the influence of R from the deformation measure, constitutive equations may be formulated in terms of the right and left Cauchy-Green deformation tensors, i.e. C=FTF=U² and B=FF^T=V², the latter of which may include the following spectral decomposition and invariants:

$\begin{array}{cc}B=\sum _{i=1}^3{\lambda }_i^2{v}^{\left(i\right)}\otimes v^{\left(i\right)}& \left(2\right)\end{array}$ $\begin{array}{cc}{I}_1=tr\left(B\right)=tr\left(B\right)={\lambda }_1^2+{\lambda }_2^2+{\lambda }_3^2& \left(3\right)\end{array}$ $\begin{array}{cc}{I}_2=1/2\left(t{r\left(B\right)}^2-tr\left(B^2\right)\right)={\lambda }_1^2{\lambda }_2^2+{\lambda }_2^2{\lambda }_3^2+{\lambda }_1^2{\lambda }_3^2& \left(4\right)\end{array}$ $\begin{array}{cc}{I}_3=\det \left(B\right)={\lambda }_1^2{\lambda }_2^2{\lambda }_3^2& \left(5\right)\end{array}$

In the case of incompressibility or isochoric motion, I₃=1. For an isotropic material, the Cauchy stress (σ) may be co-axial with B and V. The principle Cauchy stresses can be expressed relative to the principle stretches and principle vectors of V and a strain energy function W. Adding in the additional constraint of incompressibility results in the following expression:

$\begin{array}{cc}\sigma =\sum _{i=1}^3{\lambda }_i\frac{\partial W}{\partial {\lambda }_l}{v}^{\left(i\right)}\otimes v^{\left(i\right)}-pI& \left(6\right)\end{array}$

where, p is a Lagrangian multiplier that may be added to enforce the incompressiblity constraint.

Non-Limiting Example Spherical Pressure Vessel Model of the Bladder

Relating wall stresses to vessel pressure, the biological compartment may be modeled as a spherical, thin-walled pressure vessel. The classic relation for Cauchy stress along the circumference of the sphere (set to the first principle Cauchy stress σ₁) is given by:

$\begin{array}{cc}{\sigma }_1=\frac{Pr}{2t}=\frac{Pr}{4h}& \left(7\right)\end{array}$

where P represents the gauge pressure relative to the exterior vessel wall and r, t and h are the radius, wall thickness and half-wall thickness of the deformed sphere. σ₂, the principle Cauchy stress in the radial direction, takes the following values at the boundaries of the vessel walls

$\begin{array}{cc}{\sigma }_2=\{\begin{array}{c}-P\mathrm{at}\mathrm{the}\mathrm{interior}\mathrm{vessel}\mathrm{wall}\\ 0\mathrm{at}\mathrm{the}\mathrm{exterior}\mathrm{vessel}\mathrm{wall}\end{array}& \left(8\right)\end{array}$

If r/4h >>1, which is the case in a thin-walled vessel, the difference in the wall stresses can be approximated about the thickness of the vessel wall by

$\begin{array}{cc}{\sigma }_1-{\sigma }_2\approx P\frac{r}{4h}& \left(9\right)\end{array}$

In the case of incompressiblity, the principle stretches can be defined by the vessel radius r relative to its initial value ro. They can then, in turn, be determined from the initial and deformed volumes V and V₀ via r=(3/4V/π)^1/3. These relations may be of the forms:

$\begin{array}{cc}{\lambda }_1=\frac{r}{r_0}={\left(\frac{V}{V_0}\right)}^{1/3}& \left(10\right)\end{array}$ $\begin{array}{cc}{\lambda }_2=\frac{r_0^2}{r^2}={\left(\frac{V_0}{V}\right)}^{2/3}& \left(11\right)\end{array}$

The initial fill volume V₀ may be difficult to acquire practically. In some configurations to address this challenge, V₀ may be determined from deformed and undeformed tissue wall thickness as:

$V_0={\left(\frac{h}{h_0}\right)}^{3/2}V,$

where h is the deformed tissue wall thickness, and h₀ is the undeformed tissue wall thickness.

Lamb Wave Propogation in a Pre-Stressed Plate

In some configurations, a Lamb wave analysis may be performed. An expression for Lamb waves propagating in a pre-stressed, isotropic, incompressible and hyperelastic plate surrounded by fluid is provided below in eq. (12), which is a dispersion equation for axisymmetric transverse waves propogating in such a condition

$\begin{array}{cc}\gamma {s_1\left(1+s_2^2\right)}^2\tanh \left(s_{1}kh\right)-\gamma {s_2\left(1+s_1^2\right)}^2\tanh \left(s_{2}kh\right)+\frac{{\rho }^F{c}^2}{\xi }\left(s_1^2-s_2^2\right)=0& \left(12\right)\end{array}$

where h is half the thickness of the deformed hyperelastic plate (i.e. h=λ₂h₀), c is the wave speed, pF is the mass density of the fluid (assumed to be 1 g/ml). ζ is a parameter defined by

$\xi =\sqrt{1-\frac{c^2}{c_p^2}}$

where c is the speed of the transverse wave and c_p is the speed of sound in the fluid (assumed to be 1480 m/s in this current analysis). k is the wave number and it defines the frequency dependence of the wave speed through the relation k=2πf/c. s₁² and s₂² are roots of the quadratic equation defined by

$\begin{array}{cc}\gamma s^4-\left(2\beta -\rho c^2\right)s^2+\alpha -\rho c^2=0& \left(13\right)\end{array}$

Here, p represents the mass density of the solid wall, which in some configurations may be assumed to be 1 g/ml. γ, α and β are the acoustoelastic parameters for a pre-stressed hyperelastic solid. They may be defined by the components of the incremental elasticity tensor A_0piqj=F_prF_qs(∂²W/∂F_ir∂F_js). Using relationships derived between A_0piqj, W and λ_i, the following expressions may be generated

$\begin{array}{cc}\alpha =A_{01212}=\frac{{\lambda }_1^2\left({\lambda }_1\frac{\partial W}{\partial {\lambda }_1}-{\lambda }_2\frac{\partial W}{\partial {\lambda }_2}\right)}{{\lambda }_1^2-{\lambda }_2^2}& \left(14\right)\end{array}$ $\begin{array}{cc}\gamma =A_{02121}=\frac{{\lambda }_2^2\left({\lambda }_1\frac{\partial W}{\partial {\lambda }_1}-{\lambda }_2\frac{\partial W}{\partial {\lambda }_2}\right)}{{\lambda }_1^2-{\lambda }_2^2}& \left(15\right)\end{array}$ $\begin{array}{cc}2\beta =A_{01111}+A_{02222}-2A_{01122}-2A_{01221}={\lambda }_1^2\frac{{\partial }^{2}W}{\partial {\lambda }_1^2}+{\lambda }_2^2\frac{{\partial }^{2}W}{\partial {\lambda }_2^2}-2{\lambda }_1{\lambda }_2\frac{{\partial }^{2}W}{\partial {\lambda }_2{\lambda }_1}-2\frac{{\lambda }_1{\lambda }_2\left({\lambda }_2\frac{\partial W}{\partial {\lambda }_1}-{\lambda }_2\frac{\partial W}{\partial {\lambda }_2}\right)}{{\lambda }_1-{\lambda }_2}& \left(15\right)\end{array}$

Here, λ₁ may be aligned with the circumferential direction of the vessel wall while λ₂ is in the radial direction. These expressions are a consequence of the same co-axility between stress and strain measures and so material isotropy may be assumed in the expressions. A substitution of the components of the principle Cauchy stresses (e.g. σ_i=λ_i∂W/∂λ_i) from eq. (6) into eq. (14) and eq. (15) produces

$\begin{array}{cc}\alpha =\frac{{\lambda }_1^2}{{\lambda }_1^2-{\lambda }_2^2}\left({\sigma }_1-{\sigma }_2\right)& \left(16\right)\end{array}$ $\begin{array}{cc}\gamma =\frac{{\lambda }_2^2}{{\lambda }_1^2-{\lambda }_2^2}\left({\sigma }_1-{\sigma }_2\right)& \left(17\right)\end{array}$

Eq. (17) and eq. (18) are model-independent, but substitutions may not be made for β without producing derivatives of the components of σ. In some configurations, a constructive model may be used to determine the wall stress. In a non-limiting example, the Demiray-Fung constitutive model may be used to produce an expression for β of a similar form as that of eq. (17) and eq. (18). A strain energy function can be expressed in terms of the first invariant of B (e.g. eq. (3)) as follows

$\begin{array}{cc}W=\frac{{\mu }_0}{2b}\left(e^{b\left(I_1-3\right)}-1\right)& \left(18\right)\end{array}$

Where μ₀ represents the initial shear modulus of the medium and b is a unitless parameter that determines the rate of strain stiffening within the medium. Substituting eq. (19) into eq. (16) produces the following expression for ß for the non-limiting example model:

$\begin{array}{cc}2\beta =2\left({\lambda }_1^4+{\lambda }_2^4-2{\lambda }_1^2{\lambda }_2^2\right){\mu }_{0}b\left(e^{b\left(I_1-3\right)}\right)+2\left({\lambda }_1^2+{\lambda }_2^2\right){\mu }_0\left(e^{b\left(I_1-3\right)}\right)& \left(19\right)\end{array}$

An expression for stresses can be produced by substituting eq. (19) into eq. (6) and subtracting the components to eliminate p, which produces:

$\begin{array}{cc}{\sigma }_1-{\sigma }_2={\lambda }_1\frac{\partial W}{\partial {\lambda }_1}-{\lambda }_2\frac{\partial W}{\partial {\lambda }_2}=\left({\lambda }_1^2-{\lambda }_2^2\right)\frac{\partial W}{\partial I_1}=\left({\lambda }_1^2-{\lambda }_2^2\right){\mu }_0{e}^{b\left(I_{1}3\right)}& \left(20\right)\end{array}$

An expression for β in terms of stress can be produced by substituting eq. (21) into eq. (20), which results in:

$\begin{array}{cc}2\beta =\frac{\left({\lambda }_1^4+{\lambda }_2^4-2{\lambda }_1^2{\lambda }_2^2\right)b+{\lambda }_1^2+{\lambda }_2^2}{{\lambda }_1^2-{\lambda }_2^2}\left({\sigma }_1-{\sigma }_2\right)& \left(21\right)\end{array}$

The only material parameter remaining is the nonlinearity parameter, b, which may be assumed or estimated. Having γ, α, and β in terms of stress and deformation now permits eq. (12) to be defined in terms of pressure and geometric measures. The principle stretch terms may be encompassed under coefficients, A and B, as follows:

$\begin{array}{cc}A=\frac{{\lambda }_1^2}{{\lambda }_1^2-{\lambda }_2^2}=\frac{{\left(\frac{V}{V_0}\right)}^{2/3}}{{\left(\frac{V}{V_0}\right)}^{2/3}-{\left(\frac{V_0}{V}\right)}^{4/3}}& \left(22\right)\end{array}$ $B=\frac{\left({\lambda }_1^4-{\lambda }_2^4-2{\lambda }_1^2{\lambda }_2^2\right)b+{\lambda }_1^2+{\lambda }_2^2}{{\lambda }_1^2-{\lambda }_2^2}=\frac{({\left(\frac{V}{V_0}\right)}^{8/3}-{\left(\frac{V_0}{V}\right)}^{16/3}-{\left(\frac{V}{V_0}\right)}^{2/3}{\left(\frac{V_0}{V}\right)}^{4/3}}{{\left(\frac{V}{V_0}\right)}^{2/3}-{\left(\frac{V_0}{V}\right)}^{4/3}}+\frac{{\left(\frac{V}{V_0}\right)}^{2/3}+{\left(\frac{V_0}{V}\right)}^{4/3}}{{\left(\frac{V}{V_0}\right)}^{2/3}-{\left(\frac{V_0}{V}\right)}^{4/3}}$

Propagating all substitutions from eq. (17), eq. (18) and eq. (22) into eq. (12) and substituting for 01-02 with eq. (9) produces the non-limiting dispersion equation shown in eq. (25).

$\begin{array}{cc}\frac{V_0^2}{V^2}A\left(P_{4h}^r\right)\sqrt{\frac{\frac{\rho c^2}{P\frac{r}{4h}}-2B}{2\frac{V_0^2}{V^2}}+\sqrt{\frac{\left(2B\left(P\frac{r}{4h}\right)-\rho c^2\right)-4\left(\frac{V_0^2}{V^2}{A}^2\left(P\frac{r}{4h}\right)\rho c^2\right)}{2\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)}}}\times {\left(1+\frac{\frac{\rho c^2}{P\frac{r}{4h}}-2B}{2\frac{V_0^2}{V^2}A}-\frac{\sqrt{{\left(2B\left(P\frac{r}{4h}\right)-\rho c^2\right)}^2-4(\frac{V_0^2}{V^2}{A^2\left(P\frac{r}{4h}\right)}^2-\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)\rho c^2}}{2\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)}\right)}^2\times \tanh \left(\frac{\frac{\rho c^2}{P\frac{r}{4h}}-2B}{2\frac{V_0^2}{V^2}A}+\frac{\sqrt{{\left(2B\left(P\frac{r}{4h}\right)-\rho c^2\right)}^2-4(\frac{V_0^2}{V^2}{A^2\left(P\frac{r}{4h}\right)}^2-\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)\rho c^2})}{2\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)}\mathrm{kh}\right)-\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)\left(\frac{\frac{\rho c^2}{P\frac{r}{4h}}-2B}{2\frac{V_0^2}{V^2}A}-\frac{\sqrt{{\left(2B\left(P\frac{r}{4h}\right)-\rho c^2\right)}^2-4\left(\frac{V_0^2}{V^2}{A^2\left(P\frac{r}{4h}\right)}^2-\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)\rho c^2\right)}}{2\frac{V_0^2}{V^2}A\left(P_{4h}^r\right)}\right)\times {\left(1+\frac{\frac{\rho c^2}{P\frac{r}{4h}}-2B}{2\frac{V_0^2}{V^2}A}+\frac{\sqrt{{\left(2B\left(P\frac{r}{4h}\right)-\rho c^2\right)}^2-4\left(\frac{V_0^2}{V^2}{A^2\left(P\frac{r}{4h}\right)}^2-\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)\rho c^2\right)}}{2\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)}\right)}^2\times \tanh \left(\frac{\frac{\rho c^2}{P\frac{r}{4h}}-2B}{2\frac{V_0^2}{V^2}}-\frac{\sqrt{{\left(2B\left(P\frac{r}{4h}\right)-\rho c^2\right)}^2-4\left(\frac{V_0^2}{V^2}{A^2\left(P\frac{r}{4h}\right)}^2-\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)\rho c^2\right)}}{2\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)}\mathrm{kh}\right)+\frac{{\rho }^F{c}^2}{\xi }\left(2\frac{\sqrt{{\left(2B\left(P\frac{r}{4h}\right)-\rho c^2\right)}^2-4\left(\frac{V_0^2}{V^2}{A^2\left(P\frac{r}{4h}\right)}^2-\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)\rho c^2\right)}}{2\frac{V_0^2}{V^2}A\left(P\frac{r}{4h}\right)}\right)=0& \left(24\right)\end{array}$

Where h=λ₂h₀=h₀(V₀/V)^2/3 and r=(3/4V/π)^1/3. Eq. (26) is fully defined in terms of pressure (P), solid and fluid wave densities and wave speeds (ρ, ρ^F, C_p and c), and compartment geometric measurements (V₀, V and h₀).

Non-Limiting Example Ex Vivo Bladder Measurements

In a non-limiting example, UBV measurements collected from ex vivo pig bladders were used for validation. Nine ex vivo bladders were measured before and after a crosslinking treatment with formalin. One bladder was excluded from the analysis due to a low initial volume (110 ml) that precluded it from inclusion. The remaining bladders had initial filling volumes ranging from 190 to 330 ml. During measurements, a syringe was used to add water to the bladder at 10 ml increments up to a total of 280 ml of added fluid and a pressure gauge was used to measure the vessel pressure of the bladder. A programmable ultrasound system equipped with a linear array ultrasound transducer was used to collect UBV measurements. Lamb waves were excited at the wall of the bladder through acoustic radiation force by a focused ultrasound beam with 400 us tone burst. Ultra-fast ultrasound imaging at 4000 frames per second captured the resulting wave propagation. Five acquisitions were collected at each volume.

After adjusting for the curvature of the bladder wall, the axial particle velocities were averaged along the bladder wall and particle velocities were estimated using a phase-based autocorrelation technique. In order to exclude acquisitions with poor signal quality, a wave signal-to-noise calculation was adopted. Each acquisition with a SNR lower than 2 was excluded from analysis. A 2D Fourier transform-based technique was used to estimate the frequency-resolved wave speed of the propagating Lamb wave. To avoid the influence of outlier wave speeds, the median wave speed at each filling volume and frequency was used for curve fitting.

In some configurations, signal filters, such as a median filter and the like, may be applied to the particle velocity data prior to averaging. In some configurations, alternative statistical measures of central tendency, such as medial particle velocity, may be used for averaging.

Estimation of P was accomplished by optimization-based curve fitting applied to eq. (26). Fitting was applied using a simplex search method to find a minimum. If Θ_i(P) is the value of the left side of eq. (26) for the ith frequency, then P may be estimated by:

$\begin{array}{cc}{P}^{est}=\arg \min (\sum _{i=1}^N{\Theta }_{i}P*\stackrel{\_}{{\Theta }_i\left(P\right)}\}& \left(25\right)\end{array}$

The bar above Θ_i(P) denotes its complex conjugate, which may be used to ensure a real-valued quantity. The initial guess for minimization for P may be any appropriate value, in a the non-limiting example the initial guess was set to 100 mmHg. The frequency range used for fitting was 150-500 Hz. This range was selected in order to match the range used in in-vivo UBV measurements.

The nonlinearity parameter b was set to an assumed value of 5. This is the same value used for modeling vascular tissue and is comparable to the measurement of fibroglandular tissue and ex vivo liver tissue. It is a reasonable assumed value for soft tissue, even if not specialized to bladder tissue in particular. The reference volume (i.e. V₀ in eq. (10) and eq. (11)) was set to 150 ml, based on the minimum physiological volumes reported for porcine bladders upon excision.

Referring to FIGS. 6A-D, a non-limiting example application to a pressurized biological compartment is shown. FIG. 6A shows the compartment modeled as a hyperelastic spherical pressure vessel 650 with Lamb waves 652 propagating along the circumference 654. Internal pressure P 656 is shown as applying pressure pushing out on hyperelastic spherical pressure vessel 650.

FIG. 6B shows the acoustic pulse producing Lamb waves propagating forward and backward away from the pulse.

FIG. 6C depicts a 2D-FFT based method used to estimate frequency dependent phase velocities from both waves (only a single wave shown).

FIG. 6D shows curve fitting applied to the median wave speeds from the 5 acquisitions collected at each volume.

Referring to FIG. 7, non-limiting examples of vessel pressures plotted along with the vessel pressure measured by the pressure gauge as a function of bladder volume are shown. These plots illustrate how the analysis performs at recreating the underlying pressure from wave propagation and deformation (i.e. volume) measurements. Table I reports the corresponding root mean squared deviation (RMSD), mean error (ME) and relative mean error (RME) of estimated pressures relative to measured pressures for the experiments shown in FIG. 7. These quantities are defined respectively as

$\begin{array}{cc}RSMD=\sqrt{1/N\sum _{i=1}^N{\left(P_l^{\left(est\right)}-P_l^{\left(meas\right)}\right)}^2}& \left(26\right)\end{array}$ $\begin{array}{cc}\mathrm{ME}=1/N\sum _{i=1}^N\left(P_i^{\left(\mathrm{est}\right)}-P_i^{\left(\mathrm{meas}\right)}\right)& \left(27\right)\end{array}$ $\begin{array}{cc}\mathrm{RME}=100\times \left(1/N\sum _{i=1}^N\left(\frac{(P_i^{\left(\mathrm{est}\right)}-P_i^{\left(\mathrm{meas}\right)}}{P_i^{\left(\mathrm{meas}\right)}}\right)\right)& \left(28\right)\end{array}$

Where the added superscripts denote estimated (est) versus measured (meas) pressures for the ith volume and N is the number of estimates. The estimated pressures were broadly in agreement with the measured pressure (RMSD: 7.90 (mmHg); ME:−1.04 (mmHg); RME:−5.13 (%), though there is considerable variation between bladders (RMSD: 0.58-15.21 (mmHg); ME:-11.77-7.99 (mmHg); RME: −27.47-21.56 (%))

Referring to FIG. 8, a correlation plot of all measurements from the non-limiting example 8 bladders are shown. R2 is calculated as:

$\begin{array}{cc}{R}^2=1-\frac{\sum _{i=1}^N{\left({\mathrm{Pi}}_l^{\left(est\right)}-P_i^{\left(meas\right)}\right)}^2}{\sum _{i=1}^N{\left(\langle P_i^{\left(est\right)}\rangle -P_i^{\left(est\right)}\right)}^2}& \left(29\right)\end{array}$

Where < > denotes the mean of the quantity. The combined data align with an ideal one-to-one match between measured and inferred vessel pressure (R2=0:81).

TABLE I
Descriptive Statistics of the errors in estimated pressures calculated
for individual bladders and for all bladders combined.
	RSMD (mmHg)	ME (mmHg)	RME[%]
B1 (pre)	0.579	0.01	01 0.32
B1 (post)	5.86	5.31	20.53
B2 (pre)	5.13	−4.35	−22.09
B2 (post)	15.21	−7.1	−5.22
B3 (pre)	6.00	−5.63	−27.47
B3 (post)	15.68	−11.77	−20.48
B4 (pre)	4.36	−4.22	−26.44
B4 (post)	5.22	1.59	6.35
B5 (pre)	3.91	−3.78	−21.87
B5 (post)	6.03	4.07	8.78
B6 (pre)	3.42	−2.43	−16.33
B6 (post)	3.35	−1.84	−8.64
B7 (pre)	4.41	−4.15	−26.9
B7 (post)	11.56	7.99	21.56
B8 (pre)	8.95	2.96	7.83
B8 (post)	9.92	4.89	7.78
Combined (all)	7.90	−1.04	−5.13
RSMD = root mean squared deviation; ME = mean error; RME = relative mean errors; AND B# = bladder (#).

Material isotropy and a static fluid may also be assumed in some configurations, such as bladder applications. Other biological compartments, such as large blood vessels and the like, may include other assumptions. Cylindrical geometries may also be used with appropriate modifications to the non-limiting example expressions.

Effects of dynamic loading may also be included, such as by adding additional viscous dispersion and attenuation to the wave propagation. Incorporating the effects of dynamic loading may also produce large strain, viscoelastic stress relaxation.

Non-Limiting Example SNR Calculation

The signal-to-noise ratio of the propagating Lamb wave may be determined in some configurations, where signal and noise regions of-interests (ROIs) may be collected from the spatiotemporal particle velocity maps. The signal ROI may be collected by fitting a line along the peaks of the propagating waves and collecting an ROI as a window around the linear regression.

Referring to FIG. 9, non-limiting example signal and noise ROIs are shown. A robust linear fitting may be applied to the time-to-peak arrival times as a function of the lateral location. A noise region-of-interest (ROI) may contain the entire background except for the trough of the propagating wave, such that large peaks or troughs occurring anywhere in the spatiotemporal propagation domain outside of the signal ROI may be included in the noise ROI. Such large peaks act approximately as impulse functions that produce large ripples in the k-space representation of the wave. These ripples can substantially alter the wave speed estimation even when most of the wave propagation has a strong signal.

A time window may also be determined in some configurations. As a non-limiting example, a time window of 1.2 ms may be used to define the signal ROI instead of a spatial window, such as to limit the influence of wavelength variation on the SNR calculation, as wavelength will vary with wave speed at a given frequency. The SNR may then be calculated as:

$\begin{array}{cc}SNR=\frac{\langle S_i\rangle }{std\left(N_i\right)}& \left(30\right)\end{array}$

Where Si and Ni are the ith particle velocity values in the signal and noise ROIs respectively and < > and std( )denote the calculated mean and standard deviation

The present disclosure has described one or more preferred embodiments, and it should be appreciated that many equivalents, alternatives, variations, and modifications, aside from those expressly stated, are possible and within the scope of the invention.

Claims

1. A method for determining a kinetic value of a tissue volume, the method comprising: with a transducer, generating Lamb waves in a tissue wall of the tissue volume using radiation force, wherein the tissue volume is formed by the tissue wall that spatially separates a fluid material from a rigid material;forming ultrasonic echo data by detecting ultrasonic energy reflected by multiple locations along the tissue volume that is subject to the radiation force;determining frequency-resolved wave speed data for the generated Lamb waves from the ultrasonic echo data; anddetermining pressure of the tissue volume based upon the determined frequency-resolved wave speed data.

2. The method of claim 1, wherein determining the kinetic value includes determining pressure and includes optimization-based curve fitting a dispersion relation that includes the frequency-resolved wave speed data.

3. The method of claim 2, wherein the dispersion relation includes a geometry of the tissue volume.

4. The method of claim 3, wherein the geometry is at least one of a sphere or cylinder.

5. The method of claim 4, wherein an undeformed inner volume of the at least one sphere or cylinder is determined from a wall thickness of the tissue wall of the tissue volume.

6. The method of claim 1, further comprising determining axial particle velocities in the tissue wall using a phase-based autocorrelation technique.

7. The method of claim 6, wherein determining frequency-resolved wave speed data includes performing a 2D Fourier transformation of the determined axial particle velocities.

8. The method of claim 6, further comprising determining signal to noise ratios (SNRs) for axial particle velocity data and rejecting axial particle velocity data with SNRs below a determined threshold signal to noise ratio (SNR).

9. The method of claim 1, wherein the tissue volume is a bladder, and the tissue wall is a wall of the bladder.

10. The method of claim 9, further comprising adjusting a dispersion relation for a curvature of the bladder wall.

11. A system for determining a kinetic value of a tissue volume, the system comprising: a transducer configured to generate Lamb waves in a tissue wall of the tissue volume using radiation force and detect ultrasonic energy reflected by multiple locations along the tissue volume that is subject to the radiation force to form ultrasonic echo data, wherein the tissue volume is formed by the tissue wall that spatially separates a fluid material from a rigid material;a computer system configured to: i) determine frequency-resolved wave speed data from the generated Lamb waves from the ultrasonic echo data; andii) determine a kinetic value of the tissue volume based upon the determined frequency-resolved wave speed data.

12. The system of claim 11, wherein the kinetic value is a pressure and the computer system is further configured to determine the pressure using an optimization-based curve fitting of a dispersion relation that includes the frequency-resolved wave speed data.

13. The system of claim 12, wherein the dispersion relation includes a geometry of the tissue volume.

14. The system of claim 13, wherein the geometry is at least one of a sphere or cylinder.

15. The system of claim 14, wherein the computer system is further configured to determine an undeformed inner volume of the at least one sphere or cylinder from a wall thickness of the tissue wall of the tissue volume.

16. The system of claim 11, wherein the computer system is further configured to determine axial particle velocities in the tissue wall using a phase-based autocorrelation technique.

17. The system of claim 16, wherein the computer system is further configured to determine frequency-resolved wave speed data by performing a 2D Fourier transformation of the determined axial particle velocities.

18. The system of claim 16, wherein the computer system is further configured to determine signal to noise ratios (SNRs) for axial particle velocity data and reject axial particle velocity data with SNRs below a determined threshold signal to noise ratio (SNR).

19. The system of claim 11, wherein the tissue volume is a bladder, and the tissue wall is a wall of the bladder.

20. The system of claim 19, wherein the computer system is further configured to adjust a dispersion relation for a curvature of the bladder wall.