4592

A fast approach for estimation of Spark of the sensing matrix for Compressed Sensing applications.
Bhairav Bipin Mehta1, Mingrui Yang2, and Mark Alan Griswold1

1Radiology, Case Western Reserve University, Cleveland, OH, United States, 2Department of Biomedical Engineering, Cleveland Clinic, Cleveland, OH, United States

Synopsis

Compressed sensing (CS) has been extensively used with wide spread application in MRI and other signal processing fields. Spark of the sensing matrix is at the heart of the CS framework for determining the success of the signal recovery for a given designed CS system. However, estimation of Spark of the sensing matrix is a combinatorial process, thus, practically difficult to estimate for realistic sizes of sensing matrices. The purpose of this work is to present a new optimization-problem-based approach for estimation of the Spark of the sensing matrix which will overcome the existing limitations, thereby, a tool to assess and design CS framework based systems.

Purpose & Introduction:

Compressed sensing (CS)[1-3] has been extensively used with wide spread application in MRI[4] and other signal processing fields[5]. CS has been instrumental in the past decade in battling the acquisition speed limit of MRI. In order to use the CS framework for a particular application, one needs to address the following questions: What is the optimal sparsifying transform for their application? What is the best sampling pattern for a given design? What is the highest possible acceleration factor? Is the current system design achieving it? When would a given system fail? etc. For vast majority of CS applications, these questions have been answered heuristically or intuitively. However, these approaches might not suffice for certain scenarios, for instance, usage of the CS framework for clinical application or hardware implementation of CS. One of the primary reasons for the use of heuristic or intuitive approaches is that the existing assessment tools are either impractical to use for realistic systems or they provide a rough assessment of the system, which might be far away from underlying theoretical bounds. For example, Spark or Restricted isometry property (RIP) of the sensing matrix are at the heart of the CS framework for determining the success of the signal recovery for a given designed CS system[2,3,6]. However, estimation of Spark or RIP is a combinatorial process, thus, practically difficult to estimate for realistic sizes of sensing matrices. Thus, in this work we focus on developing a tool to assess the sensing system/matrix. Specifically, the purpose of this work is to present a new optimization-problem-based approach for estimating the Spark of the sensing matrix which will overcome the existing limitations.

Theory:

In general, we are looking to solve the following inverse problem

$$Sx_0=y\quad\quad\quad(1)$$

where $$$S\in\mathbb{C}^{M\times N}$$$ is the sensing matrix with $$$M<N,x_0\in\mathbb{C}^{N}$$$ is the ground truth object in the desired transform domain and $$$y\in\mathbb{C}^{M}$$$ is the sampled data. The solution space for a given $$$S$$$ and $$$y$$$ is

$$\{x |x=x_0+n,n\in null(S)\}\quad\quad\quad(2)$$

Donoho et. al.[3] defined Spark of matrix A, as the smallest possible number of columns from A that are linearly dependent. Using this definition Donoho et. al.[3] laid down the theoretical limit and thereby the sampling theorem of the CS framework. They proved that, for any $$$k$$$-sparse, $$$x_0$$$, to be the sparsest element within the solution space, the following condition should be satisfied[3]

$$\parallel\space n\space\parallel_0\space>\space2k\quad\quad\forall\space n\in\{n|Sn=\bar{0},n\neq\bar{0}\}\quad\quad\quad(3)$$

Thus, Spark is at the heart of the assessment of guaranteed signal recovery for a given system. The brute force approach to estimate Spark requires combinatorial assessment of the columns of the sensing matrix. To the best of our knowledge, currently there is no approach for estimating Spark in practically feasible time for realistic sizes of sensing matrices. Mutual coherence is a surrogate measure and could be far away from the actual value of the Spark[2,3,5]. Since Spark is by definition the $$$\ell_0$$$ norm of the sparsest element of the subspace $$$\{ n|Sn=\bar{0},n\neq \bar{0}\}$$$, we propose to estimate Spark by solving the following optimization problem.

$$\min_n\space\parallel n\parallel_1\quad\quad subject\space to\space Sn=\bar{0} \space and\space\parallel n\parallel_2=1\quad\quad\quad(4)$$

where Spark would be the $$$\ell_0$$$ norm of the solution of the optimization problem in equation (4).

Methods:

For assessment of the approach we performed Monte Carlo (10 repetitions) based simulations, where sensing matrices were generated with sizes ranging from 8 by 16 to 15 by 30. The sensing matrices were generated randomly using a uniform distribution. For each sensing matrix, the value of Spark was estimated through two approaches, namely, the brute-force approach and the proposed optimization-based approach. The steps for brute-force approach are shown in figure 1. For the optimization-based approach, equation 4 was solved using interior-point algorithm. The accuracy and computation time of the proposed approach was evaluated by comparing with the brute-force approach. Finally, the proposed approach was used to compare three sampling patterns and evaluated using CS reconstruction of the Shepp-Logan phantom.

Results & Discussions:

Figure 2 illustrates the agreement between the two approaches for different sizes of sensing matrix. Figure 3 and 4 show the proposed approach substantially reduces the estimation time of Spark. Figure 5 illustrates the estimated Spark value identifies which sampling pattern is more suitable for the given system design. Furthermore, it also accurately predicts the feasibility of the CS reconstruction given the sparsity of the underlying signal.

Conclusion:

The proposed approach provides a fast, accurate and practically feasible means for estimation of Spark for realistic sizes of sensing matrices, thereby, a tool to assess and design CS framework based systems.

Acknowledgements

The authors would like to acknowledge funding from Siemens Healthcare. This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University.

References

1. D. L. Donoho, “Compressed sensing,” IEEE Trans. Info. Theory, vol. 52, no. 4, pp. 1289–1306, Sep. 2006.

2. E. J. Candes, “Compressive sampling,” in Int. Congress of Mathematicians, Madrid, Spain, 2006, vol. 3, pp. 1433–1452.

3. D. L. Donoho and M. Elad, “Optimally sparse representation in general (nonorthogonal) dictionaries via $$$\ell_1$$$ minimization,” Proc.Nat. Acad. Sci., vol. 100, no. 5, pp. 2197–2202, Mar. 2003.

4. M. Lustig, D. L. Donoho, and J. M. Pauly, “Sparse MRI: The application of compressed sensing for rapid MR imaging,” Magn. Res. Med.,vol. 58, no. 6, pp. 1182-1195, Dec. 2008.

5. M. Duarte and Y. Eldar, “Structured compressed sensing: From theoryto applications,” IEEE Trans. Signal Process., vol. 59, no. 9, pp. 4053–4085, Sep. 2011.

6. E. J. Candes, J. K. Romberg, and T. Tao, “Stable signal recovery from incomplete and inaccurate measurements,” Comm. Pure Appl.Math., vol. 59, no. 8, pp. 1207–1223, 2006.

Figures

Figure 1: Steps used in the brute-force approach. Parallel computing power of our institutions cluster was used to reduce the computation time of the brute-force approach.

Figure 2: Spark estimation accuracy of the proposed L1 minimization approach. A: Bar plot displaying the estimated spark value using the two approaches for different matrix sizes. B: Bland-Altman plot comparing the two approaches. Both the plots illustrate the proposed L1 minimization approach accurately estimates the Spark value for different sizes of sensing matrix.

Figure 3: Comparison of computational time between the two approaches. A: Bar plot illustrating the computation time for the brute-force approach. B: Bar plot illustrating the computation time for the proposed L1 minimization approach. We see a significant reduction in computation time using the proposed approach. Please note the brute-force approach used parallel computing power of a high performance cluster while the L1 minimization approach didn’t. Also, the higher standard deviation in times for brute force approach is primarily due to variation in availability of cluster nodes.

Figure 4: Computational speed up factors for different sensing matrix sizes. The speed up factor is measured as a ratio of computation time using brute force approach divided by the corresponding time using L1 minimization approach. As anticipated, the speed up factor increases as the matrix size increases.

Figure 5: Example assessment of sampling pattern and evaluation of the feasibility of signal recovery using Spark. For comparison of the three sampling patterns, mutual coherence[4] (mc) suggests pattern 3 is most suitable while Spark suggests pattern 2 is most suitable. If we try to recover the Shepp-Logan phantom (transform sparsity=741) from the sampled data, then system with pattern 2(row 2) recovers the phantom (since Spark>2k) while system with pattern 3(row 3) fails to recover (since Spark<2k). Interestingly, if we use a modified phantom (transform sparsity=688) then the system with pattern 3(row 4) successfully recovers the data since Spark>2k.

Proc. Intl. Soc. Mag. Reson. Med. 27 (2019)
4592