Non-convex Fraction Function Penalty: Sparse Signals Recovered from Quasi-linear Systems

TL;DR

This paper introduces a non-convex fraction function penalty for quasi-linear compressed sensing, enabling effective sparse signal recovery in nonlinear systems, with an iterative thresholding algorithm outperforming existing methods.

Contribution

The paper proposes a novel non-convex fraction function penalty and an iterative thresholding algorithm tailored for quasi-linear compressed sensing problems.

Findings

The method effectively recovers sparse signals in nonlinear systems.

Numerical experiments show superior performance over state-of-the-art methods.

The algorithm's adaptability via parameter $a$ enhances recovery results.

Abstract

The goal of compressed sensing is to reconstruct a sparse signal under a few linear measurements far less than the dimension of the ambient space of the signal. However, many real-life applications in physics and biomedical sciences carry some strongly nonlinear structures, and the linear model is no longer suitable. Compared with the compressed sensing under the linear circumstance, this nonlinear compressed sensing is much more difficult, in fact also NP-hard, combinatorial problem, because of the discrete and discontinuous nature of the $\ell_{0}$-norm and the nonlinearity. In order to get a convenience for sparse signal recovery, we set most of the nonlinear models have a smooth quasi-linear nature in this paper, and study a non-convex fraction function $\rho_{a}$ in this quasi-linear compressed sensing. We propose an iterative fraction thresholding algorithm to solve the regularization problem $(QP_{a}^{\lambda})$ for all $a>0$. With the change of parameter $a>0$, our algorithm could get a promising result, which is one of the advantages for our algorithm compared with other algorithms. Numerical experiments show that our method performs much better compared with some state-of-art methods.

Equations

$$(P_0)\underset{x\in {\mathcal{R}}^n}{\min }\Vert x{\Vert }_0\mathrm{subject}\mathrm{to}Ax=b$$

$$A(x)=b.$$

$$F:{\Re }^n\to {\Re }^{m\times n}$$

$$A(x)=F(x)x$$

$$(Q{P}_0)\underset{x\in {\Re }^n}{\min }\Vert x{\Vert }_0\mathrm{subject}\mathrm{to}F(x)x=b.$$

$$(Q{P}_1)\underset{x\in {\Re }^n}{\min }\Vert x{\Vert }_1\mathrm{subject}\mathrm{to}F(x)x=b.$$

$$(QP_{1}^{\lambda})\ \ \ \ \ \ \min_{x\in\Re^{n}}\Big{\{}\|F(x)x-b\|_{2}^{2}+\lambda\|x\|_{1}\Big{\}}$$

$$P(x)=P_a(x)=\sum _{i=1}^n{\rho }_a(x_i),a>0$$

$${\rho }_a(t)=\frac{a|t|}{a|t|+1}$$

$$\lim_{a\rightarrow+\infty}\rho_{a}(x_{i})=\left\{\begin{array}[]{ll}0,&{\ \ \mathrm{if}\ x_{i}=0;}\\ 1,&{\ \ \mathrm{if}\ x_{i}\neq 0.}\end{array}\right.$$

$$(Q{P}_a)\underset{x\in {\mathcal{R}}^n}{\min }{P}_a(x)\mathrm{subject}\mathrm{to}F(x)x=b$$

$$(QP_{a}^{\lambda})\ \ \ \ \ \min_{x\in\mathcal{R}^{n}}\Big{\{}\|F(x)x-b\|_{2}^{2}+\lambda P_{a}(x)\Big{\}}.$$

$$f_{\lambda }(\beta )=(\beta -\gamma {)}^2+\lambda \cdot {\rho }_a(\beta )$$

$${\mathrm{prox}}_{a,\lambda }^{\beta }(\gamma ):=\arg \underset{\beta \in \mathcal{R}}{\min }{f}_{\lambda }(\beta ).$$

$$\mathrm{prox}_{a,\lambda}^{\beta}(\gamma)=\left\{\begin{array}[]{ll}g_{a,\lambda}(\gamma),&\ \ \mathrm{if}\ {|\gamma|>t^{\ast};}\\ 0,&\ \ \mathrm{if}\ {|\gamma|\leq t^{\ast}.}\end{array}\right.$$

$$g_{a,\lambda }(\gamma )=sign(\gamma )(\frac{\frac{1+a|\gamma |}{3}(1+2\cos (\frac{\varphi (\gamma )}{3}-\frac{\pi }{3}))-1}{a}),$$

$$\varphi (\gamma )=\arccos (\frac{27\lambda a^2}{4(1+a|\gamma |{)}^3}-1)$$

$$t_{a,\lambda}^{\ast}=\left\{\begin{array}[]{ll}t_{a,\lambda}^{1},&\ \ \mathrm{if}\ {\lambda\leq\frac{1}{a^{2}};}\\ t_{a,\lambda}^{2},&\ \ \mathrm{if}\ {\lambda>\frac{1}{a^{2}}.}\end{array}\right.$$

$$t_{a,\lambda }^1=\frac{\lambda }{2}a,t_{a,\lambda }^2=\sqrt{\lambda }-\frac{1}{2a}.$$

$$G_{\lambda ,P}(x)=({\mathrm{prox}}_{a,\lambda }^{\beta }(x_1),\cdots ,{\mathrm{prox}}_{a,\lambda }^{\beta }(x_n){)}^T$$

$$C_1(x)=\Vert F(x)x-b{\Vert }_2^2+\lambda P_a(x)$$

$$\begin{array}[]{llll}C_{2}(x,y)&=&\mu\|F(y)x-b\|_{2}^{2}+\lambda\mu P_{a}(x)\\ &&-\mu\|F(y)x-F(y)y\|_{2}^{2}+\|x-y\|_{2}^{2}.\end{array}$$

$$\min_{x\in\Re^{n}}\Big{\{}\|x-B_{\mu}(y)\|_{2}^{2}+\lambda\mu P_{a}(x)\Big{\}}$$

$$C_2(x,y)$$

$$\min_{x\in\Re^{n}}\Big{\{}\|x-B_{\mu}(y)\|_{2}^{2}+\lambda\mu P_{a}(x)\Big{\}}.$$

$$C_2(x^{\ast },x^{\ast })\le C_2(x,x^{\ast })$$

$$C_2(x,x^{\ast })$$

$$x^{\ast }=G_{a,\lambda \mu }(B_{\mu }(x^{\ast }))$$

$$x^{\ast}_{i}=\left\{\begin{array}[]{ll}g_{a,\lambda\mu}(|B_{\mu}(x^{\ast})_{i}|),&\ \ \mathrm{if}\ {|B_{\mu}(x^{\ast})_{i}|>t_{a,\lambda\mu}^{\ast};}\\ 0,&\ \ \mathrm{if}\ {|B_{\mu}(x^{\ast})_{i}|\leq t_{a,\lambda\mu}^{\ast}.}\end{array}\right.$$

$$x^{k+1}=G_{\lambda \mu ,P}(B_{\mu }(x^k))$$

Full text

I Introduction

In compressed sensing (see, e.g., [1, 2]), the problem of reconstructing a sparse signal under a few linear measurements which are far fewer than the dimension of the ambient space of the signal can be modeled into the following $\ell_{0}$-minimization:

[TABLE]

where $A\in\Re^{m\times n}$ is a $m\times n$ real matrix of full row rank with $m<n$, and $b\in\Re^{m}$ is a nonzero real vector of $m$-dimension, and $\|x\|_{0}$ is the $\ell_{0}$-norm of the real vector $x$, which counts the number of the non-zero entries in $x$ (see, e.g., [3, 4, 5]). In general, the problem $(P_{0})$ is computational and NP-hard because of the discrete and discontinuous nature of the $\ell_{0}$-norm.

However, many real-life applications in physics and biomedical sciences carry some strongly nonlinear structures [6], so that the linear model in problem $(P_{0})$ is no longer suitable. We consider a map $A:\Re^{n}\rightarrow\Re^{m}$, which is no longer necessarily linear, and reconstruct a sparse vector $x\in\Re^{n}$ from the measurements $b\in\Re^{m}$ given by

[TABLE]

In order to get a convenience for sparse signal recovery, we set most of the nonlinear models have a smooth quasi-linear nature. By this means, there exists a Lipschitz map

[TABLE]

such that

[TABLE]

for all $x\in\Re^{n}$.

The sparse signals recovered under the quasi-linear case can be mathematically viewed as the following form

[TABLE]

Similarly, the quasi-linear compressed sensing is also combinatorial and NP-hard (see, e.g., [6, 7]).

The $\ell_{1}$-norm is the most famous convex relaxation (see, e.g., [6, 7]), and the minimization for quasi-linear compressed sensing has the following form

[TABLE]

for the constrained problem and

[TABLE]

for the regularization problem, where $\|x\|_{1}=\sum_{i=1}^{n}|x_{i}|$ is the $\ell_{1}$-norm of vector $x$.

In problem $(QP_{1})$, many excellent theoretical works (see, e.g., [6, 7]) have shown that the $\ell_{1}$-norm minimization can really make an exact recovery in some specific conditions. In general, however, it may be suboptimal for recovering a sparse signal, and the regularization problem $(QP_{1}^{\lambda})$ leads to a biased estimation by shrinking all the components of the vector toward zero simultaneously, and sometimes results in over-penalization in the regularization model $(QP_{1}^{\lambda})$ as the $\ell_{1}$-norm in linear compressed sensing.

Inspired the good performance of the fraction function in image restoration and linear compressed sensing (see, e.g., [14, 17]), in this paper, we replace $\|x\|_{0}$ with a continuous sparsity promoting penalty function

[TABLE]

where

[TABLE]

is the fraction function, and it is called ”strictly non-interpolating” in [14]. Clearly, $\rho_{a}(t)$ is increasing and concave in $t\in[0,+\infty]$.

With the change of parameter $a>0$, the non-convex function $P_{a}(x)$ interpolates the $\ell_{0}$-norm

[TABLE]

By this transformation, the minimization problem $(QP_{0})$ could be translated into the following minimization problem

[TABLE]

for the constrained problem and

[TABLE]

for the regularization problem.

The rest of this paper is organized as follows. In Section II, we propose an iterative fraction thresholding algorithm to solve the regularization problem $(QP_{a}^{\lambda})$ for all $a>0$. In Section III, we present some numerical experiments to demonstrate the effectiveness of our algorithm. The concluding remarks are presented in Section IV.

II Iterative fraction thresholding algorithm (IFTA)

In this section, the iterative fraction thresholding algorithm (IFTA) is proposed to solve the regularization problem $(QP_{a}^{\lambda})$ for all $a>0$. Before we embark to this discussion, some results need to be expressed before IFTA is proposed to solve the regularization problem $(QP_{a}^{\lambda})$.

We define a function of $\beta\in\mathcal{R}$ as

[TABLE]

and

[TABLE]

Lemma 1**.**

The operator $\mathrm{prox}_{a,\lambda}^{\beta}$ defined in (14) can be expressed as

[TABLE]

where $g_{a,\lambda}(\gamma)$ is defined as

[TABLE]

[TABLE]

and the threshold value satisfies

[TABLE]

where

[TABLE]

The proof of Lemma 1 used the Cartan s root-finding formula expressed in terms of hyperbolic functions and it is a special case of the reference [18], and the detailed proof can be seen in [17].

Definition 1**.**

The iterative thresholding operator $G_{\lambda,P}$ can be defined by

[TABLE]

where $\mathrm{prox}_{a,\lambda}^{\beta}$ is defined in Lemma 1.

Nextly, we will show that the optimal solution to $(QP_{a}^{\lambda})$ could be expressed as a operation.

For any fixed positive parameters $\lambda>0$, $\mu>0$, $a>0$ and $x\in\Re^{n}$, let

[TABLE]

and

[TABLE]

Clearly, $C_{2}(x,x)=\mu C_{1}(x)$.

Theorem 1**.**

For any fixed positive parameters $\lambda>0$, $\mu>0$ and $y\in\Re^{n}$, $\displaystyle\min_{x\in\Re^{n}}C_{2}(x,y)$ is equivalent to

[TABLE]

where $B_{\mu}(y)=y+\mu F(y)^{\ast}(b-F(y)y)$.

Proof.

By the definition, $C_{2}(x,y)$ can be rewritten as

[TABLE]

which implies that $\displaystyle\min_{x\in\Re^{n}}C_{2}(x,y)$ for any fixed positive parameters $\lambda>0$, $\mu>0$ and $y\in\mathcal{R}^{n}$ is equivalent to

[TABLE]

∎

Theorem 2**.**

For any fixed positive parameter $\lambda>0$ and $0<\mu<L_{\ast}^{-1}$ with $\|F(x^{\ast})x-F(x^{\ast})x^{\ast}\|_{2}^{2}\leq L_{\ast}\|x-x^{\ast}\|_{2}^{2}$. If $x^{\ast}$ is the optimal solution of $\displaystyle\min_{x\in\Re^{n}}C_{1}(x)$, then $x^{\ast}$ is also the optimal solution of $\displaystyle\min_{x\in\Re^{n}}C_{2}(x,x^{\ast})$, that is

[TABLE]

for any $x\in\Re^{n}$.

Proof.

By the definition of $C_{2}(x,y)$, we have

[TABLE]

∎

Theorem 2 shows that $x^{\ast}$ is the optimal solution of $\displaystyle\min_{x\in\Re^{n}}C_{2}(x,x^{\ast})$ as long as $x^{\ast}$ solves $\displaystyle\min_{x\in\Re^{n}}C_{1}(x)$. Moreover, combined with Lemma 1 and Theorem 1, we can immediately conclude that the thresholding representation of $(QP_{a}^{\lambda})$ can be given by

[TABLE]

where the thresholding operator $G_{a,\lambda\mu}$ is obtained in Definition 1 by replacing $\lambda$ with $\lambda\mu$.

Corollary 1**.**

For any fixed $\lambda>0$, $\mu>0$ and vector $x^{\ast}\in\Re^{n}$, let $x^{\ast}=G_{\lambda\mu,P}(B_{\mu}(x^{\ast}))$, then

[TABLE]

where the threshold value $t_{a,\lambda\mu}^{\ast}$ is obtained in $(\ref{equ17})$ by replacing $\lambda$ with $\lambda\mu$.

With the thresholding representations (23), the IFTA for solving the regularization problem $(QP_{a}^{\lambda})$ can be naturally defined as

[TABLE]

where $B_{\mu}(x^{k})=x^{k}+\mu F(x^{k})^{\ast}(b-F(x^{k})x^{k})$.

It is fairly well known that the quantity of the solution of a regularization problem depends seriously on the setting of the regularization parameter. Here, the cross-validation method is accepted to select the proper regularization parameter. Nevertheless, when some prior information is known for a regularization problem, this selection is more reasonable and intelligent. When doing so, the IFTA will be adaptive and free from the choice of the regularization parameter.

To make this selection clear, we suppose that the vector $x^{\ast}$ of sparsity $r$ is the optimal solution of the regularization problem $(QP_{a}^{\lambda})$, without loss of generality, we suppose that

[TABLE]

By Corollary 1, the following inequalities hold

[TABLE]

[TABLE]

According to $t_{a,\lambda\mu}^{2}\leq t_{a,\lambda\mu}^{1}$, we have

[TABLE]

which implies

[TABLE]

Above estimation helps to set the optimal regularization parameter. For convenience, we denote by $\lambda_{1}$ and $\lambda_{2}$ the left and the right of above inequality respectively.

[TABLE]

A choice of $\lambda$ is

[TABLE]

In practice, we approximate $B_{\mu}(x^{\ast})_{i}$ by $B_{\mu}(x^{k})_{i}$ in (27), and we can take

[TABLE]

in applications.

One more thing needs to be mentioned here is that the threshold value

[TABLE]

Notice that (27) is valid for any $\mu>0$ satisfying $0<\mu\leq\|F(x_{k})\|_{2}^{-2}$. In general, we can take $\mu=\mu_{k}=\frac{1-\epsilon}{\|F(x_{k})\|_{2}^{2}}$ with any small $\epsilon\in(0,1)$ below.

III Numerical experiments

In the section, we carry out a series of simulations to demonstrate the performance of IFTA, and compare them with those obtained with some state-of-art methods (iterative soft thresholding algorithm (ISTA)[6,17], iterative hard thresholding algorithm (IHTA)[6,17]). For each experiment, we repeatedly perform 30 tests and present average results and take $a=1$.

In our numerical experiments, we set

[TABLE]

where $A_{1}\in\Re^{30\times 100}$ is a fixed Gaussian random matrix, $x_{0}\in\Re^{100}$ is a reference vector, $f:[0,\infty)\rightarrow\mathcal{R}$ is a positive and smooth Lipschitz continuous function with $f(t)=\ln(t+1)$, $\eta$ is a sufficiently small scaling factor (we set $\eta=0.003$), and $A_{2}\in\Re^{30\times 100}$ is a fixed matrix with every entry equals to 1. Then the form of nonlinearity considered in (29) is a quasi-linear, and the more detailed accounts of the setting in form (29) can be seen in [6,17]. By randomly generating such sufficiently sparse vectors $x_{0}$ (choosing the non-zero locations uniformly over the support in random, and their values from $N(0,1))$, we generate vectors $b$. By this way, we know the sparsest solution to $F(x_{0})x_{0}=b$, and we are able to compare this with algorithmic results.

The stopping criterion is usually as following

[TABLE]

where $x_{k}$ and $x_{k-1}$ are numerical results from two continuous iterative steps and $\mathrm{Tol}$ is a given small number. The success is measured by the computing

[TABLE]

where $x^{\ast}$ is the numerical results generated by IFTA, and $\mathrm{Re}$ is also a given small number. In all of our experiments, we set $\mathrm{Tol}=10^{-8}$ to indicate the stopping criterion, and set $\mathrm{Re}=10^{-4}$ to indicate a perfect recovery of the original sparse vector $x_{0}$.

The graphs presented in Fig.3 and Fig.4 show the performance of the ISTA, IHTA and IFTA in recovering the true (sparsest) signals. From Fig.3, we can see that IFTA performs best, and IST algorithm the second. From Fig.4, we can get that the IFTA always has the smallest relative error value with sparsity growing.

IV Conclusion

In the paper, we take the fraction function as the substitution for $\ell_{0}$-norm in quasi-linear compressed sensing. An iterative fraction thresholding algorithm is proposed to solve the regularization problem $(QP_{a}^{\lambda})$ for all $a>0$. With the change of parameter $a>0$, our algorithm could get a promising result, which is one of the advantages for our algorithm compared with other algorithms. We also provide a series of experiments to assess performance of our algorithm, and the experiment results show that, compared with some state-of-art algorithms, our algorithm performs the best in the sparse signal recovery. However, the convergence of our algorithm is not proved theoretically in this paper, and it is our future work.