<figure><img src="image_1.png" /></figure>

HAL
open science

Generalized LASSO with under-determined
regularization matrices

Junbo Duan, Charles Soussen, David Brie, Jérôme Idier, Mingxi Wan,
Yu-Ping Wang

► To cite this version:

Junbo Duan, Charles Soussen, David Brie, Jérôme Idier, Mingxi Wan, et al. Generalized LASSO with under-determined regularization matrices. Signal Processing, 2016, 127, pp.239-246. 10.1016/j.sigpro.2016.03.001. hal-01295086

HAL Id: hal-01295086

https://hal.science/hal-01295086v1

Submitted on 7 Jan 2025

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Published in final edited form as:

*Signal Processing.* 2016 October; 127: 239–246. doi:10.1016/j.sigpro.2016.03.001.

# Generalized LASSO with under-determined regularization matrices

Junbo Duan<sup>a</sup>, Charles Soussen<sup>b</sup>, David Brie<sup>b</sup>, Jérôme Idier<sup>c</sup>, Mingxi Wan<sup>a</sup>, and Yu-Ping Wang<sup>d,e</sup>

<sup>a</sup>Key Laboratory of Biomedical Information Engineering of Ministry of Education and Department of Biomedical Engineering, School of Life Science and Technology, Xi’an Jiaotong University, Xi’an 710049, China

<sup>b</sup>Université de Lorraine and the CNRS at the CRAN (Centre de Recherche en Automatique de Nancy), UMR 7039, Campus Sciences, BP 70239, 54506 Vandoeuvre-lés-Nancy Cedex, France

<sup>c</sup>L’UNAM Université, Ecole Centrale de Nantes and the CNRS at the Institut de Recherche en Communications et Cybernétique de Nantes (IRCCyN), UMR 6597, 1 rue de la Noë, BP 92101, F-44321 Nantes Cedex 3, France

<sup>d</sup>Department of Biomedical Engineering, and the Department of Biostatistics and Bioinformatics, Lindy Boggs Center, Tulane University, New Orleans, LA 70118, USA

## Abstract

This paper studies the intrinsic connection between a generalized LASSO and a basic LASSO formulation. The former is the extended version of the latter by introducing a regularization matrix to the coefficients. We show that when the regularization matrix is even- or under-determined with full rank conditions, the generalized LASSO can be transformed into the LASSO form *via* the Lagrangian framework. In addition, we show that some published results of LASSO can be extended to the generalized LASSO, and some variants of LASSO, *e.g.*, robust LASSO, can be rewritten into the generalized LASSO form and hence can be transformed into basic LASSO.

Based on this connection, many existing results concerning LASSO, *e.g.*, efficient LASSO solvers, can be used for generalized LASSO.

## Keywords

LASSO; generalized LASSO; robust LASSO; deconvolution; solution path; diagonally dominant; total variation

<footnote><sup>a</sup>Corresponding author.</footnote>
# 1. Introduction

The least absolute shrinkage and selection operator (LASSO) [1] is one of the most popular approaches for sparse linear regression in the last decade, which is usually formulated as

$$x^*(\lambda) = \arg \min_{x} \left\{ \frac{1}{2} \|y - Ax\|^2 + \lambda \|x\|_1 \right\}, \quad (\text{LASSO})$$

where $y \in \mathbb{R}^n$ gathers $n$ observed measurements; $A \in \mathbb{R}^{m \times p}$ contains $p$ predictors of dimension $\dot{x} \in \mathbb{R}^p$ contains $p$ coefficients; $\| \cdot \|$ and $\| \cdot \|_1$ stand for the $\ell_2$- and $\ell_1$-norm respectively; and $\lambda > 0$ is the regularization parameter, controlling the tradeoff between the data fidelity and the model complexity. The most attracting feature of LASSO is the use of $\ell_1$-norm regularization, which yields sparse coefficients. The $\ell_1$-norm regularization results in the piecewise linearity [2] of the solution path $\{x^*(\lambda)|\lambda \in (0, +\infty)\}$ (i.e., the set of solutions with respect to continuous change of the regularization parameter $\lambda$), allowing for efficient reconstruction of the whole solution path. Based on this property, famous path tracking algorithms such as the least angle regression (LARS) [3] and the homotopy [2, 4] have been developed. LARS is a greedy algorithm working with decreasing $\lambda$ value. At each iteration, a dictionary atom is selected and appended to the set previously selected, and the next critical $\lambda$ value is computed. Homotopy is an extension of LARS, performing both forward and backward selection of the atoms already selected.

The generalized LASSO [5] (or analysis [6, 7], least mixed norm [8]) extends the basic LASSO (or synthesis) by imposing a regularization matrix $D \in \mathbb{R}^{m \times p}$ on the coefficient vector $x^*$:

$$x^*(\lambda) = \arg \min_{x} \left\{ \frac{1}{2} \|y - Ax\|^2 + \lambda \|Dx\|_1 \right\}, \quad (\text{generalized LASSO})$$

where $D$ typically contains the prior knowledge (e.g., structure information) [9] about $x$. For example, if $x$ is expected to be a piecewise constant signal (i.e., implying that its first order derivative is sparse), then $D$ is taken to be the first order derivative operator [10]. Some variants of LASSO can be regarded as the generalized LASSO by forming a structured matrix $D$. For example, the fused LASSO proposed by Tibshirani et al. [10] imposes the $\ell_1$ regularization on both the coefficients and their first order derivatives to encourage the solution to be locally constant. When one of the two regularization parameters is fixed, a fused LASSO problem can be rewritten as a generalized LASSO by cascading a scaled identity matrix with a first order derivative matrix to form the $D$ matrix. The graph-guided fused LASSO [11, 12] incorporates the network prior information into $D$ for correlation structure. Application of generalized LASSO can be found in image restoration [8], visual recognition [13], electroencephalography (EEG) [14], bioinformatics [15], ultrosonics [16], etc.

Since the LASSO has been proposed, results concerning the recovery condition [3, 17], solution path property [18], degree of freedom [19], model selection consistency [20], efficient algorithms as well as software development [2, 4] have been widely studied. One
may want to know whether these results can be applicable to a generalized LASSO problem, e.g., solving a generalized LASSO problem with a LASSO solver. An immediate example is when $D$ is a full rank square (hence invertible) matrix. By a simple change of variables $u = Dx$, the original generalized LASSO problem can be transformed into the basic LASSO form with a predictor matrix $AD^{-1}$ and a coefficient vector $u$. Therefore it can be solved by calling the LASSO subroutine.

Although the generalized LASSO and LASSO have been studied from various aspects [5, 7, 21], their connections are not fully explored, which is the main focus of this paper. Elad et al. [6] showed that they are equivalent, but confined the discussion to the denoising case, where $A$ is an identity matrix. Tibshirani and Taylor [5] showed that the former can be transformed to the latter; however, their method needs to introduce a matrix $D_0$ (see Appendices), which brings other potential questions as discussed in the conclusion.

The paper is organized as follows: in Section II, we show that when the regularization matrix is even- or under-determined with full rank conditions, the generalized LASSO can be transformed into the LASSO form via the Lagrangian framework. Based on this formula, in Section III we show that some published results of LASSO can be extended to the generalized LASSO. In Section IV, two variants of LASSO, namely the regularized deconvolution and the robust LASSO are analyzed under the generalized LASSO framework. We conclude the paper in Section IV.

## 2. Condition and formula of transformation

The simplification of the generalized LASSO depends on the setting of $D$ [6]. For the even-determined case ($m = p$)<fnref n="1" />, if $D$ has full rank, by a simple change of variables $u = Dx$, the original generalized LASSO problem can be transformed into the basic LASSO with a predictor matrix $AD^{-1}$ and a coefficient vector $u$. Once the LASSO solution $\hat{u}(\lambda)$ is known, the original generalized LASSO solution is immediate: $x(\lambda) = D^{-1}\hat{u}(\lambda)$.

For the over-determined case ($m > p$), if we change the variable $x$ to $u$, then $u$ has a higher dimension than $x$, or higher degree of freedom. If only a solution $\hat{u}$ satisfies the constraint $DD^{\dagger}\hat{u} = \hat{u}$ where $D^{\dagger}$ is the Moore-Penrose pseudoinverse of $D$, the original solution is guaranteed to be found [6]. Therefore, the generalized LASSO under this case cannot be transformed into a basic LASSO but a LASSO with equality constraint [22].

We consider the under-determined case ($m < p$). The following Theorem 1 states that a generalized LASSO problem can be transformed into a basic LASSO form under some conditions.

<footnote n="1">¹In general, the term even- over- and under-determination are used to describe the predictor matrix $A$. In this paper, we borrow these terms to characterize the regularization matrix $D$.</footnote>
**Theorem 1**

$\Phi = \begin{bmatrix} A \\ D \end{bmatrix}$ has full column rank (implying $m+n \ge p$), and $D$ has full row rank (implying $m \le p$), then the generalized LASSO problem can be transformed into the following LASSO form:

$$u^*(\lambda) = \arg \min_{u} \left\{ \frac{1}{2} \|z - H u\|^2 + \lambda \|u\|_1 \right\}, \quad (1)$$

where

$$H = Q_1 Q_2^\dagger \quad (2)$$

$$z = (I_n - Q_1(I_p - Q_2^\dagger Q_2)Q_1^T)y$$

and $Q_1 \in \mathbb{R}^{n \times p}$ and $Q_2 \in \mathbb{R}^{n \times p}$ are defined from the (unique) QR decomposition

$$\Phi = QR = \begin{bmatrix} Q_1 \\ Q_2 \end{bmatrix} R. \quad (3)$$

**Proof**—Since $\Phi$ has full column rank, its QR decomposition (3) is unique and the square matrix $R \in \mathbb{R}^{p \times p}$ is invertible. Moreover, $Q_2 \in \mathbb{R}^{n \times p}$ is full row rank because $D$ is full row rank. (3) rereads:

$$A = Q_1 R, \quad (4)$$

$$D = Q_2 R. \quad (5)$$

From $Q^T Q = I_p$ we have

$$Q_1^T Q_1 + Q_2^T Q_2 = I_p. \quad (6)$$

The generalized LASSO problem is equivalent to the following constrained optimization problem

$$(x^*(\lambda), u^*(\lambda)) = \arg \min_{x, u} \left\{ \frac{1}{2} \|y - Ax\|^2 + \lambda \|u\|_1 \right\}$$
subject to $Dx=u$. (7)

The Lagrange function associated with (7) reads

$$L(x, u, \mu) = \frac{1}{2} \|y - Ax\|^2 + \lambda \|u\|_1 + \mu^T (u - Dx),$$

where $\mu$ gathers $m$ Lagrange multipliers. The optimality condition reads

$$\frac{\partial L(x, u, \mu)}{\partial x} = A^T(Ax - y) - D^T\mu = 0,$$

$$\frac{\partial L(x, u, \mu)}{\partial \mu} = u - Dx = 0.$$

From these two equations, we have the following system

$$\begin{bmatrix} A^T A & -D^T \\ D & 0 \end{bmatrix} \begin{bmatrix} x \\ \mu \end{bmatrix} = \begin{bmatrix} A^T y \\ u \end{bmatrix}. \quad (8)$$

If matrix

$$M \triangleq \begin{bmatrix} A^T A & -D^T \\ D & 0 \end{bmatrix} \quad (9)$$

is invertible, we can find the equivalent form of the generalized LASSO problem. By substituting $A$ and $D$ with (4) and (5) respectively, we have

$$M = \begin{bmatrix} R^T Q_1^T Q_1 R & -R^T Q_2^T \\ Q_2^T R & 0 \end{bmatrix} = \begin{bmatrix} R^T & 0 \\ 0 & I_m \end{bmatrix} \begin{bmatrix} Q_1^T Q_1 & -Q_2^T \\ Q_2 & 0 \end{bmatrix} \begin{bmatrix} R & 0 \\ 0 & I_m \end{bmatrix},$$

where $I_m \in \mathbb{R}^{m \times m}$ is the identity matrix. Since $R$ and $I_m$ are full rank, the invertibility of $M$ is equivalent to that of

$$N \triangleq \begin{bmatrix} Q_1^T Q_1 & -Q_2^T \\ Q_2 & 0 \end{bmatrix}.$$

From (6),
$$N = \begin{bmatrix} I_p - Q_2^T Q_2 & -Q_2^T \\ Q_2 & 0 \end{bmatrix} = I_{p+m} + \begin{bmatrix} Q_2^T & 0 \\ 0 & I_m \end{bmatrix} \begin{bmatrix} -I_m & -I_m \\ I_m & -I_m \end{bmatrix} \begin{bmatrix} Q_2 & 0 \\ 0 & I_m \end{bmatrix}.$$

From Woodbury's matrix identity [23, p.141], N is invertible if and only if

$$W \triangleq \begin{bmatrix} -I_m & -I_m \\ I_m & -I_m \end{bmatrix}^{-1} + \begin{bmatrix} Q_2 & 0 \\ 0 & I_m \end{bmatrix} \begin{bmatrix} Q_2^T & 0 \\ 0 & I_m \end{bmatrix} = \begin{bmatrix} Q_2 Q_2^T - \frac{1}{2} I_m & \frac{1}{2} I_m \\ -\frac{1}{2} I_m & \frac{1}{2} I_m \end{bmatrix}$$

is invertible. Since

$$\det(W) = \det(\frac{1}{2} Q_2 Q_2^T),$$

and $\det(Q_2 Q_2^T) \neq 0$ ($Q_2$ is full rank), the invertibility of $W$ is guaranteed, and therefore $M$ is invertible.

From the block matrix inversion lemma [23, p.108], we have

$$W^{-1} = \begin{bmatrix} (Q_2 Q_2^T)^{-1} & -(Q_2 Q_2^T)^{-1} \\ (Q_2 Q_2^T)^{-1} & 2 I_m - (Q_2 Q_2^T)^{-1} \end{bmatrix}.$$

Therefore,

$$N^{-1} = I_{p+m} - \begin{bmatrix} Q_2^T & 0 \\ 0 & I_m \end{bmatrix} W^{-1} \begin{bmatrix} Q_2 & 0 \\ 0 & I_m \end{bmatrix},$$

and

$$M^{-1} = \begin{bmatrix} R^{-1} & 0 \\ 0 & I_m \end{bmatrix} N^{-1} \begin{bmatrix} R^{-T} & 0 \\ 0 & I_m \end{bmatrix} = \begin{bmatrix} R^{-1}(I_p - Q_2^T(Q_2 Q_2^T)^{-1}Q_2)R^{-T} & R^{-1}Q_2^T(Q_2 Q_2^T)^{-1} \\ -(Q_2 Q_2^T)^{-1}Q_2 R^{-T} & (Q_2 Q_2^T)^{-1} - I_m \end{bmatrix} \\ = \begin{bmatrix} R^{-1}(I_p - Q_2^\dagger Q_2)R^{-T} & R^{-1}Q_2^\dagger \\ -(Q_2^\dagger)^T R^{-T} & (Q_2 Q_2^T)^{-1} - I_m \end{bmatrix},$$

where $Q_2^\dagger = Q_2^T(Q_2 Q_2^T)^{-1}$. Finally, from (8) we have

$$x = R^{-1}(I_p - Q_2^\dagger Q_2)Q_1^T y + R^{-1}Q_2^\dagger u, \quad (10)$$

and
$$ \mu = -(Q_2^\dagger)^T Q_1^T y + ((Q_2 Q_2^\dagger)^{-1} - I_m)u. $$

By substituting (10) into (7), we obtain the results (1) and (2).

Once the solution path {$\hat{u}(\lambda)|\lambda \in (0, +\infty)$} is known, the corresponding solution path {$x(\hat{\lambda})|\lambda \in (0, +\infty)$} can be calculated according to (10).

By substituting $Q_1 = AR^{-1}$ and $Q_2 = DR^{-1}$ into (2), and utilizing

$$ \Psi \triangleq \Phi^T \Phi = A^T A + D^T D = R^T R, \quad (11) $$

up to a few manipulations, an useful equivalent form of (2) without QR decomposition is obtained as

$$ H = A \Psi^{-1} D^T (D \Psi^{-1} D^T)^{-1} \quad (12) $$

$$ z = (I_n - A \Psi^{-1} A^T + H D \Psi^{-1} A^T) y $$

And from (10), $x$ can be recovered as

$$ x = \Psi^{-1} (A^T - D^T H^T) y + \Psi^{-1} D^T (D \Psi^{-1} D^T)^{-1} u \quad (13) $$

Note that $H$ has the same dimensionality as $AD^\dagger$, but normally they are not equal unless $\Psi$ is a diagonal matrix.

### Remark 1

There are many real-world examples that satisfy this condition. e.g., the total variation denoising problem where $A$ and $D$ are an identity and a first order derivative matrix, respectively. By substituting $A$ and $D$ into $M$ in (9), one can check that $M$ is invertible, and hence satisfies the condition in Theorem 1. Another two examples are the regularized deconvolution problem and the robust LASSO, which will be presented in detail in Sec. 4.

### Remark 2

Under the case $m + n = p$, $\Phi$ is a full rank (square) matrix. Therefore, the row vectors of $Q$ are orthogonal, and $Q_1 Q_2^T = 0$. From (1)-(2), we have $H = 0$ and the solution is obviously $u = 0$. According to (10), the original solution is $x = R^{-1} Q_1^T y$. This solution does not depend on $\lambda$, and the solution path {$x(\hat{\lambda})|\lambda \in (0, +\infty)$} is a singleton.

Another viewpoint is obtained by a change of variables. Defining $v = Rx$, the generalized LASSO problem rereads
$$v^*(\lambda) = \arg \min_{v} \left\{ \frac{1}{2} \|y - Q_1 v\|^2 + \lambda \|Q_2 v\|_1 \right\}. \quad (14)$$

Since the rows of $Q_1$ and $Q_2$ form an orthonormal basis, it is obvious that the minimum $\ell_2$-norm least square solution $\hat{v}=Q_1^\dagger y=Q_1^T y$ satisfies $Q_2 \hat{v} = 0$. This means that $v$ minimizes both the $\ell_2$- and $\ell_1$-term simultaneously. Therefore, the cost function in (14) is equal to 0 whenever $\lambda > 0$.

### Remark 3

Under the case $m+n<p$, $\mathrm{Null}(\Phi) \neq \emptyset$, therefore one can check that the generalized LASSO does not yield unique solution, and there are unlimited number of solutions

$$\{x^*|x^*=\Phi^\dagger \begin{bmatrix} y \\ 0 \end{bmatrix} + x_0, x_0 \in \mathrm{Null}(\Phi)\} \quad \text{with } Dx^*=0 \text{ and } Ax^*=y \text{ whatever } \lambda > 0.$$

Note that this case does not satisfy the condition in Theorem 1 since $\Phi$ does not have full column rank. This case serves to complete the discussion.

## 3. Extension of existing LASSO results

Since the LASSO has been intensively studied during the last decade, many results concerning the computational issue have been published. In this section, we extend some of them into the generalized LASSO problem.

### 3.1. Monotonic solution path property

Our earlier result [18, Theorem 2] shows that if $(H^T H)^{-1}$ is diagonally dominant, then the solution path of the corresponding LASSO problem changes monotonically with respect to $\lambda$. In other words, the solution path has at most $m$ segments. Therefore, the ‘forward’ algorithm, LARS, yields the same solution path as the ‘forward-backward’ algorithm, homotopy, and both of them can recover the complete solution path within $m$ iterations. Based on this fact, the computational complexity can thus be reduced. The following corollaries extend this result to the generalized LASSO.

**Corollary 1**—If a generalized LASSO problem satisfying condition in Theorem 1 also satisfies that $((D\Psi^{-1}D^T)^{-1} - I)^{-1}$ is diagonally dominant, then the complete solution path can be recovered within $m$ iterations of LARS and homotopy algorithm.

**Proof:** When a generalized LASSO problem satisfying condition in Theorem 1, it can be transformed into a LASSO problem with $H$ and $z$ being defined in (2). By utilizing equations (4), (5), (6) and (11), we have

$$H^T H = (Q_2 Q_2^T)^{-1} - I = (D \Psi^{-1} D^T)^{-1} - I,$$

therefore, $(H^T H)^{-1}$ is diagonally dominant. From [18], this corollary is straightforward.
There exist many matrices satisfying above condition. Obvious examples are the orthogonal dictionaries like Dirac basis or Hadamard basis. By Monte Carlo simulation, we study the probability of random matrices to satisfy the condition. p is fixed to 100, and for each given configuration (n, m), 1000 trials A and D are generated, whose entries obey i.i.d. normal distribution or Bernoulli distribution with parameter b = 0.1 (the probability for 1 is b, for 0 is 1 - b). The frequency of A and D satisfying condition in Corollary 1 is shown in Fig. B.1. From the simulation results, it is shown that random matrices satisfy the condition when n ≫ m.

When $A^T A$ forms an orthogonal design, i.e., $A^T A = I$, above condition can be further simplified to the following corollary.

**Corollary 2**—If a generalized LASSO problem satisfying condition in Theorem 1 also satisfies that A forms an orthogonal design, and $DD^T$ is diagonally dominant, then the complete solution path can be recovered within m iterations of LARS and homotopy algorithm.

**Proof:** Since D has full row rank, $DD^T$ is invertible. By utilizing Searle's set of matrix identities [24, page 151], we can verify that

$$D(I+D^T D)^{-1} D^T = (I+(DD^T)^{-1})^{-1}.$$

When $A^T A = I$, we have

$$(H^T H)^{-1} = ((D \Psi^{-1} D^T) - I)^{-1} = ((D(I+D^T D)^{-1} D^T)^{-1} - I)^{-1} = DD^T$$

Therefore, condition in Corollary 1 can be simplified into checking the diagonal dominance of $DD^T$.

Corollaries 1 and 2 show that in the best case, the computational complexity of homotopy can reach as low as m iterations to compute the complete regularization path on $\lambda \in (0, +\infty)$. It is shown [25] that in the worst case, the computation complexity can reach as high as $(3^m - 1)/2$. Therefore, the computation complexity highly varies depending on the specific problem. Instead of using a generic LASSO solver such as LARS, fast and efficient algorithms can be developed for some specific problems. e.g., the algorithm in [26] permits fast computation and low storage for the regularized deconvolution problem introduced later in Sec. 4.1 with $\alpha = 0$ and $\beta = -1$. The computational advantage is significant when the dimension is extremely large.

## 3.2. Decoupling into 1-dimensional subproblems

Efron et al. [3] and Friedman et al. [27] have shown that if the dictionary A forms an orthogonal design (or uncorrelated design), i.e., $A^T A = I$, a LASSO problem can be decoupled into a set of subproblems, and then solved separately.
To be more specific, when $A^T A = I$ the LASSO criterion can be rewritten as

$$ \frac{1}{2} \|A^T y - x\|^2 + \lambda \|x\|_1 + y^T y - y^T A A^T y. $$

This shows that the p-dimensional LASSO problem can be decoupled into p' 1-dimensional subproblems, and each subproblem can be solved separately. The closed-form solution reads

$$ x^* = \eta_{\lambda}(A^T y), \quad (15) $$

where $\eta_{\lambda}(\cdot)$ is the componentwise soft thresholding function [28]:

$$ \eta_{\lambda}(x) = \begin{cases} x - \lambda, & x > \lambda; \\ 0, & |x| \le \lambda; \\ x + \lambda, & x < -\lambda. \end{cases} $$

Similarly, for generalized LASSO problem the following corollary provides a sufficient condition to guarantee decoupling.

**Corollary 3—If a generalized LASSO problem satisfying condition in Theorem 1 also satisfies**

$$ D\Psi^{-1}D^T = \frac{1}{2}I, \quad (16) $$

then it can be decoupled into p 1-dimensional LASSO problems, with the following closed-form solution

$$ u^* = \eta_{\lambda}(2D\Psi^{-1}A^T y). \quad (17) $$

**Proof:** Since A and D satisfy condition in Theorem 1, this generalized LASSO problem can be transformed into a LASSO problem, with H, z being defined in (2), and matrices $Q_1$, $Q_2$ and R being defined in the proof of Theorem 1.

By substituting equations (4), (5) and (11) into (16), an equivalent condition reads

$$ Q_2 Q_2^T = \frac{1}{2} I. \quad (18) $$

Therefore, $Q_2^\dagger = 2Q_2^T$.
From (2) and by utilizing (6), it is easy to prove that the Gramian matrix of $\mathbf{H}$ is an identity matrix, i.e., $\mathbf{H}^T\mathbf{H} = \mathbf{I}$, so $\mathbf{H}$ forms an orthogonal design, indicating the transformed LASSO problem can be decoupled.

$$ \left[ \begin{array}{ccccccc} 1 & \beta - \alpha & -\alpha\beta & 0 & \dots & 0 \\ \beta - \alpha & 1+(\beta - \alpha)^2 & (\beta - \alpha)(1 - \alpha\beta) & -\alpha\beta & \dots & \vdots \\ -\alpha\beta & (\beta - \alpha)(1 - \alpha\beta) & 1+\alpha^2\beta^2+(\beta - \alpha)^2 & (\beta - \alpha)(1 - \alpha\beta) & \dots & 0 \\ 0 & \dots & \dots & \dots & \dots & -\alpha\beta \\ \vdots & \dots & -\alpha\beta & (\beta - \alpha)(1 - \alpha\beta) & 1+\alpha^2\beta^2+(\beta - \alpha)^2 & (\beta - \alpha)(1 - \alpha\beta) \\ 0 & \dots & 0 & -\alpha\beta & (\beta - \alpha)(1 - \alpha\beta) & 1+\alpha^2\beta^2+(\beta - \alpha)^2 \end{array} \right] \quad (20) $$

Further more, with the help of (18), we can also prove

$$ \mathbf{H}^T \mathbf{z} = 2\mathbf{Q}_2 \mathbf{Q}_1^T \mathbf{y} = 2\mathbf{D} \mathbf{\Psi}^{-1} \mathbf{A}^T \mathbf{y}. \quad (19) $$

Finally, by replacing $\mathbf{A}$ and $\mathbf{y}$ in (15) with $\mathbf{H}$ and $\mathbf{z}$ respectively, and substituting (19), the closed-form solution reads in (17).

**Remark 4**—If $\mathbf{A}^T\mathbf{A}$ is invertible, by utilizing Searle's set of matrix identities [24, page 151], an alternative condition of (16) reads

$$ \mathbf{D}(\mathbf{A}^T \mathbf{A})^{-1} \mathbf{D}^T = \mathbf{I} $$

Note that this condition is stronger than (16) since for under-determined $\mathbf{A}$, $\mathbf{A}^T\mathbf{A}$ is not invertible.

One can check that condition in Corollary 3 fulfills when $\mathbf{D}$ is a square matrix with full rank and $\mathbf{A}^T\mathbf{A} = \mathbf{D}^T\mathbf{D}$.

## 4. Two examples on the analysis of LASSO variants

Several variants of LASSO can be unified under the generalized LASSO framework, such as the total variation regularized deconvolution for signal and image recovery [29], the robust LASSO for face recognition and sensor network [30], the adaptive LASSO for variable selection [31], the fused LASSO for gene expression data analysis [10], the $l_1$ trend filtering for time series analysis [32], and the adaptive generalized fused LASSO for road-safety data analysis [33]. In this section, the first two will be discussed.
## 4.1. Regularized deconvolution

Deconvolution is a challenging problem, which can be formulated as a LASSO problem with a Toeplitz predictor matrix. Total variation [34] is a regularization tool to restore piecewise constant signals. This leads to the generalized LASSO form in which A is a Toeplitz matrix containing the convolution kernel and D is also a Toeplitz matrix containing the total variation kernel. For high dimensional problems, the computation burden is a major factor. The following theorem sheds some lights on this issue.

**Corollary 4**—Let A and D be lower square triangle Toeplitz matrices with the first column [1, α, α², α³, ..., αⁿ⁻¹]ᵀ and [1, β, 0, 0, ..., 0]ᵀ, respectively. If α and β fall within the region displayed in Fig. B.2, then the generalized LASSO problem can be solved for all λ by the LARS and homotopy algorithms in at most n iterations.

The total variation regularization corresponds to the special case when $β = -1$.

**Proof:** According to Theorem 1, since D is invertible, this generalized LASSO problem is equivalent to the basic LASSO problem (1) with $H = AD^{-1}$. Since A is squared and invertible, the Gramian matrix inverse reads $(H^T H)^{-1} = DA^{-1}A^{-T}D^T$. So in the following we will show that this matrix is diagonally dominant.

One can check that $A^{-1}$ is also a lower triangle Toeplitz matrix with the first column as [1, -α, 0, 0, ..., 0]ᵀ. Then the Gramian matrix $DA^{-1}A^{-T}D^T$ reads as (20). In order to have the diagonal dominance property, the following three inequalities are required, which correspond to the first three rows in matrix (20) (The conditions related to the other rows are obviously implied by the third inequality):

$$1 \geq |\beta - \alpha| + |\alpha\beta|,$$

$$1+(\beta - \alpha)^2 \geq |\beta - \alpha|+|\beta - \alpha||1 - \alpha\beta|+|\alpha\beta|,$$

$$1+\alpha^2\beta^2+(\beta - \alpha)^2 \geq 2|\alpha\beta|+2|\beta - \alpha||1 - \alpha\beta|.$$

Those inequalities yield three regions in the 2D plane (α, β) whose intersection is displayed on Fig. B.2.

Corollary 4 in this paper, Theorem 2 in [18], Theorem 1 and Corollary 2 in [17] jointly show that for deconvolution problem with total variation regularization ($β = -1$), if the convolution kernel is high pass ($-1 \le \alpha < 0$) or all pass ($α = 0$), but not low pass ($0 < \alpha \le 1$), the true solution can be recovered via the generalized LASSO. If there is no total variation constraint ($β = 0$), then any convolution kernel having exponential attenuation ([1, α, α², α³, ..., αⁿ⁻¹]ᵀ, -1 ≤ α ≤ 1) admits the perfect recovery.
## 4.2. Robust LASSO

LASSO with the first order total variation method can recover piecewise constant signals. However, if there are outliers in the signal $y$, this approach tends to introduce false detections. A second component signal $s \in \mathbb{R}^n$ is introduced to cope with the outliers in the observation [30]. Since outliers are usually due to burst error, which are assumed to obey the Laplace distribution, a second $\ell_1$ regression term is incorporated into the objective function, yielding the robust generalized LASSO as follows,

$$ (x^*, s^*) = \arg \min_{x,s} \left\{ \frac{1}{2} \|y - Ax - s\|^2 + \lambda_1 \|Dx\|_1 + \lambda_2 \|s\|_1 \right\}. $$

This problem can be rewritten into a generalized LASSO form with

$$ A_{\text{robust}} = [A, I], D_{\text{robust}} = \begin{bmatrix} D & 0 \\ 0 & \tau I \end{bmatrix}, \text{ and } x_{\text{robust}} = \begin{bmatrix} x \\ s \end{bmatrix}, \text{ where } \tau = \lambda_2/\lambda_1. \text{ When} $$

$\Phi = \begin{bmatrix} A \\ D \end{bmatrix}$ has full column rank, $D$ has full column rank, and $\tau > 0$, one can check that

$$ \Phi_{\text{robust}} = \begin{bmatrix} A_{\text{robust}} \\ D_{\text{robust}} \end{bmatrix} \text{ has full column rank, } D_{\text{robust}} \text{ has full row rank, and } M \text{ is invertible.} $$

Therefore, from Theorem 1 the robust LASSO problem can be transformed into a basic LASSO form.

## 5. Conclusion

This paper discusses the simplification of a generalized LASSO problem into the basic LASSO form. When the regularization matrix $D$ is even- or under-determined ($m \le p$), we showed that this simplification is possible. Otherwise, there is no guarantee that this simplification can be done. In the former case, optimization tools dedicated to LASSO can be straightforwardly applied to the generalized LASSO.

Tibshirani and Taylor [5] gave a simple way to transform a generalized LASSO to the basic LASSO form when $D$ is not a square matrix ($m < p$). As shown in Appendix Appendix A,

they introduced a matrix $D_0$ to form $\tilde{D} = \begin{bmatrix} D \\ D_0 \end{bmatrix}$, a square matrix $\tilde{D}$, which is invertible. This poses the question on whether their results depend on $D_0$ or not. Appendix Appendix B shows that their method yields the equivalent results as ours, and this indicates that $H$ and $z$ in their formula do not depend on $D_0$. In addition, an improperly introduced $D_0$ may incur potential numerical error.

The proposed formula reposes directly on $A$ and $D$, therefore gives an insight into their interaction in the generalized LASSO problem. Based on the proposed formula, it is shown that existing results related to LASSO can be extended to the generalized LASSO. Since some variants of LASSO can be unified under the generalized LASSO framework, they can be transformed into the basic LASSO, and hence efficient LASSO solvers can be applied for
the solution. Furthermore, under this framework, many different types of regression formulations such as the trend filtering in a recent study [35] can be unified.

Acknowledgments

This study was partially supported by National Science Foundation of China (No. 61401352), China Postdoctoral Science Foundation (No. 2014M560786), Shaanxi Postdoctoral Science Foundation, Fundamental Research Funds for the Central Universities (No. xjj2014060), National Science Foundation (No. 1539067), and National Institute of Health (No. R01MH104680, R01MH107354, and R01GM109068).

Appendix A

Tibshirani and Taylor's transformation method

Tibshirani and Taylor [5] proposed the following method to transform a generalized LASSO when rank(D) = m < p.

First a matrix $D_0 \in \mathbb{R}(p-m) \times p$ is found, whose rows are orthogonal to these in D. Then

matrix $\tilde{D} = \begin{bmatrix} D \\ D_0 \end{bmatrix} \in \mathbb{R}^{p \times p}$ is a square matrix with full rank, and thus invertible.

For further use, let us define $\tilde{D}x = \begin{bmatrix} Dx \\ D_0x \end{bmatrix} \triangleq \begin{bmatrix} u \\ v \end{bmatrix}$ with u and v are column vectors of

length m and p-m respectively, therefore $x = \tilde{D}^{-1} \begin{bmatrix} u \\ v \end{bmatrix}$, and the generalized LASSO criterion rereads

$$ \frac{1}{2} \|y - A\tilde{D}^{-1} \begin{bmatrix} u \\ v \end{bmatrix} \|^2 + \lambda \|u\|_1. $$

Let's also define

$$ A\tilde{D}^{-1} = [A_1 \quad A_2], $$

where $A_1$ and $A_2$ are matrices of size $n \times m$ and $n \times (p-m)$ respectively. Above criterion rereads

$$ \frac{1}{2} \|y - A_1 u - A_2 v\|^2 + \lambda \|u\|_1. \quad (A.1) $$

which is a least square problem with respect to v, and the solution is

$$ v^* = A_2^\dagger (y - A_1 u), $$

where $A_2^\dagger = (A_2^T A_2)^{-1} A_2^T$.
Finally, by substituting $v^*$ into (A.1), the generalized LASSO can be transformed into form (1) with

$$H=A_1 - A_2 A_2^\dagger A_1 \quad (A.2)$$

$$z=(I_n - A_2 A_2^\dagger)y$$

# Appendix B

## The equivalence between Tibshirani and Taylor's method and the proposed method in the current paper

In this appendix, we prove that $H$'s and $z$'s in form (2) and (A.2) are equal.

### Appendix B.1. Equivalence of H's

First, for further use let's define $Q_0 \triangleq D_0 R^{-1}$ and

$$\begin{bmatrix} Q_2 \\ Q_0 \end{bmatrix}^{-1} \triangleq [P_2 \quad P_0], \quad (B.1)$$

where $P_2$ and $P_0$ are of size $p \times m$ and $p \times (p-m)$, respectively.

It is easy to verify that

$$P_0^T Q_2^T = 0, \quad (B.2)$$

$$P_2 Q_2 + P_0 Q_0 = I, \quad (B.3)$$

$$Q_1 P_2 = A_1,$$

$$Q_1 P_0 = A_2. \quad (B.4)$$

From (B.4), (6) and (B.2), we have

$$A_2^T Q_1 Q_2^T = P_0^T Q_1^T Q_1 Q_2^T = P_0^T (I - Q_2^T Q_2) Q_2^T = 0. \quad (B.5)$$
Left multiply both sides of (B.3) with $Q_1$, and substitute resultant equation into (B.5), we have

$$A_2^T (A_1 Q_2 + A_2 Q_0) Q_2^T = 0.$$

Move the second part to the right side of equals sign, and left multiply both sides with $(A_2^T A_2)^{-1}$, right multiply with $(Q_2 Q_2^T)^{-1}$, we have

$$A_2^\dagger A_1 = -Q_0 Q_2^\dagger.$$

Left multiply both sides of above equation with $-A_2$, and add with $A_1$, we have

$$A_1 - A_2 A_2^\dagger A_1 = A_1 + A_2 Q_0 Q_2^\dagger.$$

The left side of above equation is $H$ in (A.2), and the right side is equal to $H$ in (2), since left multiply both sides of (B.3) with $Q_1$, and right multiply with $Q_2^\dagger$ yields

$$A_1 + A_2 Q_0 Q_2^\dagger = Q_1 Q_2^\dagger.$$

## Appendix B.2. Equivalence of z's

To prove that $z$'s in (2) and (A.2) are equal, we have to prove

$$Q_1 Q_1^T - Q_1 Q_2 Q_2^\dagger Q_1^T = A_2 A_2^\dagger.$$

From (B.1), one can verify that

$$I - Q_2^\dagger Q_2 = P_0 P_0^\dagger.$$

Left multiply both sides of above equation with $Q_1$, and right multiply with $Q_1^T$, we have

$$Q_1 Q_1^T - Q_1 Q_2 Q_2^\dagger Q_1^T$$

$$= Q_1 P_0 P_0^\dagger Q_1^T$$

$$= A_2 (P_0^T P_0)^{-1} A_2^T \quad \text{using Eq. (B.4)}$$
$$=A_2(P_0^T P_0 - P_0^T Q_2^T Q_2 P_0)^{-1} A_2^T \quad \text{using Eq. (B.2)}$$

$$=A_2(P_0^T(I - Q_2^T Q_2)P_0)^{-1} A_2^T$$

$$=A_2(P_0^T(Q_1^T Q_1)P_0)^{-1} A_2^T \quad \text{using Eq. (6)}$$

$$=A_2 A_2^\dagger \quad \text{using Eq. (B.4)}$$

## References

1. Tibshirani R. Regression shrinkage and selection via the Lasso. J. R Statist. Soc. B. 1996; 58(1): 267–288.

2. Osborne MR, Presnell B, Turlach BA. A new approach to variable selection in least squares problems. IMA Journal of Numerical Analysis. 2000; 20(3):389–403.

3. Efron B, Hastie T, Johnstone I, Tibshirani R. Least angle regression. Annals Statist. 2004; 32(2): 407–499.

4. Malioutov, DM.; Cetin, M.; Willsky, AS. Proc. IEEE ICASSP. Vol. V. Philadelphia, PA: 2005 Mar. Homotopy continuation for sparse signal representation; p. 733-736.

5. Tibshirani RJ, Taylor J. The solution path of the generalized lasso. Annals Statist. 2011; 39(3):1335–1371.

6. Elad M, Milanfar P, Rubinstein R. Analysis versus synthesis in signal priors. Inverse Problems. 2007; 23(3):947.

7. Vaiter S, Peyré G, Dossal C, Fadili J. Robust sparse analysis regularization. IEEE Trans. Inf. Theory. 2013; 59(4):2001–2016.

8. Fu H, Ng MK, Nikolova M, Barlow JL. Efficient minimization methods of mixed $ℓ_2-ℓ_1$ and $ℓ_1-ℓ_1$ norms for image restoration. SIAM J. Sci. Comput. 2006; 27(6):1881–1902.

9. Eldar, YC.; Kutyniok, G., editors. Compressed Sensing: Theory and Applications. Cambridge University Press; 2012 May.

10. Tibshirani R, Saunders M, Rosset S, Zhu J, Knight K. Sparsity and smoothness via the fused lasso. Journal Of The Royal Statistical Society Series B. 2005; 67(1):91–108.

11. Kim S, Sohn K-A, Xing EP. A multivariate regression approach to association analysis of a quantitative trait network. Bioinformatics. 2009; 25:i204–i212. ISMB 2009. [PubMed: 19477989]

12. Chen X, Lin Q, Kim S, Carbonell JG, Xing EP. Smoothing proximal gradient method for general structured sparse regression. The Annals of Applied Statistics. 2012; 6(2):719–752.

13. Morioka N, Satoh S. Generalized lasso based approximation of sparse coding for visual recognition. Advances in Neural Information Processing Systems. 2011:181–189.

14. Vega-Hernández M, Martínez-Montes E, Sánchez-Bornot JM, Lage-Castellanos A, Valdés-Sosa PA. Penalized least squares methods for solving the EEG inverse problem. Statistica Sinica. 2008; 18:1535–1551.

15. Tibshirani R, Wang P. Spatial smoothing and hot spot detection for CGH data using the fused lasso. Biostatistics. 2008; 9(1):18–29. [PubMed: 17513312]

16. Yu C, Zhang C, Xie L. A blind deconvolution approach to ultrasound imaging. IEEE Trans. Ultrasonics Ferroelectrics Frequency Control. 2012; 59(2):271–280.
17. Donoho DL, Tsaig Y. Fast solution of l₁-norm minimization problems when the solution may be sparse. IEEE Trans. Inf. Theory. 2008 Nov.54(11):4789-4812.

18. Duan J, Soussen C, Brie D, Idier J, Wang Y-P. On LARS/homotopy equivalence conditions for over-determined LASSO. IEEE Signal Processing Letters. 2012; 19(12):894-897.

19. Zou H, Hastie T, Tibshirani R. On the "degrees of freedom" of the lasso. Annals Statist. 2007; 35(5):2173-2192.

20. Zhao P, Yu B. On model selection consistency of lasso. Journal of Machine Learning Research. 2006; 7:2541-2563.

21. She Y. Sparse regression with exact clustering. Electronic Journal of Statistics. 2010; 4:1055-1096.

22. James, GM.; Paulson, C.; Rusmevichientong, P. Tech. Rep. University of Southern California; 2013. The constrained lasso. http://www-bcf.usc.edu/~gareth/research/CLassoFinal.pdf

23. Bernste, D. Matrix Mathematics. Princeton University Press; 2009.

24. Searle, SR., editor. Matrix Algebra Useful for Statistics. John Wiley & Sons, Inc.; 1982.

25. Mairal, J.; Yu, B. International Conference on Machine Learning. Edinburgh, Scotland, UK: 2012. Complexity analysis of the lasso regularization path; p. 353-360.

26. Duan, J.; Zhang, J-G.; Lefante, J.; Deng, H-W.; Wang, Y-P. IEEE international conference on bioinformatics and biomedicine workshops. Atlanta, GA, USA: 2011. Detection of copy number variation from next generation sequencing data with total variation penalized least square optimization; p. 3-12.

27. Friedman J, Hastie T, Höfling H, Tibshirani R. Pathwise co-ordinate optimization. The Annals of Applied Statistics. 2007; 1(2):302-332.

28. Donoho DL, Johnstone IM. Adapting to unknown smoothness via wavelet shrinkage. J. Amer. Statist. Assoc. 1995; 90(432):1200-1224.

29. Chambolle A, Lions P-L. Image recovery via total variation minimization and related problems. Numer. Math. 1997; 76:167-188.

30. Nguyen NH, Tran TD. Robust lasso with missing and grossly corrupted observations. IEEE Trans. Inf. Theory. 2013; 59(4):2036-2058.

31. Zou H. The adaptive lasso and its oracle properties. J. Amer. Statist. Assoc. 2006; 101(476):1418-1429.

32. Kim S-J, Koh K, Boyd S, Gorinevsky D. $l_1$ trend filtering. SIAM review. 2009; 51(2):339-360.

33. Viallon, V.; Lambert-Lacroix, S.; Hoeing, H.; Picard, F. Tech. Rep. Université de Lyon; 2013. Adaptive generalized fused-lasso: Asymptotic properties and applications. http://www.statistics.gov.hk/wsc/IPS012-P1-S.pdf

34. Rudin L, Osher S, Fatemi C. Nonlinear total variation based noise removal algorithm. Physica D. 1992; 60:259-268.

35. Tibshirani RJ. Adaptive piecewise polynomial estimation via trend filtering. Annals Statist. 2014; 42(1):285-323.
<figure><img src="image_3.png" /><figcaption>Figure B.1.</figcaption></figure>

The frequency (in percentage) of random matrices *A* and *D* satisfying condition in Corrolary 1 with *p* = 100. (a) Entries in *A* and *D* obeys normal distribution and Bernoulli distribution with parameter *b* = 0.1, respectively. (b) The distributions of *A* and *D* are exchanged.
<figure><img src="image_4.png" /><figcaption>Figure B.2.</figcaption></figure>

The feasible region that admits matrix (12) to be diagonal dominant. The region is mirror symmetric with respect to line $α = β$, and the curve in the south-east quadrant has expression

$$ \beta = \frac{-(\alpha^2+1)+2\sqrt{\alpha(\alpha^2-\alpha+1)}}{(\alpha-1)^2}. $$