Introduction

Suppose you want to find a sparse solution to an underdetermined system of equations — fewer unknowns on are nonzero than the total number of variables. This is the central challenge in compressed sensing, medical imaging and signal processing. The natural objective adds a smooth data-fit term and a nonsmooth L1 penalty that promotes sparsity:

F(x)=12Axb2f(x) smooth+λx1g(x) nonsmoothF(x) = \underbrace{\tfrac{1}{2}\|Ax - b\|^2}_{f(x)\text{ smooth}} + \underbrace{\lambda \|x\|_1}_{g(x)\text{ nonsmooth}}

Gradient descent works beautifully on the smooth part ff, but it cannot handle gg because the L1 norm has corners — points where no gradient exists. Proximal gradient methods split the problem cleanly: take a gradient step in ff, then solve a tiny standalone minimization problem for gg. That second step is called the proximal operator, and for the L1 norm it has a simple closed form: soft-thresholding.

The basic algorithm is ISTA (Iterative Shrinkage-Thresholding Algorithm). Its convergence rate is O(1/k)O(1/k). In 2009, Beck and Teboulle showed that adding a single momentum term — extrapolating the current iterate using the previous one — yields FISTA (Fast ISTA) with rate O(1/k2)O(1/k^2), at essentially zero extra cost.

Try It: Sparse Recovery

Below is a 1-D sparse recovery problem: the algorithm must recover a sparse vector xx from noisy measurements. Use the sliders to change the L1 penalty λ\lambda and the step size tt, then press Run to watch ISTA and FISTA converge. The left panel tracks the objective value; the right panel shows the current estimate vs. the true sparse signal.

<!-- {{c_intro}} -->
<div class="controls">
  <label>{{lbl_lambda}} <span id="lval">0.10</span>
    <input type="range" id="lambda" min="0.01" max="0.5" step="0.01" value="0.10">
  </label>
  <label>{{lbl_step}} <span id="tval">0.10</span>
    <input type="range" id="step" min="0.01" max="0.19" step="0.01" value="0.10">
  </label>
  <label>{{lbl_iters}} <span id="kval">60</span>
    <input type="range" id="iters" min="10" max="150" step="5" value="60">
  </label>
</div>
<div class="btn-row">
  <button id="run">{{btn_run}}</button>
  <button id="reset" class="ghost">{{btn_reset}}</button>
</div>
<div class="charts">
  <canvas id="objChart" width="280" height="200" title="{{chart_obj_title}}"></canvas>
  <canvas id="sigChart" width="280" height="200" title="{{chart_sig_title}}"></canvas>
</div>
<div id="status" class="status"></div>
/* {{c_css_comment}} */
* { box-sizing: border-box; }
body { font-family: system-ui, sans-serif; color: #222; margin: 0; padding: 4px; }
.controls { display: flex; flex-direction: column; gap: .4rem; margin-bottom: .6rem; }
label { font-size: .85rem; display: flex; align-items: center; gap: .5rem; }
input[type=range] { flex: 1; accent-color: #1d3557; }
.btn-row { display: flex; gap: .5rem; margin-bottom: .7rem; }
button { font: 600 14px system-ui; padding: .4rem .9rem;
         border: 1px solid #1d3557; background: #1d3557; color: #fff;
         border-radius: 8px; cursor: pointer; }
button.ghost { background: #fff; color: #1d3557; }
.charts { display: flex; gap: 8px; flex-wrap: wrap; }
canvas { border: 1px solid #d0d7de; border-radius: 6px; background: #fafbfc; }
.status { font-size: .88rem; font-weight: 600; margin-top: .5rem; min-height: 1.3em; color: #0a7d33; }
// Code not found

Notice how FISTA (orange) consistently reaches low objective values in far fewer steps than ISTA (blue). The gap widens when the problem is poorly conditioned. Soft-thresholding is what drives the sparsity: components smaller than λt\lambda t are zeroed out at every iteration, and the rest are shrunk toward zero by the same amount.

The Real Complexity

How fast do these algorithms converge?

  • ISTA performs, at each iteration kk, two cheap operations: a gradient step xytf(y)x \leftarrow y - t \nabla f(y), then a soft-threshold xisign(xi)max(xiλt,0)x_i \leftarrow \text{sign}(x_i)\max(|x_i| - \lambda t, 0). After kk steps the objective gap satisfies F(xk)FO(1/k)F(x_k) - F^* \leq O(1/k).
  • FISTA adds a momentum coefficient θk\theta_k and an extrapolated point yk+1=xk+θk1θk+1(xkxk1)y_{k+1} = x_k + \tfrac{\theta_k - 1}{\theta_{k+1}}(x_k - x_{k-1}), where θk=1+1+4θk122\theta_k = \tfrac{1 + \sqrt{1 + 4\theta_{k-1}^2}}{2}. This lifts the rate to O(1/k2)O(1/k^2) with no extra gradient evaluations.
  • Optimal in the first-order oracle model. Nesterov's lower bound (1983) shows that no algorithm using only gradient information can do better than O(1/k2)O(1/k^2) for smooth convex objectives. FISTA matches this bound even though gg is nonsmooth, because the proximal step counts as exact minimization of gg.
  • Step size matters. The step size tt must satisfy t1/Lt \leq 1/L, where LL is the Lipschitz constant of f\nabla f. Too large a step causes divergence; too small slows convergence. Line-search variants remove the need to estimate LL.

The elegance of the method is that O(1/k2)O(1/k^2) acceleration comes entirely from reusing the previous iterate — not from second-order information or a more expensive oracle. Compare to interior-point methods for linear programming, which need O(n)O(\sqrt{n}) iterations but each costs O(n3)O(n^3) per step.

Where It Matters

Proximal gradient methods are the workhorses of modern large-scale optimization wherever sparse structure is expected or desired:

  • LASSO regression: the canonical example. Minimizing 12Axb2+λx1\tfrac{1}{2}\|Ax-b\|^2 + \lambda\|x\|_1 selects a sparse subset of features automatically. ISTA/FISTA solve it efficiently even for millions of variables.
  • Compressed sensing and MRI: accelerated MRI acquires far fewer measurements than the Nyquist rate, then reconstructs the image by promoting sparsity in a wavelet domain — solved by FISTA.
  • Image denoising (TV regularization): total-variation denoising is a proximal-gradient problem; the proximal step for the TV norm is a projection onto an \ell_\infty ball, solved by Chambolle's algorithm.
  • Neural network pruning: weight sparsity during training can be encouraged by proximal-gradient updates on the L1 penalty, yielding smaller and faster deployed models.
  • Portfolio optimization: sparse portfolio selection — holding few assets — maps directly to an L1-penalized quadratic, solvable by the same two-step recipe.

The method extends naturally to matrix variables (nuclear-norm minimization for low-rank matrix completion, used in recommendation systems) and to block-separable penalties where each coordinate group has its own proximal operator. Related algorithms like ADMM and the Douglas–Rachford splitting address the same composite structure with different decompositions — see non-convex optimization for what happens when ff loses convexity.

Conclusion

Proximal gradient methods reduce a hard composite problem to a simple two-step recipe: slide down the smooth slope, then apply a closed-form operator that handles the nonsmooth penalty in one shot. For L1 regularization that operator is soft-thresholding — among the simplest formulas in all of optimization.

FISTA shows that momentum is not a heuristic trick but a mathematically certified accelerator: one extra subtraction per step squares the convergence rate from O(1/k)O(1/k) to O(1/k2)O(1/k^2), matching the theoretical lower bound for first-order methods. The result is that compressed sensing, LASSO regression, and image reconstruction problems with millions of variables can now be solved in seconds on a laptop.

The deeper lesson is structural: many practical objectives are smooth-plus-nonsmooth composites, and splitting them along that seam — rather than trying to treat the whole thing with one tool — unlocks both theoretical guarantees and practical speed. That splitting idea ripples through all of modern non-convex optimization and beyond.

Share this article

Pick a channel — or use your device's native share sheet.

Comments

Loading comments...

https://www.kipuhub.com/en/article/proximal-gradient/Content licensed under CC BY-NC 4.0.