If you don't remember your password, you can reset it by entering your email address and clicking the Reset Password button. You will then receive an email that contains a secure link for resetting your password
If the address matches a valid account an email will be sent to __email__ with instructions for resetting your password
Arterial pulse waves contain clinically useful information about cardiac performance, arterial stiffness and vessel tone. Here we describe a novel method for non-invasively assessing wave properties, based on measuring changes in blood flow velocity and arterial wall diameter during the cardiac cycle. Velocity and diameter were determined by tracking speckles in successive B-mode images acquired with an ultrafast scanner and plane-wave transmission. Blood speckle was separated from tissue by singular value decomposition and processed to correct biases in ultrasound imaging velocimetry. Results obtained in the rabbit aorta were compared with a conventional analysis based on blood velocity and pressure, employing measurements obtained with a clinical intra-arterial catheter system. This system had a poorer frequency response and greater lags but the pattern of net forward-traveling and backward-traveling waves was consistent between the two methods. Errors in wave speed were also similar in magnitude, and comparable reductions in wave intensity and delays in wave arrival were detected during ventricular dysfunction. The non-invasive method was applied to the carotid artery of a healthy human participant and gave a wave speed and patterns of wave intensity consistent with earlier measurements. The new system may have clinical utility in screening for heart failure.
found that the energy of the first wave was markedly reduced in stable compensated systolic heart failure and that the reduction increased with increasing severity on the New York Heart Association scale. The second wave was unaffected.
similarly found that the first wave was reduced by approximately 50% in patients with dilated cardiomyopathy (i.e., impaired cardiac contractility), without significant change in the second wave. They also found that the second wave was reduced by approximately 50% in patients with hypertrophic cardiomyopathy (i.e., impaired diastolic function), without significant change in the first.
Both studies employed the methods of analysis developed by
in which wave intensity, , is calculated as the product of and —respectively the change in blood pressure () and velocity () over a short interval in the artery of interest. If the wave speed is known (and it can be derived from and ), the waves can additionally be divided into their forward-traveling and backward-traveling components (
). Measurements of and should be temporally and spatially coincident, and they need to be made with high temporal resolution to capture the rapid changes in wave intensity that occur over the cardiac cycle.
These requirements, and particularly the need for frequent pressure measurements, are problematic.
used Doppler ultrasound to obtain U and tonometry to obtain P, both in the carotid artery. Tonometry requires complex calibration to give true pressures, and the measurements of and cannot be made simultaneously.
used Doppler ultrasound to measure and M-mode ultrasound to measure arterial diameter , and assumed that is proportional to . This ignores the well-established non-linearity of arterial stress–strain curves. A technique that avoids these problems is to measure and with an intra-arterial catheter that has both a pressure transducer and a Doppler probe at its tip (
Modeling of Congenital Hearts Alliance (MOCHA) Collaborative Group. Reduced ascending aorta distensibility relates to adverse ventricular mechanics in patients with hypoplastic left heart syndrome: Noninvasive study using wave intensity analysis.
), the difference in optimal beam angles for Doppler and B-mode or M-mode imaging and the high cost and long acquisition time of MRI.
In the present work, we describe a non-invasive method in which and are both obtained from the same B-mode ultrasound images; singular value decomposition (SVD) is used to separate the weak blood and strong tissue signals, cross-correlation between images to track the moving wall and blood speckles and an ultrafast scanner with plane-wave transmission to adequately resolve the rapid acceleration and deceleration of blood during systole (
). The method is validated by comparison with data obtained using the invasive catheter-based method and is found to be capable of detecting ventricular dysfunction in rabbits. Preliminary data indicate that the technique can be translated to human participants.
All experiments complied with the Animals (Scientific Procedures) Act 1986 and were approved by the Animal Welfare and Ethical Review Body of Imperial College London. Twelve New Zealand White rabbits (HSDIF strain, Envigo, Huntingdon, UK; age 81 ± 13 d, weight 2.39 ± 0.38 kg, mean ± SD) were maintained on a 12:12 h light:dark cycle at 18°C and fed a normal laboratory diet and tap water ad libitum. Each animal was pre-medicated with acepromazine (0.5 mg/kg i.m., CVet, London, UK), anaesthetized with fentanyl fluanisone (Hypnorm, 0.3 mL/kg i.m., Janssen, High Wycombe, UK) and midazolam (Hypnovel, 0.1 mL/kg i.v., Roche, Welwyn Garden City, UK), tracheotomized and ventilated at 40 breaths/min with 45 cm H2O peak inflation pressure (Harvard Small Animal Ventilator, Harvard Apparatus, Cambourne, UK). Anaesthesia was maintained with fentanyl fluanisone (0.1 mL/kg) and midazolam (0.1 mL/kg) approximately every 45 min. Body temperature was monitored using a rectal probe and maintained with a heated pad. Blood oxygen saturation was monitored by pulse oximetry.
A Vantage 64 LE ultrasound research platform (Verasonics, Kirkland, WA, USA) equipped with an L11-4v probe was used to collect high-frame-rate images of the abdominal aorta distal to the renal arteries. The rabbit was tipped from its supine position to elevate its left side, and the probe was angled in the opposite direction and clamped to give a stable, approximately coronal view through the aortic centerline. Sixty-four elements located at the center of the probe were used to transmit and receive, giving a lateral field of view of 19.2 mm, along the axis of the aorta.
B-Mode imaging was performed using a coherent plane wave compounding scheme (
). A pulse sequence consisting of three plane waves (8 MHz, 1 cycle), steered from –5° to 5° in three 5° steps, was transmitted at a pulse repetition frequency of 47.6 kHz for an imaging depth of 20 mm. Successive pulse sequences were separated to give an effective frame rate of 1 kHz after coherent summation of the three angled plane waves. Images were acquired for 2 s. The mechanical index (MI) was 0.1 MPa. Backscattered radiofrequency (RF) signals were sampled at four times the transmit frequency; echoes of each transmission were recorded for post-processing offline.
A three-lead electrocardiogram (ECG) signal was recorded concurrently via a PowerLab 26T data acquisition system (AD instruments, Oxford, UK) connected to the Verasonics host computer running LabChart software. An analog trigger signal sent at the start of the image acquisition sequence from the Verasonics to the PowerLab enabled later time alignment of the arterial waveforms and ECG trace.
The proximal right femoral artery was exposed, lidocaine was administered topically (typically three sites, <1 mL total) and the femoral nerve was severed.
A Volcano ComboWire XT wire with pressure and Doppler sensors at the tip, connected to a ComboMap system (Phillips-Volcano, San Diego, CA, USA), was inserted into the artery and advanced into the abdominal aorta, just distal to the imaging site. The wire was rotated until the strongest Doppler signal was obtained. The wall filter was set at 400 Hz and the signal-to-noise threshold was adjusted manually to optimise the real-time Doppler envelope tracking.
The analog pressure and velocity output ports were connected to two analog inputs on the PowerLab data acquisition system (1 V = 100 mm Hg for pressure, 0.5 m/s for velocity, sample rate 1 kHz). The Doppler sampling rate on the ComboMap system is 120 Hz.
At least three data sets were acquired with the image-based (Verasonics) and catheter-based (Volcano) systems. The catheter-mounted flow sensor was disconnected during imaging as it caused interference, and the two techniques were therefore used alternately. Catheter-derived pressure data could, however, be recorded during the ultrasound imaging.
Following these measurements, the effects of ventricular dysfunction were assessed. Esmolol was administered to 10 rabbits; each received three boluses of increasing dose (1.5, 3.0 and 6.0 mg/kg, i.v.) at 10-min intervals. Successive data sets were acquired with both systems every minute, starting just before the first dose. Two rabbits were given an equivalent volume of saline intravenously as a control.
Post-processing of invasive pressure data
All post-processing was performed in MATLAB R2020b (The MathWorks, Natick, MA, USA).
Signal processing within the Volcano Combowire pressure-measuring system introduces a delay of 20 ms. This was determined in an elastic tube model in vitro by aligning the peak pressure recorded using the Volcano system with that simultaneously recorded using a high-fidelity, zero-delay catheter system (Fig. S1, online only).
Even after peak alignment, the foot of the Volcano-derived pressure waveform was less acute and occurred ∼5 ms earlier than the foot in the high-fidelity, catheter-derived waveform when both were employed in the in vitro system (Fig. S1). This difference was also observed when comparing Volcano-derived pressure with diameter measured by the non-invasive ultrasound method. Analysis of the power spectra revealed low-order, frequency-based filtering in the Volcano system. A filter with comparable characteristics was constructed; it was applied to measurements of and when processing the data with PU and ln(D)P methods. This adjustment makes data from the different sources directly comparable and avoids anomalous behaviour such as non-linear slopes of the ln(D)P and PU loops in early systole; it does not have a substantial effect on mean wave speeds or intensities (see, e.g., Fig. S2, online only).
In vivo, there was an additional time delay because of displacement of the Combowire tip from the center of the B-mode field of view. The delay varied between rabbits (from 0.5 to 5 ms) because of differences in both the displacement and the wave speed. This rabbit-specific delay and the system-specific 20-ms lag described previously were corrected together by assuming that reflections and viscoelastic effects dissipate to negligible levels prior to systolic ejection, meaning that and should initially rise together. Following
, the linearity of the systolic portion of the ln(D)P loop was optimized by applying arbitrary lags to and finding the lag that gave maximum correlation. ( had the aforementioned frequency filtering applied first.)
Post-processing of non-invasive ultrasound data
Parts of the following are based on earlier studies and validations by
). The image filtering and particle tracking methods are illustrated in Fig. S3, online only.
Beamforming and compounding
The RF channel data were reconstructed using delay-and-sum beamforming (assuming a 1540 m/s speed of sound for delay calculations) and then Hilbert transformed. Echoes from the three angled plane waves were coherently summed, resulting in a final ensemble image size of Nx = 128, Nz = 416 and Nt = 2000 (lateral pixel size = 0.15 mm, axial pixel size = 0.048 mm).
Singular value decomposition filtering
Singular value decomposition (SVD) filtering was used to separate the blood signal from the tissue signal and noise (
). The image was cropped so that its axial dimension was slightly larger than the vessel diameter. Then the ensemble image was reshaped into a 2-D spatiotemporal representation (size Nx·Nz × Nt) and decomposed as
where are the ordered singular values of , and and are the corresponding spatial and temporal singular vectors. Tissue, blood and noise have different spatiotemporal characteristics—the first (and largest) singular values typically correspond to tissue, the next to blood and the smallest to noise. The decomposition can be written as
where and are singular value thresholds. These were selected manually, based on the visual appearance of and the frequency of the temporal vectors (
). Figure S4 (online only) illustrates how the thresholds were selected and that their precise value is not critical to the derived velocity waveforms. Finally, images were envelope detected and re-sampled using MATLAB's griddedInterpolant() function to obtain isotropic pixel dimensions.
Velocity was measured by tracking the movement of blood speckles between frames in . Briefly, an “interrogation window” was defined in the vessel lumen in one image, and the cross-correlation of its speckle pattern with the pattern in a larger window in the next image of the sequence was calculated. The highest correlation identified the most likely location to which the group of scatterers had moved. Velocities were obtained from the displacement and the imaging frame rate.
Correlation was Fourier-based for computational efficiency. Consider a template window and search window where is the frame number. Cross-correlation in the spatial domain is equivalent to convolution in the frequency domain eqn (4)
where and are the Fourier transforms of the interrogation windows, and are pixel shifts in the lateral and axial directions, * is the complex conjugate and is the inverse Fourier transform. Conventionally, the displacement () is given by the location of the single Gaussian-shaped peak in the correlation function (arg max(|r|)); for this to be evident, velocity gradients within an interrogation window must be negligible.
The formulations for wave intensity assume velocity is one dimensional; we were consequently not interested in obtaining the velocity profile across the vessel, and a single interrogation window was therefore used. It was preferred to many smaller windows because it increased computational speed. The window was sized to span the entire lumen unless otherwise stated.
A distribution of velocities within a window causes the correlation peak to become broadened and skewed (
). The correlation coefficient, , represents the convolution of the “particle image” (autocorrelation of an image) with the probability density function of the velocity field within the window; the location of the maximum thus represents the modal displacement rather than the mean (
Assuming negligible mean displacement of blood in the direction, the centroid in the direction was computed as
where for a window size of , the summation ranges from to and . However, the correlation plane contains noise, and thus, integration limits had to be set for eqn (5); for large displacements, the broadened peak is shifted toward the edge of the plane, introducing bias. First the particle image diameter () was measured from the autocorrelation map (24 pixels, as measured using the 1/e2 rule [
]). For non-reversing flows, we assumed that flow at the wall is zero (no-slip condition) and increases toward the center. Thus, the integration limits were defined as pixels to the maximum displacement () . To determine , the modal displacement was found. Finally, velocity was calculated as Δxc/dt, where is the time between compounded frames (1 ms).
Pixel locking is a known problem of centroid calculations (
) but is not an issue if is >3 pixels, as is the case here. Precision of the centroid calculation was increased by interpolating the correlation map; this was done by zero padding in the frequency domain before inverse transformation.
Interrogation windows were centered on the part of the vessel where the diameter measurements were made. The window width was set to 80 pixels in the first image, and the search window to 80 plus twice the expected maximum pixel displacement in the second image (we anticipated velocities would not exceed 1 m/s [
]); this ensured the tails of the broadened peak were captured for large displacements.
Cross-sectional mean velocity estimation
The methods described in the previous section provided a measure of the mean velocity () in a 2-D lengthwise slice of the lumen that included the centerline of the vessel. Although the maximum velocity in the 2-D slice is the same as the maximum in the 3-D vessel, that is not necessarily true of the mean velocity. For Poiseuille flow in a cylindrical tube, where the velocity profile is a paraboloid, it is well established that the mean velocity is 0.5 times the maximum. However, in a thin 2-D slice through the centerline of such flow, for which the velocity profile is a parabola, the mean is 0.75 times the maximum. (There is comparatively less of the slow-moving, near-wall fluid.) Real arterial flow is pulsatile and more nearly of the Womersley type than the Poiseuille type.
Wave intensities derived using the 2-D velocity and Womersley-based 3-D velocity were nearly linearly related (Fig. S5, online only). The former was used in calculations of wave intensity, because only relative values were required to compare P-U and D-U methods or ventricular function before and after drug administration. However, absolute numbers were required when calculating wave speeds, to permit comparison with previous values, and was therefore used. It was calculated as described in Appendix S1 (online only).
The artery was assumed to change dimension only in the radial direction. Diameter waveforms were computed by 1-D cross-correlation of successive frames in . First a region of interest (ROI) with a width equal to 5 A-lines and a height of 30 pixels was selected on the anterior wall of the aorta. The signal values within the ROI were averaged laterally, and the location of the maximum corresponding to the inner layers of the wall was identified from the envelope. This was further refined by subpixel Gaussian fitting.
One-dimensional cross-correlation functions were computed for each lateral location, according to equation 4 with and vectors of equal size. The resulting correlation functions were averaged, and the maximum was identified to give the axial displacement, again using subpixel Gaussian fitting. This was repeated for the remaining image pairs, with the window offset by the previous displacement each time. The wall motion waveform was given by cumulatively summing the displacements over time. The process was repeated for the posterior artery wall. The diameter change waveform was then given by the difference between the anterior and posterior waveforms; addition of the initial diameter measured in frame 1 gave absolute diameter D.
Wave intensity analysis
Waveforms were smoothed using a second-order Savitsky–Golay filter with window length of 19 ms, equivalent to 19 samples for all data sets. The length was optimized by trial and error. Successive beats in the pressure, diameter and velocity waveforms were ensemble averaged using the R wave of the ECG as a fiducial marker; only complete cardiac cycles were included.
Net wave intensities were calculated as
where the prefix “n” refers to the non-invasive diameter formulation, and is the time between two consecutive samples. was normalized by to make the magnitude of independent of the sampling frequency (
The peak magnitudes, area under the curve (“wave energy”) and timings with respect to the R wave of the ECG were calculated for the dominant waves. The start and end of the waves were identified by the zero crossing of the wave intensity.
Wave speed was calculated using single-point “loop methods.” In plots of ensemble-averaged (or ) against , there is a linear portion of the loop in early systole whose gradient can be used to determine according to
was also calculated when and were recorded simultaneously. Not only is this value free of potential velocity errors, but the derivation of the equation does not assume unidirectional wave travel.
Values of and were obtained by non-invasive ultrasound of the carotid artery of a healthy young subject, using the parameters described previously except that the MI was increased to 0.4. Data were aligned and ensemble averaged using the R wave of the ECG as a datum.
The experiment was part of a trial approved by the National Health Service (IRAS ID 248724, HRA REC Reference 18/SC/0563) and the Imperial College Healthcare NHS Trust, and registered as ISRCTN 41232. Informed consent was given.
Values of p were obtained using Student's paired t-test unless otherwise stated. Limits of agreement in Bland–Altman plots were computed using the SimplyAgree package in R (The R Project for Statistical Computing) to account for repeated measures.
In the following, where waveforms are plotted against time, t = 0 corresponds to the R wave of the ECG. Individual animals are identified by four-letter codes.
Figure 1a and 1b illustrate, respectively, invasive pressure measurements and non-invasive diameter measurements for one animal (POEK); data for all animals are provided in Figure S6a and S6b (online only). The two types of waveform clearly have the same shape, but the relation between them is not expected to be perfectly linear: the elasticity of the vessel decreases as diameter increases (“strain stiffening”), and viscous effects cause a lag between changes in pressure and changes in diameter. The non-linearity and hysteresis are visible when is plotted against (Fig. 1c and Figure S6c [online only]); data for were filtered as described previously. Rabbit PRFH had high pressure, anomalous waveforms and a stiff vessel.
Figure 2a provides the corresponding invasive and non-invasive velocity measurements; data for all rabbits are given in Figure S7a (online only). For this comparison, the velocimetry interrogation window was reduced to a small box at the center of the vessel, as the Doppler system identifies the peak velocity, . (For the same reason, Doppler-derived velocities can be unsuitable when calculating absolute wave speeds; see
). This reduced the spread of velocities within the window, so pixel displacements could be determined using simple Gaussian peak fitting.
It is apparent in Figure S7a that the invasive and non-invasive velocity measurements, while in broad agreement, show a systematic difference: the invasive curves resemble smoothed versions of the non-invasive curves. This is particularly apparent at the systolic foot: the change in gradient is visibly sharper in the non-invasive results for many rabbits. The sampling frequency of the Volcano Combowire system and the additional low-pass filtering mean that the system cannot capture the fast early-systolic accelerations in the rabbit. There additionally appears to be a scaling error in rabbits PPEX, PPFD, PRDT and PRFH, which was most likely caused by the difficulty of perfectly aligning the orientation of the catheter with the mean flow direction: the Doppler probe at the catheter tip obtains velocity information only along the beam axis. Because of this issue, B-mode–derived velocity data are used in all subsequent analyses.
Figure 2b and Figure S7b (online only) additionally illustrate velocity waveforms calculated in three ways from an interrogation window spanning the entire vessel: was obtained assuming 2-D Poiseuille flow, by assuming a 3-D Womersley flow (in both cases the displacement being estimated from the centroid of the cross-correlation peak), and by using a Gaussian fit to find the cross-correlation peak.
Figure 2c demonstrates a linear relation between velocities derived from catheter-based Doppler and B-mode ultrasound measurements for all rabbits.
Figure 3a–c illustrates calculations of wave speed, , by three different “loop” methods—the PU loop, ln(D)U loop and ln(D)P loop—for a single rabbit (POEK); data for all rabbits are given in Figure S8 (online only). In Figure 3a and 3c, and , obtained from the non-invasive ultrasound methods, were filtered as described previously. The wave speed given for each rabbit is the mean of the wave speed obtained in each repeat, not the wave speed computed from the mean loop. All methods have anomalous results for rabbit PRFH, which had a stiff vessel.
Figure 3d–h illustrates the relation between the values obtained by the different methods for all rabbits. Compared with the ln(D)P method, which is considered a gold standard, the PU method consistently underestimated wave speed, and the ln(D)U method consistently overestimated it. However, the mean errors were small: –0.55 and +0.36 m/s, respectively, when compared with a best estimate of 5.02 m/s (the mean value obtained with the ln(D)P method). The errors do not appear to scale with wave speed but remain relatively constant for all values obtained with the ln(D)P method.
Figure 4a and 4b illustrate wave intensities calculated by both the PU and DU methods for a single rabbit (POEK). In Figure 4a, is filtered as described previously. Equivalent data for all rabbits are illustrated in Figure S9a and S9b (online only). Note that these intensities are the sum of forward-traveling and backward-traveling wavelets (with positive and negative intensities, respectively) at each time point. Both methods exhibit three waves during the cardiac cycle. There is a large net forward-traveling compression wave (“W1”) followed immediately by a smaller net backward-traveling wave (“R”) and then a further net forward-traveling expansion wave (“W2”) before a return to baseline. The backward-traveling wave comprises reflections from distal sites.
Figure 4c and 4d show comparisons between the two methods for peak intensity and wave energy (the area under the intensity curve). There is a clear positive relation between the two methods for each wave. Rabbit PRFH is again an outlier.
Figure 5 illustrates the effect on wave intensity of increasing doses of esmolol, a short-acting cardioselective β1-adrenergic blocker that reduces the force and rate of cardiac contraction. Wave intensity was assessed by the invasive and non-invasive methods; illustrative data, plotted in both 2-D and 3-D, are given for rabbit PFLH. Transient dips in wave intensity are evident as a color change in the 2-D plot and as a color change and decreased height in the 3-D plot, shortly after each dotted white line representing administration of esmolol. The dips increase in size with increasing dose. A delay in the arrival of the W1 and W2 waves is also evident after the administration of each bolus—in the lefthand panels, there is a rightward shift in the peaks below each white dotted line, with the largest delay occurring at the highest dose.
In Figure 6a–d are the average peak intensities and energies of the W1 and W2 waves, calculated by both the PU and DU methods. Mean data for nine rabbits are presented ±1 SD to illustrate variability within the sample. (Rabbit PPBT was not included in the analysis as it required an anesthetic top-up during the period of esmolol administration, rendering the results unreliable.) DU-derived intensities and energies are also given for two control rabbits administered vehicle rather than esmolol.
Dips in intensity and energy are evident after each esmolol bolus for both waves and both methods. No such trends are seen in the control animals. Statistical tests were performed for the experimental group by comparing the first measurement after each dose of esmolol with the preceding measurement. The drop in W1 peak intensity and energy was significant at all doses with both methods. The drop in W2 wave energy was also significant in all these cases. For W2 peak intensity, the effect was not significant (p > 0.05) for the lowest esmolol dose assessed with the PU method and for all doses with the DU method. The fact that statistical significance was obtained for wave energies but less consistently with peak intensities is likely due to alignment issues (see the Discussion).
Average variation in the time of arrival of the W1 wave, relative to the R wave of the ECG, is illustrated in Figure 7a. (Data for rabbits PPEX and PRDT at t = 25 min were omitted because of a poor ECG signal.) A clear delay is seen after each bolus, and the effect was significant at all doses. Delays of similar magnitude, again all significant, were seen in the interval between the W1 and W2 waves (Fig. 7b). No such changes were seen in the control animals.
No significant differences in wave speed were detected after each esmolol bolus with either the PU or ln(D)U loop method (Fig. S10 [online only], p > 0.05). Mean pressure was also not significantly altered (data not shown).
Figure 8 illustrates the feasibility of using the non-invasive method in the carotid artery of human participants.
An advantage of calculating wave intensity from arterial diameter and blood velocity is that both can be measured using non-invasive techniques. We developed methods based on ultrasound because its low cost and real-time imaging capability increase clinical utility. Wave intensities are computed from derivatives, which increases the requirement for measurement accuracy and precision. We avoided the use of Doppler ultrasound because of its fundamental limitations in quantifying blood velocity and because accuracy can be further compromised by the difficulty of employing the optimal, orthogonal beam angles for the Doppler and diameter measurements. Instead, changes in diameter and velocity were obtained from the same, successive B-mode images. This precludes the alignment issues apparent when using the current “gold standard” invasive system, which relies on separate pressure and flow transducers.
The new method presents several technical challenges. Blood is a poor ultrasound scatterer at the low frequencies required to penetrate tissue to the depth of conduit arteries. That necessitated the use of SVD to obtain separation of the blood signal from the tissue signal and noise. SVD is computationally intensive, and the manual selection of thresholds was slow. Randomized SVD is faster and has been realized in real-time applications (
). We found that precise selection is not critical: the lower threshold () could be changed between 115 and 150 and the upper threshold () could be omitted entirely without visible effect on the resulting velocity waveform (Fig. S4c, online only).
A further technical challenge is capturing the fast-moving arterial blood with sufficient temporal resolution. That was achieved by using an ultrafast, plane-wave scanner capable of frame repetition rates of ≈10 kHz. In such systems, the beam is not focused, which accelerates scanning but degrades the image. Image quality is recovered by averaging frames and by compounding scans at different angles. The net frame rate still exceeds that of conventional, focused scanners.
Velocity estimation is also complex. The use of B-mode imaging means that scattering from red blood cells can theoretically be tracked in two dimensions. In conventional ultrasound imaging velocimetry, a small field of scatterers is identified within the vessel in one image, and then a search is made for a field in the succeeding image that has the best-fit pattern of scatterers; the displacement of the field divided by the time interval between the two scans gives the velocity for that region of blood. Repeating the cross-correlation procedure for other fields across the diameter of the artery gives the velocity profile. We modified this procedure: a single field covering the whole diameter of the vessel was used. This reduces noise to a practicable level and increases the speed of computation. It is also consistent with the 1-D formulation of wave intensity. However, the resulting increase in the spread of velocities within each field introduces complexity. Velocity gradients are higher near the wall than near the center of the vessel, and the pattern of scatterers will therefore be more disrupted near the wall as the region of blood moves down the vessel. As a result, the cross-correlation algorithm is biased toward tracking scatterers near the center. Because absolute velocities are higher at the center, this leads to an overestimate of the mean flow velocity. In our method, this effect was ameliorated by using the centroid rather than the peak of the distribution of correlation coefficients
An additional issue is that the ultrasound image approximates a 2-D longitudinal slice through the center of the vessel. (In reality, this slice has finite thickness but that is ignored here.) Considering the simple case of fully developed flow, the velocity profile would be a parabola in the 2-D slice but a paraboloid in the 3-D vessel. The mean velocity for the parabola is three-fourths of the maximum velocity, whereas the mean for the paraboloid is half the maximum velocity. In the real vessel, where Womersley-type flow occurs, the velocity profile is more complex and varies over the cardiac cycle. For calculations of wave speed, where absolute rather than relative magnitudes matter, we used the inverse Womersley solution; elsewhere, to minimise computational cost, the 2-D values were used.
The assessment of arterial diameter over the cardiac cycle also uses the cross-correlation technique but its application in this case is less complex because the wall moves more like a rigid object, with a single velocity at each time step. The method is routinely employed (e.g.,
). The movement of the wall can be less than a single pixel; again, this problem arises from the need to use relatively low ultrasound frequencies, but subpixel resolution is readily obtained by fitting a curve to the Gaussian cross-correlation function.
The calculation of wave speed from the gradient of the PU loop or the ln(D)U loop assumes that all waves are forward traveling. This condition is most likely to be satisfied during the period when pressure and velocity rise at the start of ventricular ejection. Rabbits and people have broadly similar wave speeds, but the condition is more likely to be violated in the rabbit because of its higher heart rate, meaning there is less time for reflected wave energy to dissipate, and its shorter arteries, which reduce the time before reflections arrive at the measuring site. For these reasons, the present study used the first 10 ms of ejection in rabbits and the first 50 ms in the human participant.
Plots of aortic pressure, obtained with an intra-arterial catheter, and aortic diameter, obtained using the B-mode cross-correlation method (Fig. 1; Fig. S6, online only), both revealed in all animals the expected shape of a rapid systolic rise, followed by a slower decay that was clearly biphasic; that is, there was a secondary, diastolic peak or at least an inflection in the gradient of the decay, before the next systole, and not a continuous exponential decay as might be expected in the absence of wave reflections. The relation between pressure and diameter revealed the expected non-linearity and hysteresis caused by strain-stiffening and viscous effects.
Plots illustrating blood flow velocities obtained both by the Doppler ultrasound probe mounted at the tip of the intra-arterial catheter and by non-invasive velocimetry (Fig. 2; Fig. S7, online only) gave the expected triphasic aortic waveform in all animals and with both techniques: there was a period of forward-going flow in early systole, followed by a period of nearly zero flow and, finally a period of low, forward-going flow in diastole. Furthermore, the shapes of the curves obtained by the two methods were similar. A visible discrepancy between the methods was smoothing of the peaks and troughs in the Doppler measurements. The intra-arterial flow measurements were made with a device designed for human use; its sampling rate is probably insufficient to capture fine features of the shorter cardiac cycle in the rabbit. It also has inbuilt, low-pass filtering. A second discrepancy was a difference in the absolute magnitude in four rabbits: the invasive measurements gave lower peak velocities than were obtained with the non-invasive method. We attribute this to the well-known difficulty in aligning the Doppler probe with the predominant flow direction (
). With this exception, the scatterplot of the two data sets revealed a strong relation, the data lying close to the line of identity.
Wave speeds derived from the first 10 ms of the ln(D)U loop were higher than those derived from the equivalent period of the PU loop (Fig. 3; Fig. S8) but the difference was relatively small: averages were, respectively, 5.36 and 4.48 m/s. Reflections may have been present during this part of ventricular ejection. They introduce errors of opposite direction in the two methods: wave speed is overestimated when using the ln(D)U loop but underestimated when using the PU loop (
)—supports this view as the latter results lie between the other two sets of values (Fig. 3). A similar effect would have been produced by overestimating U. Note that wave speed is not required when assessing net wave intensity or net wave energy, but only for separating waves into their forward and backward components and for investigating vessel wall stiffness.
The PU and DU methods gave comparable patterns of wave intensity (Fig. 4). In all rabbits, both indicated a large net forward-traveling (“W1”) wave in early systole that was followed by a small net backward-traveling (“R”) wave and then a more complex and dispersed net forward-traveling (“W2”) wave. These are thought dominantly to correspond to the initial compression wave caused by ventricular contraction, its reflection and the forward-going expansion wave that derives from an interaction between the inertia of the blood and the reduced and then reversed contraction of the ventricle. The more disperse and complex nature of the W2 wave could be explained at least in part by the longer period between it and the R-wave of the ECG, used to time align individual heart beats when creating an ensemble average. This, coupled with beat-to-beat variation in the duration of the cardiac cycle would have the effect of smearing out a wave even it were concentrated at the same relative location in each individual beat. Additionally, the W2 wave may consist of more components than the other waves. Figure 4 reveals a strong correspondence between wave intensities measured by the two methods and also between wave energies (i.e., the integrals of the intensity curves for each wave). The relation in both cases is curvilinear, as expected from the strain-stiffening behaviour of the wall.
Ventricular dysfunction was transiently induced with esmolol, a short-acting cardioselective β1-adrenergic receptor blocker that decreases the force and rate of ventricular contraction. In addition to its clinical use, the drug has been employed in studies investigating the potential performance of left ventricular assist devices in heart failure (e.g.,
). Both the PU and DU methods were unequivocally able to detect decreased intensity, decreased energy and a delayed arrival of the W1 wave at all the esmolol doses employed (Fig. 5, Fig. 6, Fig. 7). Both could also detect similar effects on the energy and arrival time of the W2 wave (Fig. 5, Fig. 6, Fig. 7). The DU method was less sensitive than the PU method at detecting the esmolol-induced decrease in peak intensity of the more spread out W2 wave (Fig. 6). Neither method achieved a significant result at the lowest dose, suggesting that intensity is less consistent than the energy for this wave.
The new technique, using the diameter–velocity formulation of wave intensity analysis and based on B-mode rather than separate Doppler and M-mode measurements, is sufficiently accurate and non-invasive that it may have clinical utility. The rabbit aorta was chosen as a test bed because its diameter, depth and blood velocity are similar to those of human peripheral arteries, Furthermore, rabbits have wave patterns similar to those of people, unlike small rodents, for example. A proof-of-concept experiment indicated that the non-invasive method can be applied to the human carotid artery. The pattern of wave intensities and the calculated wave speed (Fig. 8) were consistent with data from other non-invasive techniques (
). Potential uses include screening for heart failure, tracking its progression and improving the accuracy of prognoses. Although the algorithms required to achieve accurate results have proven complex, they are not difficult to implement now that they have been developed. A recent modeling study (
Modeling of Congenital Hearts Alliance (MOCHA) Collaborative Group. Reduced ascending aorta distensibility relates to adverse ventricular mechanics in patients with hypoplastic left heart syndrome: Noninvasive study using wave intensity analysis.