Experimental performance assessment of the sub-band minimum variance beamformer for ultrasound imaging

Recent progress in adaptive beamforming techniques for medical ultrasound has shown that current resolution limits can be surpassed. One method of obtaining improved lateral resolution is the Minimum Variance (MV) beamformer. The frequency domain implementation of this method effectively divides the broadband ultrasound signals into sub-bands (MVS) to conform with the narrow-band assumption of the original MV theory. This approach is investigated here using experimental Synthetic Aperture (SA) data from wire and cyst phantoms. A 7 MHz linear array transducer is used with the SARUS experimental ultrasound scanner for the data acquisition. The lateral resolution and the contrast obtained, are evaluated and compared with those from the conventional Delay-and-Sum (DAS) beamformer and the MV temporal implementation (MVT). From the wire phantom the Full-Width-at-Half-Maximum (FWHM) measured at a depth of 52 mm, is 16.7 µ m (0.08 λ ) for both MV methods, while the corresponding values for the DAS case are at least 24 times higher. The measured Peak-Side-lobe-Level (PSL) may reach -41 dB using the MVS approach, while the values from the DAS and MVT beamforming are above -24 dB and -33 dB, respec-∗

used with the SARUS experimental ultrasound scanner for the data acquisition.The lateral resolution and the contrast obtained, are evaluated and compared with those from the conventional Delay-and-Sum (DAS) beamformer and the MV temporal implementation (MVT).From the wire phantom the Full-Width-at-Half-Maximum (FWHM) measured at a depth of 52 mm, is 16.7 µm (0.08λ) for both MV methods, while the corresponding values for the DAS case are at least 24 times higher.The measured Peak-Side-lobe-Level (PSL) may reach -41 dB using the MVS approach, while the values from the DAS and MVT beamforming are above -24 dB and -33 dB, respec-

Introduction
Adaptive beamforming techniques have been used for decades in numerous applications of array processing [1][2][3][4] in fields such as sonar, radar, and seismology.The commercial use of such techniques is mainly related to military applications [5] or telecommunications [6].In general, adaptive beamformers aim to maximize the signal strength from a particular location and suppress signals from all other locations.This is accomplished by processing the received responses of an array to obtain constructive and destructive interference respectively.Improved transducers, reduced costs, and the availability of processing with Field-Programmable Gate Arrays (FPGAs) or Graphics Processing Units (GPUs) makes it possible to introduce similar real-time adaptive processing to medical ultrasound imaging [7,8].Initial results indicated that increased resolution and contrast can be achieved.Such research includes the linearly constrained adaptive beamformer [9,10], the adaptive beamformers suggested by Viola and Walker [11], and the Minimum Variance (MV) beamformer [12][13][14][15].The latter was originally developed by Capon [16] for use with seismic arrays with the objective of localizing earthquakes with greater precision.From a theoretical perspective, the MV beamformer is intended to provide unit gain in a selected direction and minimize the signal power for all other directions that are normally contributions from side-lobes.
The MV method has been extended unmodified to broadband ultrasound imaging, in the time-domain [17], or in the frequency domain [18] where division of transducer element signals into frequency sub-bands (MVS) precedes the processing.The fre-quency division ensures that the original narrow-band condition of the beamformer is met as laid out by Capon [16].As a result the MVS is expected to achieve improved resolution compared to the temporal implementation (MVT).In medical ultrasound imaging, the MVS was first introduced by Holfort et al. [18] with a quantitative evaluation on simulated data showing, by some measures, one order of magnitude higher image resolution compared to the conventional Delay-And-Sum (DAS) beamformer.
Applying lateral oversampling in simulated ultrasound data during the receive processing, resulted in further resolution gains [19].Particularly, the main-lobe width of a point target located at a depth of 40 mm was found to be 22 times narrower with MVS beamforming when compared to that achieved by DAS beamformers.A −13 dB sidelobe reduction was also noticed in favor of the adaptive approach.A 10-fold resolution improvement was maintained for point targets located at greater depths, up to 80 mm.
Further results from a circular cyst phantom showed that the MVS yielded 3 dB higher contrast compared to the best DAS beamformer, which also distorted the initial cyst shape.
The above simulation studies on MVS motivate the experimental validation.In an experimental setting, the cancellation of unwanted signals becomes less reliable, and the interference of adjacent targets is likely to compromise the accuracy of the method.Thus, in this work the MVS was applied to real ultrasound data from a wire-target and a cyst phantom.The MVS was combined with a Forward-Backward (FB) spatial smoothing technique [20], as it has been shown to increase the robustness of the timedomain MV beamformer implementation [21,22].Quantitative resolution and contrast metrics were used to evaluate the MVS performance and to compare it with the MVT, and the DAS beamformer, which is widely used in commercial ultrasound systems.

Methods
The standard way to process the signals received by a transducer array [23] is the DAS beamformer.The channel signals are time-delayed, weighted, and finally summed to form the beamformer output.The apodization weights depend on depth with a fixed F-number rather than on the data, and therefore expand with increasing depth.The MVS method, the experiment, and the quantitative analysis are described below.

Sub-band Minimum Variance Beamforming
The MVS method calculates a set of data-dependent apodization weights.This dependence on the acquired Radio Frequency (RF) data renders the beamformer adaptive.
where w(ω) = [w 0 (ω), w 1 (ω), ..., w M−1 (ω)] H is the complex weights vector, H is the vector of the Fourier Transform of the segmented channel signals, and {.} H denotes conjugate transpose.The MVS minimizes the power of each b(ω, r p ) corresponding to a single frequency bin, while preserving the signal from the position r p .The power is given by: (2) where E{.} denotes the expectation value and R(ω) is the covariance matrix given by: The MV objective can be expressed as: where, e is the time-delay vector that is only a vector of ones, since the time delays already have been applied to the signals.Lagrangian multiplier theory [24] yields an analytical solution to this constrained optimization problem.Given that R −1 exists, the MV weights are calculated by: The minimization goal is expressed for each frequency band, and the constraint refers to the distortionless response (unity gain) from the focus point [25,26].The MVS weight calculation is followed by the summation of the individual sub-band responses.
For K sub-bands, the final beamformer output B(ω, r p ) averaged over a number of N emissions, is: An important aspect of frequency domain implementation of the MV beamformer is the ability to calculate different weights for each sub-band and each point as seen from ( 6), which averages the processed channel data.After the MVS weight calculation, the inverse DFT is employed to derive the broadband MVS response, which for r p is centred around t = 0.

Forward-Backward Sub-array Averaging
A simple substitution of w into (1) would result in the calculation of the output of the MVS beamformer.While increased aperture size provides improved resolution, the increased number of channel data may result in inaccurate covariance matrix estimation, and thus incorrect weight calculation [17].To reduce the correlations between the received signals, the transducer array is divided into a number of overlapping subarrays, and the covariance matrix is replaced by the sample covariance matrix, which is estimated from several samples instead of the whole array.The sample covariance matrix may be derived by samples starting from the left of the array and moving to the right in Forward averaging ( RF ), or by the average of RF and RB , where RB is the averaging starting from the opposite direction (Backward averaging).The RF for a single frequency component, can be expressed as: where L is the sub-array length, and G l is the set of signals from the lth sub-array, in in [22], where J is the exchange matrix.In the Forward-Backward (FB) averaging technique the sample covariance matrix, RFB is given by: The FB averaging allows RFB to be inverted for larger L values than Forward only averaging does, making it possible to use larger sub-apertures during the processing.
The latter naturally increases the resolution limits.Once the optimized apodization weights, w, are calculated, with the use of the RFB , the beamformer output for each frequency bin, can be given by:

Experimental Setup and Data Analysis
The measurements were performed using the 1024 channel experimental ultrasound scanner SARUS [27], and all the parameters of the scans are summarized in Table 1.A 7 MHz, 192 element, linear array transducer with ≈ λ pitch was used to scan two phantoms containing wires and cysts respectively.In the first phantom, wires of a diameter of 0.07 mm were separated by 10 mm axially starting at a depth of 42 mm and reaching up to 122 mm.The speed of sound, c was measured to 1484 m/s based on the water temperature [28], resulting in a wavelength λ = c/ f 0 equal to 212 µm.
The cyst phantom contained a collection of different sized cylinders with diameters of In transmit, the active aperture consisted of 128 elements emitting a focused field.
The virtual source [29,30] was placed at a depth of 53.2 mm resulting in a F-number equal to 2, and Hanning transmit apodization was also used to reduce edge waves [31].
The lateral co-ordinate of the aperture centre was moved by a distance equal to one pitch between successive emissions, starting from the position of element #64 and ending to the position of element #128.RF data from 65 emissions in total were acquired from all 192 channels individually in receive, and were combined to provide a final high-resolution image as in standard Synthetic Aperture (SA) imaging [32].The MVS method was used to beamform a full image, by calculating an apodization weight for each image pixel.Synthetic aperture images using the MVT [17] as well as fixed Boxcar and Hanning [33] apodization weights with receive F-number equal to 1.5, were also formed as a standard for comparison.Adaptive apodization weights with L values ranging from 32 (= M/6) to 128 (= 2M/3) were extracted from the wire and cyst phantom data.For the wire phantom, areas of 6.4 mm in the lateral and 3.3 mm in the axial direction were beamformed separately.The selected areas included only one wire to avoid interference between neighboring scatterers and evaluate the effect of the beamformers on the side-lobes.The number of pixels in each image was 491×33, with small pixel lateral dimension of 13 µm (= pitch/16).The smaller pixel size increases the number of pixels and thus weights to be calculated, and was found to improve the lateral resolution when MV beamforming is used with point scatter data [19].Further decrease than the selected pitch/16 value in the lateral pixel size, did not result in additional lateral resolution improvements.For the cyst phantom, received data from all 65 emissions were used to form a complete image with dimensions 30 mm × 60 mm, with the same pixel size as in the wire phantom case.

Performance Assessment
Quantitative measurements on the acquired images were employed to evaluate the performance of the different beamformers.The lateral Full-Width-at-Half-Maximum (FWHM) and the Peak-Side-lobe-Level (PSL) were measured from the Point Spread Function (PSF) of an isolated wire.The lateral FWHM measures the width of the PSF main-lobe with narrower main-lobes indicating better resolution.The PSL quantifies the side-lobe suppression with lower values indicating contrast improvement.From the cyst phantom, the power ratio (PR), the contrast-to-noise ratio (CNR) and the speckle signal-to-noise ratio (sSNR) were used to assess the contrast resolution.The power ratio was calculated using the envelope detected image data by [18,34]: where P c is the mean power of a circular area inside an anechoic region (cyst) and P s the mean power of a circular area from the uniform scattering medium (speckle) of similar size.The CNR was calculated using the following equation [35,36]: where µ c and µ s are the mean intensity of a cyst and speckle at the same depth, and σ c and σ s are their corresponding intensity standard deviations.The sSNR was defined as µ/σ where µ is the mean value of the speckle amplitude and σ its standard deviation [36,37].The power in dB (y-axis) across the lateral beam width (x-axis) at a 52 mm depth is shown for all methods in Fig. 2. The values of the lateral FWHM and the PSL associated with this figure are displayed in Table 2.The lateral FWHM and PSL variation in respect to the different L values are shown in Fig. 3 for the wire-target located at a depth of 52 mm.For L = 32, the MVS results are comparable to those of the DAS beamformers (Table 2).The lateral FWHM varied between ≈ 0.3 mm and ≈ 0.02 mm, taking lower values at increasing L (Fig. 3(a)).The smallest value, and thus, best performance, was found for the largest L (= 128).The PSL was relatively constant around −20 dB for all L values up to 112 (Fig. 3(b)).The side-lobes dropped significantly to −41 dB only for L = 128, demonstrating, as in the FWHM case, the best image quality for L = 128.Further L increase resulted in noise-only images, from which no FWHM or PSL could be measured.The MVT results (Fig. 3 ,Table 2 showed no significant differences compared to the MVS, apart from a small difference in the PSL for the larger L values, where the MVS was at best 8 dB improved.The MVS implemented with L = 48 showed overall very similar performance to DAS Hanning with a 2 dB average difference, in favor of the adaptive approach.

Cyst Phantom
In this study the DAS and the MV methods were used to beamform an entire image instead of the isolated targets of the previous subsection.In Fig. 5   Visually the first 3 beamformed responses of the cyst phantom in Fig. 5, appear very similar, which was confirmed quantitatively (Fig. 6 and Table 3).At 30 mm depth, the PR was between −29 and −30 dB, the CNR between 1.94 and 2.05 and the sSNR between 2.16 and 2.27, demonstrating no significant improvement for the MVS.The three leftmost images also have two strong specular reflections at the top and bottom of the cyst.These characteristics are similar for all MV responses using L sizes between M/6 and M/2.The maximum sub-array length L = 128 used, which provided maximum resolution for the wire phantom (Fig. 1) was found to randomize the speckle appearance and therefore resulted in a varying intensity across the MVT and MVS images with alternating bright and dark vertical zones particularly at the top.Due to this intensity variation, the contrast at 30 mm was significantly reduced to    [19,36,[39][40][41].
The point scatterer results obtained using real data here, confirm the previous findings derived in a simulation environment.In this work, the MVS provided at best 24 times lower FWHM and −17 dB improved side-lobe suppression compared to DAS beamforming.These numbers are comparable to those mentioned in the simulation study (22 times and −13 dB respectively) [19].However, the experimental results have been acquired by deploying an optimized processing that involved a larger subarray length value (L = 2M/3 = 128) and target isolation.The use of such a high L value was enabled by using the FB averaging technique.It is commonly accepted that the FB averaging outperforms the standard forward averaging [21], providing a more robust sample covariance matrix.The forward only averaging is usually combined with sub-array lengths that are between M/4 and M/2 [18,38], since there is a trade-off between sub-aperture size and sample covariance matrix accuracy.Importantly, for L values smaller than 128, the MVS showed some resolution gains compared to conventional beamforming (Fig. 1), but the level of improvement was significantly lower compared to the simulation results [18].Moreover, the MVS beamformer was applied to small regions centred around a single wire to ensure that the highest possible performance is achieved.The beamforming of larger structures, minimized the resolution gain of the adaptive method as was demonstrated by the cyst data processing.From Fig. 5(a)-(c) and Table 3, it is not possible to identify a significant advantage of the MVS over the DAS.The deterioration of the MVS image in Fig. 5(d)-(e) is due to the larger sub-aperture used which, given the large number of scatterers that were included in this phantom, reduces the possibility of optimal signal cancellation.The cyst phantom results are not in full accordance with the initial simulations, where the circular shape of a cyst located at 40 mm depth was preserved with the MVS compared to the distorted DAS response [18,42].Recent MV studies on cystic resolution [39,43] show that it is only towards the edges of small cysts that the MV beamforming may result in higher contrast, which is not in disagreement with the results here.
A comparison of the MVT and the MVS, did not demonstrate a clear advantage of one implementation over the other.From the wire-target experiment, there is little difference in PSL between the two adaptive approaches, as in simulation [19].This is best reflected in the PSF appearance for the wire target closest to the virtual source (at 52 mm depth), where the target is more clearly defined for the MVS derived image (Fig. 1(e), first row), while side-lobes are visible in the MVT case (Fig. 1(d), first row).From the cyst phantom, the resulting values of all contrast resolution metrics are similar for MVT and MVS, while there was a −9 dB contrast improvement for the MVS in the simulation results [19].Overall, the expected theoretical advantage of dividing the broadband ultrasound signals into sub-bands, was not confirmed experimentally.However, as noted in [44], beamforming methods such as the MVS, are in general, sub-optimal since correlations between the frequency domain channel signals of different sub-bands are not taken into account in the derivation of the broadband beamformer output.Considering the additional computational load, which is attributed to the number of matrix multiplications needed for the weight calculation as shown from ( 5), (7), and (8) (proportional to L 3 ), and to the individual processing of each frequency band, it can be concluded that there is no clear benefit in using the MVS method in structural/anatomical imaging.Both wire and cyst phantom experiments confirm that the MV efficiency depends on the relation between the number of available channel signals and the number of scatterers to be resolved [17].The MV performance is likely to be further compromised when imaging structures of the human body by the tissue induced aberration [45], mainly due to the variations in the speed of sound [18] and attenuation [46].The MV beamformer would require further development to compensate for such environments.
The applicability of the MV method remains open for B-mode imaging, and despite the limitations described above, it has been shown that it is feasible to implement MV beamforming for real-time cardiac ultrasound imaging [7] or imaging of the eye [8].
The results here show that the MVS using L = 48 is a balanced MV implementation offering 33% improved lateral resolution compared to DAS, while also maintaining similar contrast resolution with lower than 5% deviation based in all the criteria selected for the quantitative evaluation, as shown in Table 3.However, the high-sub-array length MV implementations appear particularly attractive for use in point scatter imaging.The emerging field of super-resolution ultrasound contrast imaging is an obvious example.It is well established that single microbubbles are very efficient point scatterers [47], and recent developments have utilized this fact to explore techniques available from other fields of sensing.In essence all the techniques aim to locate the centre of a particle signal and minimize side-lobes.With the aid of such contrast microbubbles, and an a priori knowledge of point source scatter, high resolution transcranial images of vascular structure have been obtained [48].This was accomplished by applying aberration correction methods based on the position estimation of individual bubbles, thus achieving resolution beyond the diffraction limit.Similarly, based on the localization of isolated signals from microbubbles, in-vivo imaging of the mouse ear microvasculature with 5-fold resolution gains was performed with the additional feature of a super-resolution blood velocity mapping [49].In other work, improved microbubble localization with ultrafast Ultrasound Localization Microscopy (uULM) applied to conventionally beamformed data, resulted in the mapping of vessels up to 10 times smaller than the ultrasound wavelength, during in-vivo measurements on anaesthetized rats [50].Whereas super-resolution imaging is mainly image-based, the MV beamformer offers a complementary method in the processing of signals.The advantage of using such a method does not only rely on the narrower main-lobe width of a PSF (FWHM), but also in improved side-lobe suppression (PSL).This suggests the potential for reduced variability of the PSF and reduced background clutter or noise.
Both of these may improve the statistics of detecting microbubbles in an image, further improving accuracy and reproducibility of image processing, while also increasing the number of bubbles possible to use.The lack of axial resolution improvement using the MV method is not a major obstacle as the PSF has a very well defined shape, which may facilitate the image analysis implemented after the image formation.

Conclusion
The performance of the frequency domain implementation of the MV beamformer was experimentally examined for medical ultrasound imaging.The adaptive method provided up to 24-fold resolution gains and up to 17 dB improved side-lobe suppression over the conventional DAS beamformers in the lateral localization of individual point scatterers.A comparison with the time domain MV beamformer showed no difference in resolution and up to 8 dB improvement in the side-lobe suppression.These results were acquired using experimental ultrasound data from point scatterers, and confirmed previous simulation findings.Further, the adaptive method did not demonstrate its usefulness for entire images in a cyst phantom study, where the contrast resolution was at best similar to the one provided by the DAS beamformers.
The received channel data are focused as in a normal DAS beamformer to generate the input signals to the MVS algorithm.The short-time Fourier transform (STFT) is used to divide the time delayed channel signals into frequency sub-bands, and each band is thereafter processed separately.For a single focus point, r p , the Discrete Fourier Transform (DFT) is applied on segments with a period t d , hence STFT, transforming the time domain input signals into the frequency domain.The segment size depends on the excitation pulse and the 2-way impulse response of the transducer used.The mth segmented, channel signal y m (t) is given for t ∈ [−t d /2,t d /2].The beamformer output b(ω, p), for a single emission, for a transducer with M elements, that are all used in receive, and for each frequency sub-band ω, is given by:

8, 4 ,
and 2 mm at various depths starting from 10 mm to 60 mm (Dansk Fantom Service, Frederikssund, Denmark).The cyst phantom was homogeneous with a constant speed of sound equal to 1540 m/s, resulting in a wavelength equal to 220 µm.Data were initially sampled at 70 MHz, and then the sampling frequency, f s was decimated by a factor of 2 to 35 MHz.Averaging was used along with the decimation, through accumulation of successive samples, effectively implementing a rectangular filter with a sinc transfer function.

Figure 1 :
Figure 1: Responses of individual wire-targets at different depths are shown for 5 different sets of apodization weights as resulted from conventional beamforming (a) and (b), and MV adaptive beamforming (c), (d) and (e).A 40 dB dynamic range display was used.

Figure 2 :
Figure 2: Lateral variations of the beamformed responses of Fig. 1 (first row) at a depth of 52 mm.

Figure 3 :
Figure 3: Lateral FWHM and PSL variation in respect to sub-array length L, for 65 emission MVT and MVS responses.Sub-array length L values up to 2M/3 were used.

Figure 4 :
Figure 4: Lateral FWHM variation (a) and PSL variation (b) as a function of depth for 65 emission DAS and MV responses.Wire-targets between depths of 42 mm and 122 mm were imaged.
the beamformed responses of the cyst phantom are shown with a dynamic range of 60 dB.Similarly to Fig. 1 two MVS images are shown with sub-array lengths L = M/4 = 48 and L = 2M/3 = 128 and one MVT with L = 2M/3 = 128.In Fig. 6 the lateral variations at 30 mm depth are shown, and the images from the cyst at 30 mm depth are also displayed separately in Fig. 7 for more detail.The calculated contrast resolution metrics can be found in Table 3 for the 4 mm diameter cyst centred at (x, z = −1 mm, 30 mm) based on the yellow circled areas shown in Fig. 7(a).

Figure 5 :
Figure 5: Responses of the cyst phantom are shown for 5 different sets of apodization weights as resulted from conventional beamforming (a) and (b), and MV adaptive beamforming (c), (d) and (e).A 60 dB dynamic range display was used.

− 16
dB and −15 dB for MVT and MVS respectively.The corresponding CNR and
4. DiscussionA quantitative assessment of the Minimum Variance Sub-band (MVS) beamformer, using experimental ultrasound data was investigated for the first time.It was shown that such adaptive apodization weights achieve super-resolution lateral localization of point sources, with FWHM values of λ/12, while at the same time keeping the side-lobes below −40 dB.It is difficult to compare the above findings with other MV implementations due to the use of varying scan parameters, scanned object dimensions, or metrics definitions.However, to the best of the authors' knowledge such low FWHM values have never been presented in the MV beamforming literature for medical ultrasound.The MV processing (as opposed to the MVS) is mainly time domain-based and has provided λ/10 at best, for simulated data elsewhere

Table 1 :
Scan Parameters for the Wire-and Cyst-Phantom Measurements