# A Novel Iterative Thresholding Algorithm Based on Plug-and-Play Priors for Compressive Sampling
## 2 Compressive Sampling and Fast Iterative Shrinkage-Thresholding Algorithm  
The compressive sampling recovery is an underdetermined problem. More than one solution can yield the same measurements, but if \(x\) is sparse, we can recover the original signal by minimizing the following L1 problem. \[\hat{x}=\underset{x}{\text{argmin}}\frac{1}{2}\|y-\mathbf{A}x\|_{2}^{2}+ \lambda\|x\|_{1}, \tag{1}\]  
where \(\lambda\) is a regularization parameter. The above problem can be efficiently solved by various methods [22; 23; 24; 25], one of which is ISTA (iterative shrinkage-thresholding algorithm) [26]. Specifically, the general step of ISTA is summarized in Algorithm 1, where \(f=\|y-\mathbf{A}x\|_{2}^{2}/2\), and \(\
abla f(x^{k})=\mathbf{A}^{*}(\mathbf{A}x^{k}-y)\) denotes the gradient of the function \(f\) at the point \(x^{k}\); \(\mathbf{A}^{*}\) denotes the conjugate transpose of \(\mathbf{A}\); \(c\) is a step size. The undersampled Fourier matrix \(\mathbf{A}=\mathbf{S}\mathbf{F}\) is used as the measurement matrix, where \(\mathbf{F}\) is the 2D Fourier transform and \(\mathbf{S}\) is a selection matrix containing \(m\) rows of the identity matrix. Step (b) is a proximity operator. It has a close form solution, which can be expressed as:  
\[x^{k}=soft(x_{g},\lambda c)=\text{sgn}(x_{g})\cdot\text{max}(abs(x_{g})- \lambda c,0). \tag{2}\]  
FISTA [21] offers an accelerated version of ISTA by adding an accelerated step (step (c)), which is summarized in Algorithm 2. ```
Input: the CS measurements \(y\) and the measurement matrix \(\mathbf{A}\) Initialization:\(x^{0}=0\), \(\lambda\) for\(k=1\) to \(K\)do (a) \(x_{g}=x^{k}-c\
abla f(x^{k})\) (b) \(x^{k}=\underset{x}{\text{argmin}}\{\lambda\|x\|_{1}+\|x-x_{g}\|_{2}^{2}/(2c)\}\) endfor
```  
**Algorithm 1**ISTA [26]```
Input: the CS measurements \(y\) and the measurement matrix \(\mathbf{A}\) Initialization:\(x^{0}=r^{1}=0\), \(\lambda\) for\(k=1\) to \(K\)do (a) \(x_{g}=r^{k}-c\
abla f(r^{k})\) (b) \(x^{k}=\underset{x}{\operatorname{argmin}}\{\lambda\|x\|_{1}+\|x-x_{g}\|_{2}^{2 }/(2c)\}\) (c) \(t^{k+1}=\frac{(1+\sqrt{1+4(t^{k})^{2}})}{2}\); \(r^{k+1}=x^{k}+(\frac{(t^{k}-1)}{\mu^{k+1}})(x^{k}-x^{k-1})\) endfor
```

**Algorithm 2** FISTA [21]

## 3 CS via Composite Regularization and Adaptive Thresholding

### The New Composite Model

Images are generally sparse in wavelet and gradient domain, although not strictly. Moreover, there exist abundant similar image patches in natural images. These sparse priors can be utilized by using three regularizations (TV, wavelet sparse, and nonlocal wavelet regularization). The nonlocal sparse regularization can preserve image fine details, while TV and wavelet sparse regularization can maintain the consistency of the entire image. In order to gain the advantages of these regularizations, we combine the nonlocal self-similarity of patches and the global sparsity in wavelet and gradient domain to jointly enforce the solution in image CS recovery, and formulate a composite regularization problem as:

\[\hat{x}=\underset{x}{\operatorname{argmin}}\frac{1}{2}\|y-\mathbf{A}x\|_{2}^{ 2}+\alpha\|x\|_{TV}+\beta\|x\|_{Wa\!\!\!\!\!\



Steps (b)-(d) are the so-called proximal map, which is formally defined in [27], i.e., given a continuous convex function \(g(x)\) and any scalar \(\rho>0\), the proximal map associated to function \(g\) is defined as:

\[prox_{\rho}(g):=\operatorname*{argmin}_{u}(u)+\frac{1}{2\rho}\|u-x\|_{2}^{2}, \tag{5}\]

In the plug-and-play prior approach [18; 19], this proximal map is equal to a denoising operation, where \(x\) is a noisy image, \(g(u)\) is a regularization corresponding to a denoiser, and \(\rho\) is the denoising threshold. In the steps (b)-(d), we regard \(x_{g}\) as a noisy image; then the corresponding denoising methods can be adopted to obtain the solution of each sub-problem. For step (b), we use a TV-based denoising method which is based on the FISTA algorithm [27].