 Research article
 Open Access
 Published:
Deterministic compressive sampling for highquality image reconstruction of ultrasound tomography
BMC Medical Imaging volume 17, Article number: 34 (2017)
Abstract
Background
A wellknown diagnostic imaging modality, termed ultrasound tomography, was quickly developed for the detection of very small tumors whose sizes are smaller than the wavelength of the incident pressure wave without ionizing radiation, compared to the current goldstandard Xray mammography. Based on inverse scattering technique, ultrasound tomography uses some material properties such as sound contrast or attenuation to detect small targets. The Distorted Born Iterative Method (DBIM) based on firstorder Born approximation is an efficient diffraction tomography approach. One of the challenges for a high quality reconstruction is to obtain many measurements from the number of transmitters and receivers. Given the fact that biomedical images are often sparse, the compressed sensing (CS) technique could be therefore effectively applied to ultrasound tomography by reducing the number of transmitters and receivers, while maintaining a high quality of image reconstruction.
Methods
There are currently several work on CS that dispose randomly distributed locations for the measurement system. However, this random configuration is relatively difficult to implement in practice. Instead of it, we should adopt a methodology that helps determine the locations of measurement devices in a deterministic way. For this, we develop the novel DCSDBIM algorithm that is highly applicable in practice. Inspired of the exploitation of the deterministic compressed sensing technique (DCS) introduced by the authors few years ago with the image reconstruction process implemented using l _{1} regularization.
Results
Simulation results of the proposed approach have demonstrated its high performance, with the normalized error approximately 90% reduced, compared to the conventional approach, this new approach can save half of number of measurements and only uses two iterations. Universal image quality index is also evaluated in order to prove the efficiency of the proposed approach.
Conclusions
Numerical simulation results indicate that CS and DCS techniques offer equivalent image reconstruction quality with simpler practical implementation. It would be a very promising approach in practical applications of modern biomedical imaging technology.
Background
Since Wilhelm Roentgen discovered the Xray beam in 1885, there has been a big leap in the clinical diagnostic field in which more advanced technologies for Biomedical imaging applications have been developed, e.g. Magnetic Resonance Imaging (MRI), Computed Tomography (CT), Ultrasound Tomography (UT), Electron Paramagnetic Resonance (EPR), etc.… Among the newly developed techniques, ultrasound imaging techniques have become the widely used tool in the health sector due to its implementation ability for diagnostic and therapeutic as well as its many advantages such as low cost, noninvasive nature, painless test, mobility and fast diagnosis.
Ultrasound imaging which uses sound waves in the range between 20 kHz and 1 GHz is commonly used since the development of sonar in 1910. Based on the principle of sonar, one of the techniques that can be widely used is Bmode imaging [1]. This technique is used for nondestruction evaluation and biomedical imaging. Bmode image represents a qualitative change of acoustic impedance function. As a result of this change, it allows to distinguish between environments in the region of interest. However, this imaging technique, using feedback of sound waves when encountering target, only provides the qualitative information of the imaged targets. Meanwhile, ultrasound tomography, based on inverse scattering technique, provides the quantitative information of those targets. Indeed, when sound waves encounter a heterogeneous environment, some of the energy will then be scattered in all directions. The scattered data will be obtained by the receivers which are set up around the target of interest. Therefore, a set of measurements of the scattered field is obtained. Inverse scattering problem includes estimating the distribution of acoustical parameters (e.g. speed of sound, attenuation and density) to reconstruct the target of interest in the inhomogeneous environment. This technique allows a more detailed description of the imaged target. Instead of using acoustical impedance parameter in Bmode imaging, it uses one of parameters of acoustical properties. Therefore, acoustic tomograms display quantitative information of the target under examination.
Although ultrasound tomography has many advantages, this technique has not been widely applied in practice. One of the reasons is the lack of applications that can take advantage of inverse scattering techniques. Currently, the main application of this technique is only for breast imaging in women to detect cancercausing cells [2,3,4]. Another limitation of inverse scattering techniques is the lack of efficient and powerful calculation methods. Inverse scattering techniques have high computational complexity and it is also the main reason that there is only a certain number of commercialized tomography devices. Hence, stateoftheart inverse scattering techniques primarily focus on reducing the computational complexity and constantly improving the quality of imaging. Most of research work on ultrasound tomography are based on Born approximation. Born Iterative Method (BIM) and Distorted Born Iterative Method (DBIM) are wellknown for diffraction tomography [5]. The DBIM is a quantitative approach in image reconstruction of the very small target. In this method, the background medium is considered inhomogeneous and is updated with each iteration. Therefore, the equation for Green’s function and the equation for incident field are updated with each iteration.
Compressed sensing (CS), which is introduced by Candes and Tao [6] and Donoho [7] in 2006, could acquire and reconstruct sparse signals at a rate lower than that of Nyquist. Random measurement approach in the detection geometry configuration is proposed in [8, 9]. A set of measurements of the scattered field is performed using sets of receiver’s random positions. This method can reduce the computational complexity and improve the quality of the reconstruction of the sound contrast, compared to the linear measurement method. However, this method does not well denoise and the implementation of compressed sampling technique based on random sampling leads to a difficult restriction on the hardware of the ultrasound tomography system. It means that the probes will have to be randomly distributed on the measurement system. This problem is difficult to implement because in order to create a random sequence, one must use a hardware random number generator (HRNG) or a pseudorandom number generator (PRNG). The implementations of HRNG and PRNG in practice are very complex. In [10], the authors proposed to use a deterministic measurement matrix, which is deterministic, instead of random one. The elements of this deterministic measurement matrix are chosen from the sequence generated by a logistic deterministic system. Previous simulated results indicate that this approach offers a very good performance compared to the random method. Moreover, using deterministic CS system inherits a simpler hardware implementation than using the random one. The hardware implementation of the deterministic generator is just a simple nonlinear circuit which is much simpler than a random generator. In this paper, based on the compressed sensing technique, we propose a methodology which helps to determine the locations of transmitters and receivers in a deterministic way for ultrasound tomography. Then, we develop a novel DCSDBIM (Deterministic CSDBIM) algorithm that is highly applicable in practice. As a result, this approach offers a very high performance, compared to the conventional DBIM method.
Methods
Distorted born iterative method
Proposal of measurement configuration
Consider the setting of transmitters and receivers in Fig. 1.
As shown in Fig. 1, the measurement configuration is arranged circularly around the target. The pressure signal from a transmitter is propagated, scattered and measured by receivers. In principle, transmitters and receivers can be distributed uniformly, randomly, or deterministically based on the measurement configuration. Besides, the number of transmitters and receivers depends on different scenarios, in line with practical requirements. It should be noticed that the number of transmitters and receivers should not be so large because the measurement system will be complex. In addition, a large amount of computation and memory requirement are needed to store and process the data. It should be also noticed that the locations of transmitters and receivers may be either different or identical. However, with current ultrasound transducer technology, one transducer can both transmit and receive ultrasound signals. In other words, the location of transmitter and receiver can be identical.
For the proposed measurement configuration, we assume that there are M transmitters and N receivers. M transmitters are arranged around the target to obtain complete information about the target at different angles. Transmitterreceiver procedure of ultrasound wave is performed as follows. Initially, only the first transmitter transmits ultrasound signal and all receivers (N) will simultaneously collect the scattered ultrasonic signal. As a result, we get the 1st set of measurements (i.e. N measurements) at the first position of the transmitter. The second transmitter then begins operating and all receivers simultaneously obtain scattering signal at the second position of the transmitter. At the end of this procedure, we obtain the 2nd set of measurements (i.e. 2 × N measurements). The process continues until the last transmitter (i.e. M^{th}) where we will obtain M sets of measurements (i.e. M × N measurements). With these measured values, full information about the target at different angles around the target are received.
The measured data will then be brought to DBIM to estimate the sound contrast. The change of the sound speed would be utilized to detect any tissue if exists. In this paper, the transmitters and receivers are working under the condition of homogeneous background medium, i.e. water, where there is a target (strange tumor). This target has sound contrast different from the background (denoted as ∆c) which is embedded in this medium. Because DBIM is used to detect strange tumors in early stage in this paper, sound constrast compared to the background medium is usually very small. Therefore, in our scenarios, which will be discussed in section Results below, value of ∆c will be selected to a few percents.
DBIM method
Because of the wellknown properties, Bessel function [11] is usually used in the numerical simulation as a transmitted signal called Incident Wave whose frequency is f. Wavelength λ of this wave is calculated by
where c _{0} is the sound speed in background medium. For homogeneous medium, the received signal at the receiver is the Incident Wave itself. In presence of tumors, the medium becomes inhomogeneous; when the incident wave hits the target, the following two situations may occur: if the target size is much larger than the wavelength of the incident wave, it is reflected; if the target size is less than or equal to the wavelength of the incident wave, it is scattered in all directions around the target. Notice that ultrasound frequency for clinical diagnosis is in the range from 20 KHz to 12 MHz. As a result, the wavelength ranges from 6.2 μm to 74.2 mm if sound speed in the background is 1484 m/s. Born iterative method is used to find out the linear relationship of the scattered pressure difference and the sound contrast difference. The key of this method is that the scattering signal is considered to be very small in comparison with the incident signal. This is entirely consistent with the real requirements which need to detect strange tumors in early stage. Therefore, we consider to reconstruct targets with very small sound contrast, i.e. scattering signal is very small.
The wave equation of the scheme discussed above can be expressed by
where \( {\mathrm{p}}^{\mathrm{sc}}\left(\overrightarrow{r}\right) \), \( {\mathrm{p}}^{\mathrm{inc}}\left(\overrightarrow{r}\right) \), and \( \mathrm{p}\left(\overrightarrow{r}\right) \) are the scattered, incident, and total signals respectively. It can be seen that the known data arethe total signal and the incident signal. However, what we concern about is the reconstruction of the unknown target T (r) from the obtained data. This is the inverse problem.
Consider the wave numbers of the background and target mediums to be k _{0} and k (r) respectively. According to [12], inhomogeneous differential equation has the form
Green function is an effective method to solve inhomogeneous differential equation. Therefore, it is used to find out the nonlinear relationship of the scattered signal and the target based on the total and incident signals. Eq. (1) can be rewritten in details using the Green function G_{0} (·)
Notice that when the background medium is homogenerous, G_{0} is the 0th Hankel function of the first kind and described by
T (r) in Eq. (2) is the target function that needs to be estimated. It can be calculated as follows:
Equation (4) indicates that the ideal target function depends on the frequency of the incident signal (ω = 2πf) and the sound speed difference of the background medium (c _{ 0 }) as well as the target medium (c). In order to calculate in details each pixel inside the region of interest, scattering wave equation (2) needs to be discreterized. The method of moment (MoM) which uses basic sinc function is used to solve this problem.
The total pressure field in the observed mesh area (N × N points) can be expressed by
where \( \overline{C} \) is the Green matrix showing the interactions among pixels, Ī is unit matrix, and D (·) returns a square diagonal matrix of the input vector. The scattered signal in the form of N _{ t } N _{ r } × 1 vector is described by
where \( \overline{B} \) is the Green matrix showing the interaction of all pixels to the receiver. We have to determine two parameters \( \overline{p} \) and \( \overline{T} \) in Eqs. (5) and (6).
By rewritting these equations, we have [13]
where \( \overline{M}=\overline{B}\cdot D\left(\overline{p}\right) \). For a transmitter and a receiver, we formulate a matrix \( \overline{M} \) and a scalar value Δp ^{sc}. The target function \( \overline{T} \) has N ^{2} variables corresponding to the number of pixels in the region of interest. It can be estimated by:
where n and n1 are two consecutive discretetime points. \( \varDelta \overline{T} \) is estimated by using Tikhonov’s regularization [14]:
where \( \varDelta {\overline{p}}^{sc} \) is the difference between estimated and measured scattered signals whose size is (N _{ t } N _{ r } × 1). Besides, measurement results are assembled in a matrix form \( {\overline{M}}_t \) of (N _{ t } N _{ r } × N ^{2}) elements and γ is the regularization factor that needs to be carefully chosen.
The DBIM procedure is presented in Algorithm 1.
Deterministic compressive sampling
Fundamentals of compressive sampling
Compressive Sampling (CS) or known as Compressed Sensing [15] allows for exactly recovering signal v ∈ ℝ^{n} from a small number of “random measurements” u ∈ ℝ^{m} (m < < n), which are defined as follows:
where Φ is a m × n matrix called sampling basis. The columns of Φ have entries equal to 1 at random positions; the entries at other positions are null. This structure means that measurements are done randomly.
The core problem of compressive sampling is that assuming v has sparse representation in an orthonormal basis Ψ, i.e.
in which, w, also known as sparse, only has s < m ≪ n nonzero coefficients. Note that s is the sparsity degree of the signal u and the number m of required measurements is normally higher than s. Compressive sampling theory shows that this sparse property allows accurate recovery w with overwhelming probability success [16]. In particular, sensing basis must have incoherent property to the model basis Ψ [17]. This property is guaranteed by the randomness of the nonzero components in Φ. Therefore, the problem can be written as follows:
where A ∈ R ^{mxn} is fullrank, i.e. the m rows of A are linearly independent. By these settings, Eq. (12) is solved for w with wsparse constraint. Once w is obtained, v can be calculated from (11). With A satisfying “Restricted Isometry Property (RIP)” condition which is first established by Candes et al. [16], the CS problem can be solved through l _{0}minimization technique:
where l _{0} norm is defined as \( {\left\Vert w\right\Vert}_{l_0}:={ \max}_i\left{w}_i\right \). In this paper, instead of solving (13), the sparsest solution \( \widehat{v} \) of (12) can be found by solving the l _{1}minimization problem [18]:
in which, the l _{1} norm is defined as \( \left\Vert w\right\Vert {l}_1:={\displaystyle {\sum}_{i=1}^n\left{w}_i\right.} \)
Equation (12) is introduced under the assumption that the exact form of the reconstructed signal is given. This is rarely the case in practice, because the measurements are often affected by noise. To reconstruct the signal in case of noisy measurements, we have:
where e represents the noise \( {\left\Vert e\right\Vert}_{l_2}\le \varepsilon \). l _{1} problem in presence of noise can be expressed as follows [18]:
Compressive sampling using deterministic filters
In this paper, we use a deterministic basis generated by a pseudo random sequence instead of pure random basis. The advantage of this method in comparison with CS is that it offers a simpler hardware design, because of the deterministic nature of the proposed measurement system. In fact, to facilitate the practical implementation, researchers focused on the measurement system design whose properties are not completely random or completely deterministic. The first work exploited this view is introduced in [11], followed by others [19, 20]. In the work [21], the authors offer some simple criteria to design a deterministic compressed sampling system that ensures the ability to successfully reconstruct sparse signals. Overall, the deterministic compressed sampling technique inherits some advantages over random one, such as more effective recovery time, clear structure, efficient storage and tighter recovery limit [20].
Because of simpler hardware implementation of a deterministic system compared to a random one’s, in this work, the matrix Φ is constructed using the deterministic method. To perform this, we consider a dynamic deterministic system with deterministic nonlinear characteristic; as its dynamic is very sensitive to initial conditions, its output has randomlike property [22]. Note that the reconstruction accuracy of deterministic compressed sampling technique is ensured like the random one [23]. To construct Φ, we consider a Logistic map based dynamic structure that generates a very simple deterministic sequence which is transformed into a sequence that would have a Gaussianlike behavior [22, 23]. The following Logistic map sequence is used as follows
We then convert it by the Logit Transform to become Gaussianlike as
Remarks: When the control parameter ρ in (17) equals to 4, the dynamics described by (17) is deterministic; this means that the dynamic of (17) is very sensitive in regard of the initial condition q (0). A very small change in q (0) would largely change q (n) in a relatively short time. More detailed construction of Φ from q _{ G } (n) can be found in [22]. Finally, the image reconstruction is performed using sparse approximation algorithms.
The deterministic compressive sampling DBIM
Measurement configuration
Assume that there are N _{ t } × N _{ r } transducers where N _{ t } is the number of transmitters and N _{ r } is the number of receivers (that gives N _{ t } × N _{ r } measurements). We need to reconstruct the target which is divided into N pixels vertically and horizontally, i.e. N ^{2} variables. There are three cases for the arrangement of transducers and transmitters that are required for consideration.
Firstly, if N _{ t } × N _{ r } is larger than N ^{2}, the number of equations is greater than the number of variables. This is the overdetermined problem and it is called a dense scattering domain. In this configuration, the amount of scattered information obtained around the target is large, from which we can calculate and easily reconstruct the target image. However, a large number of transmitters and receivers is required. This high number of devices makes the imaging system much more complex in terms of setting up the measurement equipment and of large required processing time consuming in target image reconstruction. Therefore, this configuration is less likely used in practice.
Secondly, if N _{ t } × N _{ r } is approximately or equals to N ^{2}, the number of equations equals to the number of variables. This case is called a moderate scattering domain. For the highresolution images in this paper, i.e. large N ^{2}, we still need in this situation a relatively large amount of transducers to give enough measurement for correctly reconstructing images.
Thirdly, if N _{ t } × N _{ r } is smaller than N ^{2}, i.e. the number of equations is smaller than the number of variables, this is the underdetermined problem. When N _{ t } × N _{ r } is significantly smaller than N ^{2}, we call this a sparse scattering domain, that is, transducers is sparse on the measurement system due to the small number of transducers and the amount of scattered information is less and sparser than the case of the dense and moderate scattering domains. Define the compressed ratio r as the ratio of N _{ t } × N _{ r } and N ^{2}, i.e.
It is clear that for the sparse scattering problem, r < 1. When r is small, the measurement system is simple and the collected data is of small size. Consequently, the system setup and calculation are less complex. Henceforth, we need to focus on the case of the smallest value of r, with which, the imaging system is the simplest and it ensures the acceptable image reconstruction.
The implementation processes of the conventional method and the proposed method are shown in Fig. 2 and Fig. 3 respectively. In both figures, the Input represents the ideal target function and the Output represents the reconstructed target function. Nevertheless, in Fig. 2, the ideal target function of interest is reconstructed using Tikhonov regularizations while it is reconstructed using sparsebased approximation techniques such as l _{1} optimization or greedy Algorithms in Fig. 3. Besides, the measurement configuration of linear transmitteranddetector locations is used in conventional method and the measurement configuration of transmitters and detectors with linear transmitter locations and deterministic detector locations is used in the proposed method.
The signal accurate reconstruction of the CS is guaranteed if the sparse domain is not associated with the sampling domain [6, 7]. In [24], the authors compared the efficiency of sparse signal recovery methods, namely l _{1}LSP, MOSEK, PDCOCHOL, PDCOLSQR, l _{1}MAGIC and HOMOTOPY, and MOSEK. Their results show that l _{1}LSP method outperforms others by its fast reconstruction and low computational complexity. Moreover, this also works effectively in large dense problems. With these observations, the DBIM method can be employed for not only detecting small scale targets (like breast cancer) but also for other largerscale applications. Therefore, in this paper, we use l _{1}regularized leastsquares programs (l _{1}LSP) [24] for reconstructing the image target.
l_{1}regularized leastsquares and DCSDBIM procedure
The l _{1}LSP solves l _{1}regularized least squares problems using the truncated Newton interiorpoint method and has the form as follows:
where A is a data matrix in dense or sparse format with n columns and m rows; w is a vector of length n; and u is a vector of length m. In the DCSDBIM approach, A, w, u are \( {\overline{M}}_t \), \( \varDelta \overline{T} \), and \( \varDelta {\overline{p}}_t^{sc} \) respectively. It is noted that the size of \( {\overline{M}}_t \) is dependent on the number of transmitters (N _{ t }), receivers (N _{ r }) and variables (N ^{2}). In other words, it is changed in different scenarios. In the DCSDBIM form, l _{1} problem is expressed as
where ζ is the regularization parameter in l _{1} problem. The selection of regularization parameter ζ is crucial because it has a great influence on the image reconstruction quality. The large value of ζ will make the reconstructed image rough and the small value of ζ will result in high computational complexity. It is noted that the inverse solver matrix \( {\overline{M}}_t \) is changed at each iteration. In fact, in DBIM, the Green function is updated in each iteration; therefore, matrices B and C are changed, leading to the variation of \( {\overline{M}}_t \) at each iteration.
In this context, the regularization parameter ζ also changes at each iteration. It is chosen as a function of the forward error. For the computation process, the Rayleigh quotient is recursively iterated; the value σ_{0} of the first singular value of the inverse solver matrix \( \overline{M_t} \) is estimated; ζ is then selected as given in [25]. In the simulation scenario presented in this paper, the value of ζ for the first iteration is chosen as 10^{−15}. The DCSDBIM procedure is presented in Algorithm 2.
Results
Parameters set up
The frequency of the incident signal, f = 1 MHz, is selected as in a previous experimental work [26]. It is known that the convergence rate of DCSDBIM method depends on the tolerated reconstruction error. For large error, the convergence is fast; it slows down when the error is fixed smaller. In any case, after a few iterations, the normalized error reaches a floor. Reducing this floor value requires the reduction in the tolerated distortion of reconstructed images in the l _{1} algorithm. Furthermore, it will lead to a more complex computational procedure and as consequence, a much longer imaging time.
The propagation speed of ultrasound wave in the women breast environment is actually in the range of 1350 m/s to 1600 m/s (in the background medium of about 1484 m/s) [27], that is, the difference in the propagation speed in women breast is in the range from 0 to 15.6%. Therefore, choosing the wave propagation speed difference of 5% is reasonable.
The main limitation of the DBIM approach is that the divergence issue occurs when ∆φ > π, where \( \varDelta \varphi =2\omega \left(\frac{1}{c}\frac{1}{c_0}\right) R \) [28]. Therefore, the incident frequency must satisfy \( f<\frac{c_0}{2 d \times \%\varDelta c} \). The incident pressure for a Bessel beam of zero order in twodimensional case is described by
where J _{0} is the 0th order Bessel function and r − r _{ k } is the distance between the transmitter and the k ^{th} point in the range of interest.
Based on the above discussions, for the simulations, the following parameters are finally chosen: frequency f = 1 MHz; the total number of iterations N _{ sum } = 8; N = 21 (i.e. Number of variables is N ^{2} = 21 × 21 = 441); Target diameter = 7.3 mm; Sound contrast 5%; Signaltonoise ratio SNR = 20 dB; Distances from transmitters and receivers to the center of the target are 100 mm. The numerical simulation program used is MATLAB running on a PC with Intel core i3 processor and 2 GB RAM.
Performance evaluation of the DCSDBIM and DBIM
The ideal target function T (r) described in Eq. (4) is shown in Fig. 4 where it is placed at the center of the meshing area with X and Y indicating the pixels coordinates. The conventional linear configuration of transmitters and for the case of N _{ t } = N _{ r } = 22 is presented in Fig. 5a while histogram of linear detector locations over full circle in case of N _{ r } = 22 is shown in Fig. 5b.
We propose a deterministic undersampling configuration of detectors, in which the number of detectors is smaller than the number of detectors in the conventional configuration in Fig. 6. With a reduced number M of measurements and the less computational complexity in the iteration process, the proposed configuration sustains a quality of the reconstruction comparable to that obtained by the conventional configuration. Notice that the transmitters are still placed at equal distance as in the conventional configuration. The proposed configuration of transmitters and detectors with linear transmitter locations and detector DCS locations in case of N _{ t } = N _{ r } = 16 is depicted in Fig. 6a. Meanwhile, Fig. 6b shows the histogram of detector DCS locations over full circle in case of N _{ r } = 16.
To quantify the efficiency of the proposed approach, we acquire the target functions to obtain experimental data to be used in the iterative reconstruction of target image. Then, the error in the reconstructed image is determined and compared to the original image at each iteration. Suppose that m is a P × Q original image (i.e. ideal target function) and \( \widehat{m} \) is the reconstructed image. The normalized absolute error can be defined as:
The normalized absolute errors and runtimes of the DBIM and DCSDBIM methods through iterations with different number of transmitters (N _{ t }) and receivers (N _{ r }) are presented in Table 1. These simulation results have shown that after N _{ sum } iterations, in some cases, the runtime of the DCSDBIM method is significantly larger than that of the DBIM method. The price to pay for implementing easiness and the less complex hardware of the new proposed DCSDBIM approach is the time penalty compared to the classical technique. However, for the ultrasound imaging, the required batch processing would last few minutes so that this time penalty for DCSDBIM does not affect in any way its usefulness in practice.
For the normalized error after N _{ sum } iterations, the image reconstruction quality of the DCSDBIM method is worse than the DBIM method when r <0.5 and is significantly better than the one of the conventional method when r >0.5. In case of r <0.5, although the reconstruction quality of the proposed method is not as good as the conventional method, it can successfully reconstruct the target function when r is very small (in case of r = 0.082 and 0.145). Meanwhile, the conventional method cannot reconstruct the target function (i.e. NaN in Table 1). In case of r >0.5, the reconstruction quality of the proposed method is significantly better than conventional method. We are concerned in practice the case that offers the best performance with a small number of measurements. Thus, we are interested in the case of r = 0.735 (i.e. 324 measurements) that offers a much better performance compared to the conventional method, as shown in Fig. 7.
In general, the simulation results have demonstrated that the DCSDBIM method is a very robust tool with a very highquality reconstruction. It is really a very promising approach in practical applications of modern biomedical imaging technology.
The error performances of both the DCSDBIM method (in case of N _{ t } = N _{ r } = 16, i.e. the number of measurements = 16 × 16 = 256) and the conventional DBIM one (in case of N _{ t } = N _{ r } = 22, i.e. number of measurements = 22 × 22 = 484) are shown in Fig. 8a. Although the number of measurements of the DCSDBIM method is approximately half the one of the DBIM method, these methods offer the same image reconstruction quality after the sixth iteration step. With the same normalized error, DCSDBIM method only needs 3 iterations to complete while the DBIM method requires 6 iterations. Therefore, in this scenario, when using the proposed method, we save half of number of measurements and iterations. Furthermore, it is also shown in Fig. 8b that the proposed method still offers a very good performance (with 400 measurements), compared to the conventional method (with 900 measurements). However, as mentioned above, the runtime for the proposed method is much longer.
The reconstructed results of the DBIM and DCSDBIM approaches through iterations 1 to 8 in case of N _{ t } = N _{ r } = 16 (i.e. r = 0.581) are depicted in Figs. 9 and 10. It can be seen that the convergence rate is obtained very fast by the DCSDBIM approach just after the first few iterations and is not largly affected by noise. In contrast, the DBIM approach has a very low convergence rate and it is much affected by noise. With this observation, the proposed approach offers a much better quality than the conventional approach.
The normalized error and similar performance indexes such as Mean Normalized Absolute Error (MNAE), Root Mean Squared Error (RMSE) are used to evaluate the reconstruction of ultrasound tomography in some relative work [9, 29,30,31,32]. In this paper, we add one more performance index (Qindex) to compare our work (DCSDBIM) with DBIM. Universal image quality index (or Qindex) [33] comparison of the DBIM and DCSDBIM in accordance with different number of measurements is shown in Fig. 11. The obtained results have shown that the Qindex of the DCSDBIM method is much better in the range of small number of measurements. As shown in Fig. 11, a small number of measurements (about 200 measurements) in the DCSDBIM method can also provide the Qindex which is equivalent to a large number of measurements (about 800 measurements) in the DBIM method. Besides, when the number of measurements is small and moderate (less than or equal to N^{2}), as shown in Fig. 8, we found that the quality of the DCSDBIM method is much better than that of the DBIM method.
Reconstruction performances with different values of the sound contrast of the DBIM and DCSDBIM methods are shown in Fig. 12. It can be seen that the DCSDBIM method is quite sensitive to the level of sound contrast. As shown in Fig. 12, when the sound contrast increases, the normalized error also increases. The normalized error has the smallest value when the sound contrast is about 1%. Meanwhile, the DBIM method is relatively stable when changing sound contrast (as the sound contrast increases, the normalized error increases slightly, and almost does not change much). Because the purpose of ultrasound tomography is to detect tumors in early stage, sound contrast is very small, about a few percent. Therefore, in practice, we only consider the small sound contrast. In summary, in case of small sound contrast, the DCSDBIM method is far more effective than the DBIM method.
Performance evaluation of the DCSDBIM and CSDBIM
The performances of DCSDBIM and DBIM have been discussed and compared in the last part. We are going to compare DCSDBIM and CSDBIM performances. In Fig. 13, we have seen that the DCSDBIM offers a performance almost similar to that of the CSDBIM.
The numerical simulation results are quite consistent with previous work on deterministic compressed sampling technique [10]. Furthermore, because the implementation of compressed sampling technique is based on random sampling, this means that transducers are randomly distributed on measurement system. This random structure leads to a complex hardware implementation. If we apply the deterministic compressed sampling technique which uses a nonlinear deterministic system that apparently acts like a random system [21], its hardware implementation of the deterministic structure is simpler than the random one. Note that the exact reconstruction of the deterministic compressed sampling technique is guaranteed like the random one.
The implementation procedure of the DCSDBIM is presented in Fig. 14. Assume that number of pixels of the ideal object function (N) and total number of iterations (N _{ iter }) are given. The flowchart presented in Fig. 14 starts with the initialization of three parameters \( {\overline{O}}_n,{\overline{p}}_0 \), and n (\( {\overline{O}}_n={\overline{O}}_0;{\overline{p}}_0={\overline{p}}^{inc}; n= 0 \)). In addition, N _{ t } N _{ r } is also chosen so that it is smaller than N ^{2}, pixels of the ideal object function (for the best case, 0.5 < r <1). We set up the deterministic measurement configuration based on deterministic filter. Next, the DCSDBIM algorithm is used to reconstruct the target in N _{ sum } iterations. The output result of the procedure is the reconstructed function after N _{ iter } iterations.
Discussion
Based on inverse scattering theory, the DBIM is a wellknown quantitative imaging approach for detecting very small targets thanks to its mechanical properties. Deterministic compressive sampling technique is a promising approach for feasible hardware implementation in practical applications. This paper has successfully applied DCS technique for setting up the measurement configuration for the DBIM, and then the target is reconstructed using l _{1} least square problem in order to improve the quality of the image reconstruction. This method also offers a simpler setting compared to the others. Simulation scenarios of sound contrast reconstruction were implemented to demonstrate the very good performance of this method. These results have shown that the gain of this new approach merits practical considerations. For practical purposes, the reconstruction of threedimensional (3D) image is done by using many 2D slices at different positions of zaxis; indivual processing outcomes are finally all merged together [34]. Therefore, the core issue is that we need to acquire good 2D slices.
With current transducer array technology, one transducer can both transmit and receive ultrasound signal. Thus, when setting up the actual measurement configuration, depending on imaging quality requirement, we can arrange transducers on the measurement system so that the distance between two transducers can be 1°, 2°, etc.… If the distance between two transducers is small, we can arrange multiple transmitters and receivers on the measurement system such that we can reconstruct highresolution images (i.e.large number of pixels in the range of interest); reciprocally, if the distance between two transducers is large, the number of transmitters and receivers will be less. Therefore, we can reconstruct lowresolution images. The number of transmitters and receiverswill have to be chosen in the acceptable range in order to reconstruct goodenough image, i.e. 0.5 < r < 1. However, to be more reasonable, we should arrange the configuration such that the distance between two transducers are small, about 1°. With this arrangement, when we createthe deterministic sequence of the DCS, the indexes of this sequence correspond to the positions of transducers on the measurement system. This creates a randomlike system, and thus ensures conditions of reconstruction in compressed sampling technique [6, 7]. This setup does not make the imaging process more complex. In fact, not all transducers in the measurement system work, only transducers whose indexes coincide with the ones of the deterministic sequence. Therefore, the volume of calculation only depends on the number of active transducers on the measurement system.
Conclusion
Inspired of easier hardware implementation of deterministic CS, in this paper, we have proposed the deterministic measurements in the detection geometry configuration and the image reconstruction process has been implemented using l _{1} regularization. The simulation results of the proposed method have demonstrated the high performance of the proposed approach, where the normalized error is approximately 90% reduced, compared to the conventional DBIM approach. With the same quality, we can save half of number of measurements and only use two iterations when using the proposed method. Furthermore, numerical simulation results also indicate that CS and DCS techniques offer equivalent image reconstruction quality. However, DCS has the advantage of being able to execute much more convenient hardware. With current transducer technology and using the slicing technique that transforms 3D images to several 2D images, we can produce multiple sliced images at the same time. This will significantly reduce imaging time for patients and make ultrasound tomography imaging a realtime imaging tool. To implement this approach in practice while keeping the cost at an acceptable level, the following two possibilities might be considered: 1) We only need a 2D measurement system; when imaging different slices, the measurement system will shift along the z axis to produce images of different slices. However, this will take more time and the shifted measurement system will cause mechanical errors; 2) We set up the measurement system following women’s breast shape; the process of imaging at different slices will simultaneously take place. However, can we arrange transducers on various slices differently and sparsely? If this can be done, setup costs would significantly reduce. This is a fascinating issue for further study.
Abbreviations
 BIM:

Born iterative method
 CS:

Compressed sensing
 CT:

Computed tomography
 DBIM:

Distorted born iterative method
 DCS:

Deterministic compressed sensing
 DCSDBIM:

Deterministic compressed sensing distorted born iterative method
 EPR:

Electron paramagnetic resonance
 HRNG:

Hardware random number generator
 MRI:

Magnetic resonance imaging
 PRNG:

Pseudorandom number generator
 RIP:

Restricted Isometry property
 UT:

Ultrasound tomography
References
Schueler CF, Lee H, Wade G. Fundamentals of digital ultrasonic processing. IEEE Trans Sonics Ultrason. 1984;31:195–217.
Greenleaf J, Ylitalo J, Gisvold J. Ultrasonic computed tomography for breast examination. IEEE Eng Med Biol Mag. 1987;6:27–32.
Andre MP, Janée HS, Martin JP, et al. Highspeed data acquisition in a diffraction tomography system employing largescale toroidal arrays. Int J Imag Syst Technol. 1997;8(1):137–47.
Wiskin J, Borup DT, Johnson SA, et al. Fullwave, nonlinear, inverse scattering. Acoust Imag. 2007;28:183–94.
Devaney AJ. Inversion formula for inverse scattering within the Born approximation. Opt Lett. 1982;7:111–2.
Candès EJ, Romberg J, Tao T. Robust uncertainty principles: Exact signal reconstruction from highly incomplete frequency information. IEEE Trans Inf Theory. 2006;52(2):489–509.
Donoho DL. Compressed sensing. IEEE Trans Inf Theory. 2006;52(4):1289–306.
TranDuc T, NguyenLinh T, DoNgoc M. Modified Distorted Born Iterative Method for Ultrasound Tomography by Random Sampling. Australia: International Symposium on Communications and Information Technologies (ISCIT); 2012.
Sloun R, Pandharipande A, Mischi M, et al. Compressed Sensing for Ultrasound Computed Tomography. IEEE Trans Biomed Eng. 2015;62(6):1660–4.
NguyenLinh T, Dinh VP, Hussain MZ, et al. Compressed sensing using chaos filters. In: Australian Telecommunication Networks and Applications Conference (ATNAC2008). Australia: National Wine Centre Adelaide; 2008.
Krainov VP, Reiss RH, Smirnov MB. Appendix J: Properties of the Generalized Bessel Function. In: Radiative Processes in Atomic Physics. New York: John Wiley & Sons; 2005.
Avinash CK, Slaney M. Principles of Computerized Tomographic Imaging. New York: IEEE Press Society for Industrial and Applied Mathematics; 2001.
Lavarello RJ, Oelze M. Tomographic Reconstruction of ThreeDimensional Volumes Using the Distorted Born Iterative Method. IEEE Trans Med Imaging. 2009;28(1):1643–53.
Golub GH, Hansen PC, O’Leary DP. Tikhonov Regularization and Total Least Squares. SIAM J Matrix Anal Appl. 1999;21(1):185–94.
Candes EJ, Romberg JK, Tao T. Stable signal recovery from incomplete and inaccurate measurements. Comm Pure Appl Math. 2006;59(8):1207–23.
Candes EJ, Wakin MB. An Introduction To Compressive Sampling. IEEE Signal Process Mag. 2008;25:21–30.
Candès EJ, Romberg J. Sparsity and incoherence in compressive sampling. Inverse Problems. 2007;23:969.
Candès EJ. The restricted isometry property and its implications for compressed sensing. C R Acad Sci. 2008;346:589–92.
Duarte MF, Eldar YC. Structured compressed sensing: From theory to applications. IEEE Trans Signal Process. 2011;59(9):4053–85.
Calderbank R, Howard S, Jafarpour S. Construction of a large class of deterministic sensing matrices that satisfy a statistical isometry property. IEEE J Sel Top Sign Proces. 2010;4(2):358–74.
Yu L, Barbot JP, Zheng G, et al. Compressive sensing with chaotic sequence. IEEE Signal Process Lett. 2010;17(8):731–4.
Sprott JC. Chaos and timeseries analysis. Oxford: Oxford University Press; 2003.
Tropp J, Wakin BM, Duarte FM, et al. Random filters for compressive sampling and reconstruction. Toulouse: IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP 2006); 2006.
SeungJean K, Koh K, Lustig M, et al. An interiorpoint method for largescale l1regularized least squares. IEEE J Sel Top Sign Proces. 2007;1(4):606–17.
Lavarello R, Oelze M. A study on the reconstruction of moderate contrast targets using the distorted Born iterative method. IEEE Trans Ultrason Ferroelectr Freq Control. 2008;55(1):112–24.
Lavarello RJ. New Developments on Quantitative Imaging Using Ultrasonic Waves. Ph.D. Dissertation, University of Illinois at UrbanaChampaign; 2009. www.brl.uiuc.edu/Downloads/lavarello/Lavarello_PhD.pdf.
Greenleaf JF, Bahn RC. Clinical imaging transmissive ultrasonic computerized tomography. IEEE Trans Biomed Eng. 1981;28:177–85.
Slaney M, Kak AC, Larson LE. Limitations of imaging with first order diffraction tomography. IEEE Trans Microwave Theory Tech. 1984;32(8):860–73.
Lavarello RJ, Oelze M. Density imaging using inverse scattering. J Acoust Soc Am. 2009;125(2):793–802.
Lavarello RJ, Bond SD, et al. Regularized tomographic density imaging using multiple frequency information. Ultrasonics Symposium (IUS). 2010.
Haddadin OS, Ebbini ES. Imaging strongly scattering media using a multiple frequency distorted Born iterative method. IEEE Trans Ultrason Ferroelectr Freq Control. 1998;45(6):1485–96.
Zaiping N, Yerong Z. Hybrid Born iterative method in lowfrequency inverse scattering problem. IEEE Trans Geosci Remote Sens. 1998;36(3):749–53.
Wang Z, Student Member, IEEE, Bovik AC, Fellow, IEEE. A Universal Image Quality Index. IEEE Signal Process Lett. 2002;9(3):81–4.
Duric N, Littrup P, Roy O, et al. Clinical breast imaging with ultrasound tomography: A description of the SoftVue system. J Acoust Soc Am. 2014;135(4):2155.
Acknowledgements
Not applicable.
Funding
This work was supported by Ministry of Education and Training under Grant No. B2017SP208.
Availability of data and materials
All data generated or analyzed during this study are included in this published article.
Authors’ contributions
The last author TTD proposed the idea of using CS for ultrasound tomography and with the collaboration of the second author HHT, elaborated a detailed working plan. The first author TQH, with the help of the last author, executed the research described in the detailed working plan and wrote a first draft. The third author TTL with the collaboration of the second author reworked the draft several times; the second and the last author critically revised the final version in order to satisfy all authors personal point of view. All authors read and approved the final manuscript.
Competing interests
The authors declare that they have no competing interests.
Consent for publication
Not applicable.
Ethics approval and consent to participate
Not applicable.
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Author information
Authors and Affiliations
Corresponding author
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Cite this article
Huy, T.Q., Tue, H.H., Long, T.T. et al. Deterministic compressive sampling for highquality image reconstruction of ultrasound tomography. BMC Med Imaging 17, 34 (2017). https://doi.org/10.1186/s1288001702068
Received:
Accepted:
Published:
DOI: https://doi.org/10.1186/s1288001702068
Keywords
 Mammography
 Ultrasound tomography
 Inverse scattering
 Distorted born iterative method (DBIM)
 Deterministic compressive sampling (DCS)