Randomized Coordinate Gradient Descent

• Randomized coordinate gradient descent is a first-order optimization method that updates a randomly selected coordinate or block to reduce computational cost.
• It effectively combines partial-gradient updates with techniques like importance sampling and Bregman geometry to address smooth, strongly convex, and nonconvex problems.
• This approach is widely applied in quantum optimization, boosting, and network constraints, offering improved stability and convergence compared to full-gradient methods.

Randomized coordinate gradient descent is a family of first-order methods that minimizes an objective by updating only a randomly selected coordinate or block at each iteration, rather than applying a full-gradient step. In its modern form, it covers smooth problems $\min_x f(x)$, composite objectives of the form $f(x)+\sum_i g_i(x_i)$, and block-structured variables $x=(x_1,\dots,x_n)$; the basic iteration may be a partial-gradient step such as $w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$ or a one-block proximal-gradient subproblem on the sampled block (Wang et al., 2021, Wang et al., 2014, Qu et al., 2014). In several analyses, “randomized coordinate descent” and “randomized coordinate gradient descent” effectively coincide because the method is explicitly gradient-based, using partial derivatives or partial gradients rather than derivative-free coordinate searches (Wang et al., 2014, Wang et al., 2021).

1. Canonical formulation and algorithmic core

The standard setting partitions the variable into coordinates or blocks and updates only one sampled part at a time. For a composite objective with block-separable regularizer,

$$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$

a representative randomized block-coordinate proximal-gradient step is

$$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$

with all other blocks unchanged; in the purely deterministic full-data version, $\nabla_{j_t} f_{i_t}(x^t)$ is replaced by $\nabla_{j_t} f(x^t)$ (Wang et al., 2014). In the smooth unconstrained case, the update reduces to the familiar single-coordinate gradient step

$$\theta_{n+1}=\theta_n-a_n e_{i_n}\,\partial_{i_n}f(\theta_n),$$

with $i_n$ sampled uniformly from $f(x)+\sum_i g_i(x_i)$0 (Ding et al., 2023).

Quadratic least-squares problems provide the classical archetype. For

$f(x)+\sum_i g_i(x_i)$1

randomized coordinate descent chooses column $f(x)+\sum_i g_i(x_i)$2 with probability proportional to $f(x)+\sum_i g_i(x_i)$3 and performs

$f(x)+\sum_i g_i(x_i)$4

which is simultaneously an exact one-dimensional minimization and a coordinate-gradient step scaled by the inverse coordinate curvature (Ramdas, 2014). This column-wise viewpoint is central in linear solvers, where coordinate methods act on features, in contrast with row-based randomized Kaczmarz methods that act on data points (Ramdas, 2014).

The basic rationale is computational. When only one coordinate or block is touched, per-iteration work can be far smaller than in full-gradient methods. This advantage is strongest when gradients are expensive to obtain in all directions at once, when the objective is anisotropic across coordinates, or when data and variable dimensions are both large (Wang et al., 2014, Ding et al., 2023).

2. Sampling models, block smoothness, and coordinate geometry

Classical analyses rely on coordinate-wise smoothness constants. In the single-coordinate smooth setting,

$f(x)+\sum_i g_i(x_i)$5

and one defines

$f(x)+\sum_i g_i(x_i)$6

with the inequalities

$f(x)+\sum_i g_i(x_i)$7

These anisotropic constants explain why coordinate methods can outperform full-gradient descent when the global smoothness constant $f(x)+\sum_i g_i(x_i)$8 is much larger than the average coordinate curvature $f(x)+\sum_i g_i(x_i)$9 (Ding et al., 2023).

A more general framework replaces serial uniform sampling by arbitrary random subsets of blocks. ALPHA formalizes a sampling $x=(x_1,\dots,x_n)$0 with marginal probabilities $x=(x_1,\dots,x_n)$1 and analyzes it through the Expected Separable Overapproximation condition

$x=(x_1,\dots,x_n)$2

which yields a unified theory for serial, minibatch, parallel, distributed, uniform, and importance-sampled coordinate methods (Qu et al., 2014). In the serial case, the ESO vector $x=(x_1,\dots,x_n)$3 reduces to the usual block Lipschitz constants.

Regime	Sampling and update pattern	Representative source
Serial coordinate/block	One coordinate or block per iteration; partial-gradient or proximal step	(Wang et al., 2021, Wang et al., 2014)
Arbitrary or importance sampling	Random subset $x=(x_1,\dots,x_n)$4 with marginals $x=(x_1,\dots,x_n)$5; ESO-based analysis	(Qu et al., 2014)
Online or stochastic block sampling	Random sample or mini-batch index and random coordinate block in the same iteration	(Wang et al., 2014)
Feasibility-preserving constrained sampling	Random pair or path of $x=(x_1,\dots,x_n)$6 blocks so that updates remain feasible	(Necoara et al., 2015, Fan et al., 2017)

Importance sampling follows naturally from this geometry. In ALPHA’s serial smooth setting, the bound suggests the optimal but generally impractical sampling law

$x=(x_1,\dots,x_n)$7

which shows that both curvature and distance-to-solution information matter (Qu et al., 2014). This suggests that coordinate descent is best understood not as a single algorithm, but as a sampling-and-geometry design space.

3. Convergence theory in convex, strongly convex, and nonconvex regimes

For smooth convex problems, randomized coordinate methods admit the same qualitative rate hierarchy as full-gradient methods. ALPHA yields $x=(x_1,\dots,x_n)$8 convergence in the nonaccelerated case and $x=(x_1,\dots,x_n)$9 convergence in the accelerated case under arbitrary sampling, recovering gradient descent, coordinate descent, parallel coordinate descent, and accelerated coordinate descent as special cases (Qu et al., 2014). In block-structured strongly convex settings, linear convergence follows from standard coordinate smoothness and strong convexity assumptions; for quadratic SPD systems, classical randomized coordinate descent has the rate

$w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$0

which makes explicit the dependence on the spectrum of the coordinate geometry (Lok et al., 11 Jun 2025).

A central refinement concerns heterogeneous blocks. Alternating Randomized Block Coordinate Descent updates one randomly sampled block among $w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$1 by a block-gradient step and then exactly minimizes over a distinguished block $w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$2. Its rates depend only on $w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$3 and are independent of the least smooth block $w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$4, provided that block can be optimized exactly; the accelerated version AAR-BCD achieves the corresponding $w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$5 behavior (Diakonikolas et al., 2018). This shows that the special role usually associated with two-block alternating minimization is not limited to the case $w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$6.

For nonconvex problems, the modern theory centers on stationarity rather than objective optimality. In the Bregman-relative-smoothness setting, randomized Bregman coordinate descent defines a full Bregman proximal map $w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$7 and proves the exact characterization

$w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$8

Under blockwise relative smoothness, every limit point is stationary, the expected stationarity measure satisfies

$w_{t+1}=w_t-\eta_t \nabla_{i_t}F_S(w_t)e_{i_t}$9

and, when $$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$0 is strongly convex in Euclidean norm, an $$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$1-stationary point can be found in $$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$2 iterations; the convex rate improves to $$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$3, and relative strong convexity yields global linear convergence (Gao et al., 2020). For structured nonsmooth nonconvex DC objectives $$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$4, randomized coordinate subgradient descent achieves an $$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$5 rate in the expected squared norm of a block composite subgradient mapping, while all limit points are critical (Deng et al., 2019).

Strict-saddle behavior has also been analyzed. Randomized coordinate gradient descent with random coordinate choices and random stepsizes almost surely does not converge to strict saddle points; if all saddle points are strict and critical points are isolated, the iterates converge almost surely to local minima (Chen et al., 2021). A later fixed-stepsize analysis establishes the same qualitative strict-saddle escape statement through a center-stable manifold argument (Chen et al., 11 Aug 2025). This suggests that the saddle-avoidance picture known for full gradient descent extends to genuinely randomized coordinate dynamics.

4. Beyond the Euclidean Lipschitz model

The classical Euclidean model assumes coordinate-wise Lipschitz partial gradients. A major extension replaces that assumption by Bregman geometry and relative smoothness. In the composite problem

$$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$6

the randomized Bregman coordinate descent update is

$$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$7

where $$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$8 is a block Bregman distance (Gao et al., 2020). If one chooses the Euclidean reference $$\min_x\; f(x)+g(x)=\frac{1}{I}\sum_{i=1}^I f_i(x)+\sum_{j=1}^J g_j(x_j),$$9, then $$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$0, relative smoothness becomes ordinary coordinate-wise Lipschitz continuity, and the method reduces exactly to proximal randomized coordinate descent; when $$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$1, it reduces further to standard randomized coordinate gradient descent (Gao et al., 2020).

The gain is that coordinate updates remain analyzable in models where partial gradients are not globally Lipschitz. The paper highlights Poisson inverse problems, KL-type objectives, matrix and tensor factorization, and relative-entropy regression as motivating examples, using non-Euclidean references such as Burg entropy and Boltzmann–Shannon entropy (Gao et al., 2020). This directly contradicts the common assumption that coordinate descent is intrinsically tied to Euclidean quadratic upper models.

A still broader generalization transports the coordinate idea to measure space. Random Wasserstein Coordinate Descent updates a probability measure by restricting the Wasserstein gradient to one ambient axis,

$$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$2

with coordinatewise smoothness constants $$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$3 and total curvature $$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$4 (Xu et al., 2 Apr 2026). The resulting complexity mirrors the Euclidean theory: nonconvex and geodesically convex regimes have $$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$5 complexity, while PL or geodesically strongly convex regimes yield $$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$6 convergence (Xu et al., 2 Apr 2026). This suggests that “coordinate descent” can be abstracted from vector coordinates to geometry-adapted directions in nonlinear spaces.

A different structural extension enriches the direction set itself. Stochastic Spectral and Conjugate Descent shows that merely changing coordinate sampling probabilities is fundamentally limited on SPD quadratics, whereas augmenting coordinates with a few spectral or conjugate directions can interpolate between the classical randomized coordinate rate and a condition-number-independent $$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$7 regime (Kovalev et al., 2018). The paper’s SSCD scheme removes the influence of the smallest eigenvalues by sampling both coordinates and selected eigendirections. This is no longer pure coordinate descent, but it arises directly from an analysis of what classical randomized coordinate gradient descent cannot do by reweighting coordinates alone.

5. Stochastic, online, constrained, and linear-system variants

When both the number of samples and the number of features are large, coordinate randomization can be combined with stochastic sampling over data. ORBCD samples a mini-batch or data index $$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$8, then samples a coordinate block $$x_{j_t}^{t+1}=\arg\min_{x_{j_t}} \left\langle \nabla_{j_t} f_{i_t}(x^t),x_{j_t}\right\rangle + g_{j_t}(x_{j_t}) + \frac{\eta_t}{2}\|x_{j_t}-x_{j_t}^t\|_2^2,$$9, and performs a proximal step using only $\nabla_{j_t} f_{i_t}(x^t)$0 (Wang et al., 2014). In the convex stochastic case it achieves an $\nabla_{j_t} f_{i_t}(x^t)$1 rate, in the strongly convex case an $\nabla_{j_t} f_{i_t}(x^t)$2 rate, and with SVRG-style variance reduction—ORBCDVD—it converges geometrically in expectation for strongly convex objectives (Wang et al., 2014). This places randomized coordinate gradient descent squarely at the intersection of SGD and block-coordinate proximal methods.

Linear coupling constraints require a different update geometry, because a single-block move typically breaks feasibility. For network problems with the constraint $\nabla_{j_t} f_{i_t}(x^t)$3, randomized block coordinate descent therefore updates $\nabla_{j_t} f_{i_t}(x^t)$4 variables at a time along a sampled path or edge, maintaining the zero-sum constraint exactly (Necoara et al., 2015). The expected complexity is $\nabla_{j_t} f_{i_t}(x^t)$5 for smooth convex objectives, with linear convergence in the strongly convex case (Necoara et al., 2015). A later refinement for the two-block case sharpened the expected rate by extending Nesterov’s second technique to the linearly constrained setting and by constructing a basis for the feasible subspace (Fan et al., 2017).

Linear systems and least squares provide another major specialization. For overdetermined inconsistent systems, randomized coordinate descent on the columns of $\nabla_{j_t} f_{i_t}(x^t)$6 converges linearly to the least-squares solution, while randomized Kaczmarz on rows does not; for underdetermined consistent systems, the preference reverses if one wants the minimum-norm solution (Ramdas, 2014). More recent work modifies classical coordinate descent rather than the problem class. Subspace-constrained randomized coordinate descent restricts the iterates to an affine space induced by a Nyström approximation and obtains a rate controlled by the tail spectrum $\nabla_{j_t} f_{i_t}(x^t)$7, thereby avoiding degradation from large spectral outliers (Lok et al., 11 Jun 2025). Randomized conjugate gradient least squares goes further by replacing the full gradient in a CGLS-type direction-correction scheme with a randomized coordinate gradient $\nabla_{j_t} f_{i_t}(x^t)$8; the resulting RCGLS method has a linear convergence bound at least as good as standard randomized coordinate descent and admits a lightweight sparse implementation (Zeng et al., 24 May 2026).

These developments show that “randomized coordinate gradient descent” is not a single update rule but a design principle: choose a low-cost partial first-order surrogate, preserve whatever structure the problem requires, and analyze the resulting dynamics in the natural geometry of the model.

6. Applications, statistical interpretation, and common misconceptions

The application range is broad because many modern objectives make full gradients expensive or unnatural. In parameterized quantum circuits, each partial derivative is typically obtained by a separate parameter-shift measurement, so the full gradient costs roughly $\nabla_{j_t} f_{i_t}(x^t)$9 times one derivative evaluation; under a local PL condition and bounded unbiased noise, random coordinate descent has the same stochastic stability order as gradient descent and can be much cheaper in total partial-derivative cost when $\nabla_{j_t} f(x^t)$0 (Ding et al., 2023). The reported experiments include nearly $\nabla_{j_t} f(x^t)$1 fewer partial derivative evaluations than gradient descent on a transverse-field Ising VQE instance and about $\nabla_{j_t} f(x^t)$2 fewer on a distinct $\nabla_{j_t} f(x^t)$3-qubit, $\nabla_{j_t} f(x^t)$4-parameter ansatz (Ding et al., 2023).

Boosting provides a different reinterpretation. Randomized Gradient Boosting Machine views weak-learner coefficients as coordinates in a convex optimization problem and shows that classical GBM is greedy coordinate descent, while RGBM is a random-then-greedy hybrid that samples a subset of weak learners and then chooses the largest-magnitude coordinate gradient within that subset (Lu et al., 2018). This does not fit the plain “sample one coordinate and update it immediately” template, but it demonstrates how randomized coordinate-gradient ideas extend to dictionary search in supervised learning.

The statistical role of coordinate randomization is not limited to optimization speed. Stability and Generalization for Randomized Coordinate Descent proves that, in convex empirical-risk minimization, the $\nabla_{j_t} f(x^t)$5 on-average argument stability scales as

$\nabla_{j_t} f(x^t)$6

and in the strongly convex case as $\nabla_{j_t} f(x^t)$7, leading to $\nabla_{j_t} f(x^t)$8 and $\nabla_{j_t} f(x^t)$9 excess-risk bounds with suitable early stopping (Wang et al., 2021). The same work argues that RCD is more stable than SGD because only one coordinate is perturbed at each step (Wang et al., 2021). This reframes randomized coordinate gradient descent as a learning algorithm with distinctive generalization behavior, not merely a computational primitive.

Several recurrent misconceptions are clarified by the literature. First, randomized coordinate descent is not inherently derivative-free: in ORBCD, ALPHA, and the stability analysis, it is explicitly a partial-gradient or partial proximal-gradient method (Wang et al., 2014, Qu et al., 2014, Wang et al., 2021). Second, lower per-iteration cost does not by itself imply lower total complexity. In Langevin Monte Carlo, blindly replacing the full gradient by a one-coordinate unbiased estimator produces no computational gain, because the variance increases the required number of iterations; only after incorporating variance reduction through SAGA or SVRG does the coordinate strategy improve total cost, and in the underdamped case it can match the full-gradient iteration count while using one directional derivative per step (Ding et al., 2020). Third, the usual coordinate-wise Lipschitz assumption is not the only viable analytical framework: Bregman relative smoothness and Wasserstein coordinate smoothness show that coordinate methods can be formulated in genuinely non-Euclidean geometries (Gao et al., 2020, Xu et al., 2 Apr 2026).

Taken together, these results depict randomized coordinate gradient descent as a large methodological family unified by one principle: replace a full first-order update by a randomized partial one that is cheaper, structurally compatible with the model, and analyzable in an appropriate geometry. Its modern theory spans convex and nonconvex optimization, exact and proximal updates, Euclidean and Bregman geometries, stochastic and online learning, constrained network problems, linear solvers, boosting, quantum optimization, and even Wasserstein spaces.